# Improved Firefly Algorithm: A Novel Method for Optimal Operation of Thermal Generating Units

Authors: Thang Trung Nguyen; Nguyen Vu Quynh; Le Van Dai
Journal: Complexity
 (2018)

---

## Abstract

This paper presents a novel improved firefly algorithm (IFA) to deal the problem of the optimal operation of thermal generating units (OOTGU) with the purpose of reducing the total electricity generation fuel cost. The proposed IFA is developed based on combining three improvements. The first is to be based on the radius between two solutions, the second is updated step size for each considered solution based on different new equations, and the third is to slightly modify a formula producing new solutions by using normally distributed random numbers and canceling uniform random numbers of conventional firefly algorithm (FA). The effect of each proposed improvement on IFA is investigated by executing five benchmark functions and two different systems. The performance of IFA is investigated on six other study cases consisting of different types of objective function and complex level of constraints. The objective function considers single fuel with quadratic form and nonconvex form, and multifuels with the sum of several quadratic and nonconvex functions while a set of constraints taken into account are power loss, prohibited zone, ramp rate limit, spinning reserve, and all constraints in transmission power networks. The obtained results indicate the proposed improvements in terms of high optimal solution quality, stabilization of search ability, and fast convergence compared with FA. In addition, the comparisons with other methods also lead to a conclusion that the proposed method is a very promising optimization tool for systems with quadratic fuel cost function and with complicated constraints.


---

## Body

## 1. Introduction

The optimal operation of thermal generating units (OOTGU) has been widely concerned in the power system operation field due to its significant important role. In fact, thermal units are using a huge amount of fossil fuel for generating electricity while primary fuels are expensive and will be exhausted in the near future. The objective of the OOTGU problem is to determine the generated active power of thermal generating units so that the total fuel cost can be reduced as much as possible; while all the constraints of thermal generating units and of transmission power networks are taken into account. The problem becomes more and more complicated since more and more complicated constraints of generator are considered such as ramp rate, prohibited zone, voltage limitations, and generation limits. In addition, constraints related to transmission power network are not simple conditions to be exactly met such as the active power spinning reserve, active and reactive power balance, limitation of transmission lines, limitations of transformer tap, and limitation of capacitor bank’s reactive power.There have been a huge number of studies considering OOTGU as a main problem. The studies have applied conventional methods such as deterministic algorithms and new methods. Most deterministic methods are based on gradient, derivative, Lagrange optimization function while the new methods are based on metaheuristic algorithms, modified versions of metaheuristic algorithms, and combination of two different metaheuristic algorithms. In addition, the metaheuristic methods integrated to other methods dealing with constraints of problem are also effectively combined optimization tools.Deterministic methods have been applied for the problem such as Lagrange optimization function-based method (LM) [1–3], dynamic programming (DP) [1], linear programming method (LPM) [4], Hierarchical method (HM) [5], Hopfield neural networks (HNN) [3, 6], improved Hopfield neural network method (IHNNM) [7], augmented Lagrange Hopfield network (ALHN) [8], and enhanced augmented Lagrange Hopfield network (EALHN) [9, 10]. Among the methods, LM was the first optimization tool for solving the problem of OOTGU with quadratic fuel cost functions and the balance constraint of generated power and required power such as load demand and power losses. The obtained results from the method were acceptable since it could deal with the constraint of active power balance constraint and maximum error was very small. DP also achieved results as good as those from LM. However, DP required a huge number of iterations for search process, especially for large-scale systems with large number of units. LPM and HM have applied approximation techniques to simplify the complex of constraints and nonsmooth fuel cost functions and then LM was acted as an optimization search algorithm. Consequently, the applicability of the two methods was restricted for complicated systems with many nonlinear constraints and nonsmooth fuel cost functions. HNN and its improved versions in [7–9] consisting of IHNNM, ALHN, and EALHN have become more effective than LM, DP, LPM, and HM in dealing with large-scale systems with faster search process and better optimal solutions. The method group is mainly based on energy function and Hopfield neural network. However, HNN variants did not use Lagrange optimization function but they established energy function at the beginning and then Hopfield neural network was applied. On the contrary, ALHN variants have constructed an augmented Lagrange optimization function first and then the function was converted into energy function with the present of inverse sigmoid function. The ALHN variants could overcome several drawbacks of HNN variants such as high oscillation, high error, and low optimal solutions for complicated systems. In general, these methods in the deterministic group have the same limitations of application for systems considering nonconvex fuel cost function and nondifferentiable functions.New methods based on metaheuristic algorithms have been widely and successfully applied for the considered OOTGU problem even for very complicated systems with nonconvex fuel cost functions and constraints taking prohibited zones, active power spinning reserve, and ramp rate into account. Differential evolution variants have been applied for solving the problem in [11–13] consisting of conventional differential evolution (DE) [11], hybrid dynamic programming integer-coded (HDP-ICDE) [12], and colonial competitive differential evolution (CCDE) [13]. DE together with GA and PSO has been implemented for different cases with different constraints with the superiority of DE over genetic algorithm (GA) and particle swarm optimization (PSO) as a result. HDP-ICDE has used integer-coded differential evolution (ICDE) as an optimization tool for searching solutions meanwhile dynamic programming has been used as an evaluation tool for computing fitness function value of found solutions. The method could not overcome drawbacks of ICDE such as high number of solution evaluation and falling into local search. CCDE has changed the mutation operation by using different models based on so-far best solutions and replacement of worst solutions by more promising solutions that aim to improve the balance between local search and global search as well as to retain the best solutions. This method could show its superiority over DE and other versions of DE, but the test systems have been used and the result comparisons from the tests could not show its better performance than other methods because its cost and other methods’ cost were equal. The comparisons of search speed have not been carried out. The performance of cuckoo search algorithm (CSA) [14, 15] and improved cuckoo search algorithm (ICSA) [16] have been shown clearly since different types of constraints were considered and wide comparisons were implemented. However, the speed comparisons have been made only via execution time, leading to not good comparison criteria for conclusion. Different test systems considering constraints of transmission power network have been used to test the efficiency of PSO variants in [17–19]. PSO with pseudo-gradient and constriction factor (PG-CF-PSO) [17] has determined the best direction for updating velocity by using pseudo-gradient theory and used constriction factor to find good search zone. It has been demonstrated to be more effective than PSO via three IEEE systems with 30, 57, and 118 buses. New adaptive particle swarm optimization (NAPSO) [18] has integrated both mutation operation and adaptive PSO for avoiding local optimal solution. The method has used different fuzzy and self-adaptive techniques for updating parameter and solutions. Each proposed technique has been investigated by testing different PSO methods with each modification such as NAPSO without mutation (NAPSO1), NAPSO without fuzzy (NAPSO2), and NAPSO without self-adaptive (NAPSO3). These methods have had worse results than NAPSO but better results than standard PSO. Combination of particle swarm optimization and tabu search algorithm (PSO-TSA) [19] has used the ability of PSO for local search based on three different phases and TSA for global search by adjusting obtained solutions of PSO. The method has been compared to three difference implemented methods such as GA, PSO, and TSA, and it can be shown better objective function. In spite of the superiority over the three methods, PSO-TSA has been more complicated, included more control parameters, and used high number of produced new solutions. Different GA variants have been proposed for dealing with different set of constraints of the considered problem such as real-coded genetic algorithm (RCGA) [20], RCGA with arithmetic-average-bound crossover and wavelet (NRCGA) [20], hybrid real-coded genetic algorithm (HRCGA) [21], improved GA (IGA) [22], and IGA with updated multiplier (IGAMU) [22]. These methods have been developed by improving GA, using RCGA and other improved crossover operation, and combining improved version of GA and multiplier updating technique. Generally, more complicated modifications lead to more improvement of GA variants. However, there have been more drawbacks such as higher number of control parameters, difficult task of selection of such factors, and complicated implementation when combining two or more methods. The different versions of krill herd algorithm (KHA) [23] have been used for solving different cases of the considered problem such as KHA without using genetic operation (KHA1), KHA using crossover operation (KHA2), KHA using mutation operation (KHA3), and KHA using crossover operation and mutation operation (KHA4). Several disadvantages of these KHA methods have been pointed out in [24] such as low convergence speed and low success rate. Thus, opposition-based krill herd algorithm (OBKHA) [24] has been developed by applying learning technique based on opposition and it has been compared to KHA methods. As derived from the comparisons, OBKHA could be considered better than KHA methods about optimal solutions reflected via less minimum cost and stable search ability reflected via less average cost, less maximum cost, and less standard deviation; however, the convergence speed and the success rate have not been investigated since there were no comparisons of used population and used iterations. Biogeography-based optimization algorithm (BBOA) has been employed for minimizing cost and emission in [25] and for minimizing cost in [26]. BBOA has some characteristics in common with other metaheuristic algorithms consisting of mutation and selection like GA, DE, and evolution programming, but this method has much more number of control parameters such as the probability of habitat modification and mutation, the rate of mutation and immigration, and lower and upper boundaries of probability. The results reported in [26] could show good performance of BBOA via comparisons with several existing methods on four test cases, but the keywords for determining the best control parameters were missed. Clearly, the application of BBOA to real power system with complicated constraints in transmission power network was limited and studies about BBOA for complicated systems have not been carried so far. Grey wolf optimization algorithm (GWOA) [27, 28] has been applied for simple systems without complicated constraints and small number of units. Thus, the demonstration of the methods’ performance has not been persuasive. On the contrary, crisscross optimization algorithm (CSOA) [29] and exchange market algorithm (EMA) [30] have had more result comparisons in the implementation of large-scale systems and complicated constraints. The ant lion optimization algorithm (ALOA) [31] was also a new method with few applications, and the method in the study has not shown its potential search ability persuasively due to simple test employment. In addition, many new methods have been developed for solving the problem such as modified symbiotic organisms search algorithm (MSOSA) [32], mine blast algorithm (MBA) [33], clonal algorithm (CA) [34], mathematical programming algorithm (MPA) [35], improved quantum-inspired evolutionary algorithm (IQIEA) [36], cuckoo optimization algorithm (COA) [37], improved colliding bodies optimization algorithm (ICBOA) [38], flower pollination algorithm (FPA) [39], natural updated harmony search (NUHS) [40], lightning flash algorithm (LFA) [41, 42], moth swarm algorithm (MSA) [43], and orthogonal learning competitive swarm optimization algorithm (OLCSOA) [44]. Among the remaining methods, LFA, MSA, ICBOA, and OLCSOA were developed in the last two years and they were considered better than most existing methods for considered test cases. LFA and OLCSOA have been implemented for systems with different complex levels of thermal generating units such as nonsmooth fuel cost functions with valve-point loading effects, multifuel types, and prohibited zone meanwhile challenges of transmission power network constraints have not been taken into account. On the contrary, MSA has been focused on dealing with the constraint complex rather than the challenges of fuel cost function of thermal generating units. In general, metaheuristic algorithms have been widely and successfully applied for different test systems and their superiorities have been pointed out mainly based on fuel cost function comparisons. The performance of ICBOA has been measured by testing on different study cases with different types of fuel cost function and the consideration of all constraints of transmission power networks. ICBOA together with other methods such as artificial bee colony (ABC), DE, PSO, BBO, standard colliding bodies optimization algorithm (CBO), and enhanced CBO have been measured by testing on different study cases with different types of fuel cost function and the consideration of all constraints of transmission power networks. ICBO has been pointed out obtaining better optimal solutions than all methods, but the further comparisons of convergence speed have not been carried out and discussed.The firefly algorithm (FA) is also a population-based metaheuristic algorithm similar to PSO, DE, CSA, and so on built by Yang in 2008 for solving optimization problems [45]. The configuration of FA consists of three procedures of updating the distance between two considered fireflies, updating a step size, and updating new solutions. FA has shown its superiority over other traditional algorithms consisting of GA and PSO [46] by the comparisons of different benchmark optimization functions. In 2011, Yang and his coworkers applied the algorithm for solving the considered OOTGU problem [47] and they have stated that the method has been efficient for the problem by carrying out result comparisons via testing on three different systems. However, the performance of the method has not been demonstrated for large-scale systems with large number of generators and highly complicated constraints in transmission networks. Thus, there were several improvement versions of the method for dealing with the problem and the comparisons with it were also done for evaluation. An improved firefly algorithm (IFA) has been proposed by Kazemzadeh-Parsi in [48] by applying three different modifications. The first modification was to use the better solutions of the previous iteration for replacing the worst solutions in the current iteration. The second modification was to replace low solutions by randomly produced solutions, and the third modification was to move the worse solutions only to one representative solution instead of all better solutions. The representative solution is determined by using the average solution of all better solutions. The memetic firefly algorithm (MFA) in [49] has focused on the balance of exploration acting as global search and exploitation acting as local search by using adaptive attractiveness β and adaptive step parameter α with respect to the change of iteration in formula updating new solutions for each lower quality solution. Another improved FA with using adaptive step parameter α (ASPFA) was proposed in [50] by suggesting an adaptive formula for updating step parameter α based on current iteration and the maximum number of iterations. The chaotic firefly algorithm (CFA) was developed in [51] for solving the considered OOTGU problem with different test cases of fuel cost function and constraints. The method has used chaotic distribution to produce the values of attractiveness β and step parameter α instead of using random distribution in FA. Moustafa et al. [52] have implemented conventional FA, IFA, MFA, and ASPFA on 3-unit system and 6-unit system for comparison. The result comparisons have indicated that IFA with three modifications was the best one and FA was the worst one among four FA variants. However, there was no comparison between these methods with others. On the contrary, CFA in [51] has been tested on many tests with more complicated objective functions considering valve-point effects and multifuels, and more complicated constraints considering POZ, ramp rate, and active power spinning reserve. CFA has been demonstrated to be more effective than many methods such as PSO, ABC, DE, and GA variants. Another modified firefly algorithm (α-MFA) based on the modification on distance between two fireflies and the modification on generation of step parameter α was proposed in [53] for dealing with OOTGU problem with 3-, 6-, 13-, and 15-unit systems considering valve-point loading effects, ramp rate constraints, and POZ constraints. The five improved versions of FA have been demonstrated to be better than conventional FA and other methods. However, the superiority of the methods was not shown clearly since it has not been implemented for complicated systems taking all constraints of transmission power networks and generators into consideration.In the paper, we propose three modifications on the FA in order to tackle several disadvantages of FA such as premature convergence to local optimum solution and impossibility of jumping out of the search zone with many local optimum solutions. In the first modification, we propose a new formula to update radius between a considered fireflyXi (a solution) to another firefly Xj (another solution) with lower fitness function than the considered solution. The proposed radius based on Xi and the best solution XGbest become more effective better than that based on Xi and Xj used in FA. In the second modification, we propose a new algorithm for producing new solutions of an old solution by suggesting two models for the updated step size. A large or a smaller updated step size will be used to diversify new solutions and to avoid converging to a local optimum and trapping into search zone with many local optimums. In the third modification, uniform distribution is replaced with normal distribution aiming at diversifying the search zone. As a result, the improvement of the new algorithm is highly considerable compared to FA. The application of each modification and three modifications is evaluated by testing on difference of four systems with other six cases. These test cases consider thermal generating units using single fuel and multifuels, convex and nonconvex objective function, and all constraints in real power systems. The main contributions of the paper can be summarized as follows:
(i)
Carry out three modifications on the new solution produced phase of conventional FA and point out the specific impact of each modification on the obtained results.
(ii)
The proposed method has few control parameters, and the parameters are easily tuned. The proposed method is more effective than FA, but the implementation of the proposed method is as simple as that of FA.
(iii)
Thoroughly present the application of the proposed IFA method for solving different systems with different constraints in optimal operation of power systems. The exact selection of control variables can help the proposed method satisfy all considered constraints.
(iv)
Develop basic application of the proposed method for optimization problem and lead to conclusion that the proposed method is effective for systems with complicated constraints but it works less effectively for systems considering nonconvex objective functions.The rest of the paper is organized as follows. Section2 describes the problem formulation. Section 3 presents the firefly algorithms. Section 4 shows the implementation of IFA for OOTGU problem. Section 5 presents the results of numerical simulations. Finally, the conclusions and future works are given in Section 6.


## 2. Problem Formulation

### 2.1. Objective Function

In OOTGU with single fuel, the fuel cost of each generating unit is expressed as a quadratic function of its power output. The objective of the problem is to minimize the total fuel cost ofN available units and can be written as follows [1]:
(1)minF=∑s=1NFsPs.where Fs is the fuel cost function of thermal unit s and it can be represented in a single quadratic form or piecewise form depending on the number of fuel cost options that thermal plants use. If the thermal units use only one fuel, its form is the second order equation as shown in (2) while the form of Fs when using multiple fuels is a piecewise curve consisting of at least two second order equations. For the case of multifuel options, Fs is shown in (3) [1].
(2)FsPs=as+bsPs+csPs2,(3)FsPs=as1+bs1Ps+cs1Ps2,fuel1,Ps,min≤Ps≤Ps1,max,asm+bsmPs+csmPs2,fuelm,Psm,min≤Ps≤Psm,max,asMs+bsMsPs+csMsPs2,fuelMs,PsMs,min≤Ps≤Ps,max.In addition, once valve-point loading effects (VPLE) of thermal units are considered during the process of increase or decrease of power output, (2) and (3) will become more complicated models as follows [20].
(4)FsPs=as+bsPs+csPs2+ds×sinesPs,min−Ps.FiPi=as1+bs1Ps+cs1Ps2+ds1×sines1Ps,min−Ps,fuel1,Ps,min≤Ps≤Ps1,max,asm+bsmPs+csmPs2+dsm×sinesmPsm,min−Ps,fuelm,Psm,min≤Ps≤Psm,max,asMS+bsMSPs+csMSPs2+dsMS×sinesMSPsMS,min−Ps,fuelMS,PsMS,min≤Ps≤PsMS,max.For easy understanding, the case of single fuel option and the case of multifuel options neglecting and considering VPLE are depicted in Figures1 and 2. It is clear that the curves with VPLE are much more complicated.Figure 1
Fuel cost function for the case of single fuel option.Figure 2
Fuel cost function for the case of multifuel options.


### 2.2. Set of Constraints

#### 2.2.1. Generator Constraints

Active and reactive power and voltage of generators are restricted by the following inequality [17].
(5)Ps,min≤Ps≤Ps,max,s=1,…,N,(6)Qs,min≤Qs≤Qs,max,s=1,…,N,(7)Vs,min≤Vs≤Vs,max,s=1,…,N,where Qs,min and Qs,max are the lower and upper limitations of reactive power of generator s, respectively; Vs,min and Vs,max are the lower and upper limitations of voltage of generator s, respectively; and Qs and Vs are the operating values of reactive power and voltage of generator s, respectively.It is noted that the constraint (5) is always taken into account in the OOTGU problem while the constraints (6) and (7) are considered only when real transmission power networks are being considered.


#### 2.2.2. Prohibited Operating Zone (POZ) Constraints

In addition to the active power limitation constraint in (5), thermal units are also constrained by prohibited operating zones due to the security requirement. Thus, power output of each thermal unit should be within its lower limit and upper limit and outside the POZ. An example with two prohibited operating zones of thermal units is shown in Figure 3, and a typical thermal unit with prohibited operating zones is mathematically modeled as follows [20].
(8)Ps∈Ps,min≤Ps≤Ps,1lower,Ps,hupper≤Ps≤Ps,h+1lower,h=1,…,NPOZ,s,Ps,NPOZ,supper≤Ps≤Ps,max.Figure 3
Prohibited operating zones of thermal units.


#### 2.2.3. Spinning Reserve Constraint

In order to avoid lack of energy in case that the largest power unit shutdowns, total active power reserve of all thermal units should be required to be equal to or higher than the largest power unit. The requirement is as below [20].
(9)∑s=1NPSRs≥PSR.where PSR is the total active power that the power systems require for spinning reserve and PSRs is the active power spinning reserve of thermal generating unit s.


#### 2.2.4. Ramp Rate Constraint

During the operation of thermal units, generation increase and generation decrease are constrained by the ramp up and down the rate limit, respectively. Thus, initial power output of units and the ramp rate limits will determine the lowest and the highest power output of unit at the moment. The description can be presented in the following inequalities [20].
(10)Ps−Ps0≤RURs,ifhighergenerationisrequired,Ps,0−Ps≤RDRs,iflowergenerationisrequired.


#### 2.2.5. Real Power Balance

The total real power output of generating units satisfies total load demand plus system power losses [1].
(11)∑s=1NPs=PLD+PTL,and the total power loss is calculated using Kron’s (1):
(12)PLD=∑s=1N∑h=1NPsBshPh+∑s=1NB0sPs+B00.In addition to (11), real power balance is also presented by the following equation since real transmission power networks are considered [17].
(13)Pi−PLDi=Vi∑j=1NbVjGijcosδi−δj+Bijsinδi−δj,i=1,…,Nbus,where Gij and Bij are the real part and the imaginary part of admittance of a transmission line connecting buses i and j.Similarly, reactive power balance is also modeled as follows [17]:
(14)Qi+Qci−QLDi=Vi∑j=1NbVjGijsinδi−δj−Bijcosδi−δj,i=1,…,Nbus,where Qci is the generated reactive power of capacitors at bus i and constrained by the following inequality [17]:
(15)Qci,min≤Qci≤Qci,max,i=1,…,Nc,where Qci,min and Qci,max are, respectively, the lowest and the highest reactive power of the capacitor bank at bus i.


#### 2.2.6. Transformer Tap Constraints

In a real transmission power network, transformers are usually used to step down voltage and supply electricity to lower voltage transmission lines. The voltage of the lower voltage lines can be adjusted by setting the tap of the transformers, so that the following rule is always exactly met [17].
(16)Ti,min≤Ti≤Ti,max,i=1,…,Ntap,where Ti,min and Ti,max are the lower limit and the upper limit of transformer’s tap at bus i.


#### 2.2.7. Stably Working Power Network Constraints

In order to operate electric power network continuously, load bus voltage and apparent power flow of transmission lines must follow the following limitations [17].
(17)VLi,min≤VLi≤VLi,max,i=1,…,Nload,MVAi≤MVAi,max,i=1,…,Nline,where VLi,min and VLi,max are, respectively, lower voltage and upper voltage of load at bus i; MVAi,max is the apparent power flow capacity of line i; and MVAi is the current apparent power flow of line i.


## 2.1. Objective Function

In OOTGU with single fuel, the fuel cost of each generating unit is expressed as a quadratic function of its power output. The objective of the problem is to minimize the total fuel cost ofN available units and can be written as follows [1]:
(1)minF=∑s=1NFsPs.where Fs is the fuel cost function of thermal unit s and it can be represented in a single quadratic form or piecewise form depending on the number of fuel cost options that thermal plants use. If the thermal units use only one fuel, its form is the second order equation as shown in (2) while the form of Fs when using multiple fuels is a piecewise curve consisting of at least two second order equations. For the case of multifuel options, Fs is shown in (3) [1].
(2)FsPs=as+bsPs+csPs2,(3)FsPs=as1+bs1Ps+cs1Ps2,fuel1,Ps,min≤Ps≤Ps1,max,asm+bsmPs+csmPs2,fuelm,Psm,min≤Ps≤Psm,max,asMs+bsMsPs+csMsPs2,fuelMs,PsMs,min≤Ps≤Ps,max.In addition, once valve-point loading effects (VPLE) of thermal units are considered during the process of increase or decrease of power output, (2) and (3) will become more complicated models as follows [20].
(4)FsPs=as+bsPs+csPs2+ds×sinesPs,min−Ps.FiPi=as1+bs1Ps+cs1Ps2+ds1×sines1Ps,min−Ps,fuel1,Ps,min≤Ps≤Ps1,max,asm+bsmPs+csmPs2+dsm×sinesmPsm,min−Ps,fuelm,Psm,min≤Ps≤Psm,max,asMS+bsMSPs+csMSPs2+dsMS×sinesMSPsMS,min−Ps,fuelMS,PsMS,min≤Ps≤PsMS,max.For easy understanding, the case of single fuel option and the case of multifuel options neglecting and considering VPLE are depicted in Figures1 and 2. It is clear that the curves with VPLE are much more complicated.Figure 1
Fuel cost function for the case of single fuel option.Figure 2
Fuel cost function for the case of multifuel options.


## 2.2. Set of Constraints

### 2.2.1. Generator Constraints

Active and reactive power and voltage of generators are restricted by the following inequality [17].
(5)Ps,min≤Ps≤Ps,max,s=1,…,N,(6)Qs,min≤Qs≤Qs,max,s=1,…,N,(7)Vs,min≤Vs≤Vs,max,s=1,…,N,where Qs,min and Qs,max are the lower and upper limitations of reactive power of generator s, respectively; Vs,min and Vs,max are the lower and upper limitations of voltage of generator s, respectively; and Qs and Vs are the operating values of reactive power and voltage of generator s, respectively.It is noted that the constraint (5) is always taken into account in the OOTGU problem while the constraints (6) and (7) are considered only when real transmission power networks are being considered.


### 2.2.2. Prohibited Operating Zone (POZ) Constraints

In addition to the active power limitation constraint in (5), thermal units are also constrained by prohibited operating zones due to the security requirement. Thus, power output of each thermal unit should be within its lower limit and upper limit and outside the POZ. An example with two prohibited operating zones of thermal units is shown in Figure 3, and a typical thermal unit with prohibited operating zones is mathematically modeled as follows [20].
(8)Ps∈Ps,min≤Ps≤Ps,1lower,Ps,hupper≤Ps≤Ps,h+1lower,h=1,…,NPOZ,s,Ps,NPOZ,supper≤Ps≤Ps,max.Figure 3
Prohibited operating zones of thermal units.


### 2.2.3. Spinning Reserve Constraint

In order to avoid lack of energy in case that the largest power unit shutdowns, total active power reserve of all thermal units should be required to be equal to or higher than the largest power unit. The requirement is as below [20].
(9)∑s=1NPSRs≥PSR.where PSR is the total active power that the power systems require for spinning reserve and PSRs is the active power spinning reserve of thermal generating unit s.


### 2.2.4. Ramp Rate Constraint

During the operation of thermal units, generation increase and generation decrease are constrained by the ramp up and down the rate limit, respectively. Thus, initial power output of units and the ramp rate limits will determine the lowest and the highest power output of unit at the moment. The description can be presented in the following inequalities [20].
(10)Ps−Ps0≤RURs,ifhighergenerationisrequired,Ps,0−Ps≤RDRs,iflowergenerationisrequired.


### 2.2.5. Real Power Balance

The total real power output of generating units satisfies total load demand plus system power losses [1].
(11)∑s=1NPs=PLD+PTL,and the total power loss is calculated using Kron’s (1):
(12)PLD=∑s=1N∑h=1NPsBshPh+∑s=1NB0sPs+B00.In addition to (11), real power balance is also presented by the following equation since real transmission power networks are considered [17].
(13)Pi−PLDi=Vi∑j=1NbVjGijcosδi−δj+Bijsinδi−δj,i=1,…,Nbus,where Gij and Bij are the real part and the imaginary part of admittance of a transmission line connecting buses i and j.Similarly, reactive power balance is also modeled as follows [17]:
(14)Qi+Qci−QLDi=Vi∑j=1NbVjGijsinδi−δj−Bijcosδi−δj,i=1,…,Nbus,where Qci is the generated reactive power of capacitors at bus i and constrained by the following inequality [17]:
(15)Qci,min≤Qci≤Qci,max,i=1,…,Nc,where Qci,min and Qci,max are, respectively, the lowest and the highest reactive power of the capacitor bank at bus i.


### 2.2.6. Transformer Tap Constraints

In a real transmission power network, transformers are usually used to step down voltage and supply electricity to lower voltage transmission lines. The voltage of the lower voltage lines can be adjusted by setting the tap of the transformers, so that the following rule is always exactly met [17].
(16)Ti,min≤Ti≤Ti,max,i=1,…,Ntap,where Ti,min and Ti,max are the lower limit and the upper limit of transformer’s tap at bus i.


### 2.2.7. Stably Working Power Network Constraints

In order to operate electric power network continuously, load bus voltage and apparent power flow of transmission lines must follow the following limitations [17].
(17)VLi,min≤VLi≤VLi,max,i=1,…,Nload,MVAi≤MVAi,max,i=1,…,Nline,where VLi,min and VLi,max are, respectively, lower voltage and upper voltage of load at bus i; MVAi,max is the apparent power flow capacity of line i; and MVAi is the current apparent power flow of line i.


## 2.2.1. Generator Constraints

Active and reactive power and voltage of generators are restricted by the following inequality [17].
(5)Ps,min≤Ps≤Ps,max,s=1,…,N,(6)Qs,min≤Qs≤Qs,max,s=1,…,N,(7)Vs,min≤Vs≤Vs,max,s=1,…,N,where Qs,min and Qs,max are the lower and upper limitations of reactive power of generator s, respectively; Vs,min and Vs,max are the lower and upper limitations of voltage of generator s, respectively; and Qs and Vs are the operating values of reactive power and voltage of generator s, respectively.It is noted that the constraint (5) is always taken into account in the OOTGU problem while the constraints (6) and (7) are considered only when real transmission power networks are being considered.


## 2.2.2. Prohibited Operating Zone (POZ) Constraints

In addition to the active power limitation constraint in (5), thermal units are also constrained by prohibited operating zones due to the security requirement. Thus, power output of each thermal unit should be within its lower limit and upper limit and outside the POZ. An example with two prohibited operating zones of thermal units is shown in Figure 3, and a typical thermal unit with prohibited operating zones is mathematically modeled as follows [20].
(8)Ps∈Ps,min≤Ps≤Ps,1lower,Ps,hupper≤Ps≤Ps,h+1lower,h=1,…,NPOZ,s,Ps,NPOZ,supper≤Ps≤Ps,max.Figure 3
Prohibited operating zones of thermal units.


## 2.2.3. Spinning Reserve Constraint

In order to avoid lack of energy in case that the largest power unit shutdowns, total active power reserve of all thermal units should be required to be equal to or higher than the largest power unit. The requirement is as below [20].
(9)∑s=1NPSRs≥PSR.where PSR is the total active power that the power systems require for spinning reserve and PSRs is the active power spinning reserve of thermal generating unit s.


## 2.2.4. Ramp Rate Constraint

During the operation of thermal units, generation increase and generation decrease are constrained by the ramp up and down the rate limit, respectively. Thus, initial power output of units and the ramp rate limits will determine the lowest and the highest power output of unit at the moment. The description can be presented in the following inequalities [20].
(10)Ps−Ps0≤RURs,ifhighergenerationisrequired,Ps,0−Ps≤RDRs,iflowergenerationisrequired.


## 2.2.5. Real Power Balance

The total real power output of generating units satisfies total load demand plus system power losses [1].
(11)∑s=1NPs=PLD+PTL,and the total power loss is calculated using Kron’s (1):
(12)PLD=∑s=1N∑h=1NPsBshPh+∑s=1NB0sPs+B00.In addition to (11), real power balance is also presented by the following equation since real transmission power networks are considered [17].
(13)Pi−PLDi=Vi∑j=1NbVjGijcosδi−δj+Bijsinδi−δj,i=1,…,Nbus,where Gij and Bij are the real part and the imaginary part of admittance of a transmission line connecting buses i and j.Similarly, reactive power balance is also modeled as follows [17]:
(14)Qi+Qci−QLDi=Vi∑j=1NbVjGijsinδi−δj−Bijcosδi−δj,i=1,…,Nbus,where Qci is the generated reactive power of capacitors at bus i and constrained by the following inequality [17]:
(15)Qci,min≤Qci≤Qci,max,i=1,…,Nc,where Qci,min and Qci,max are, respectively, the lowest and the highest reactive power of the capacitor bank at bus i.


## 2.2.6. Transformer Tap Constraints

In a real transmission power network, transformers are usually used to step down voltage and supply electricity to lower voltage transmission lines. The voltage of the lower voltage lines can be adjusted by setting the tap of the transformers, so that the following rule is always exactly met [17].
(16)Ti,min≤Ti≤Ti,max,i=1,…,Ntap,where Ti,min and Ti,max are the lower limit and the upper limit of transformer’s tap at bus i.


## 2.2.7. Stably Working Power Network Constraints

In order to operate electric power network continuously, load bus voltage and apparent power flow of transmission lines must follow the following limitations [17].
(17)VLi,min≤VLi≤VLi,max,i=1,…,Nload,MVAi≤MVAi,max,i=1,…,Nline,where VLi,min and VLi,max are, respectively, lower voltage and upper voltage of load at bus i; MVAi,max is the apparent power flow capacity of line i; and MVAi is the current apparent power flow of line i.


## 3. Firefly Algorithms

### 3.1. The Conventional Firefly Algorithm

Firefly algorithm, a metaheuristic algorithm, was first built by Yang in 2008 for solving optimization problems [45]. The configuration of FA was developed by the three following ideas:
(i)
Each firefly attracts another one by its brightness.
(ii)
The fireflies with higher brightness have higher attractiveness level to other ones.
(iii)
The fireflies with lower brightness level move to other ones with higher brightness level.The three behaviors of nature fireflies have inspired Yang in building an optimization algorithm called firefly algorithm. There is a mutual connection between the behaviors of fireflies and the construction of FA. In fact, each firefly corresponding to each optimal solution will own its brightness corresponding to the fitness function of the optimal solutions. The action that the fireflies with darker brightness will look for and get to other fireflies producing higher brightness level is similar to the newly produced solutions based on old solutions with a better fitness function. Consequently, in FA, each old solution can be newly produced several times depending on the comparison of its brightness with other ones. As a result, only one new solution of each old solution is kept based on the comparison of fitness function.Suppose that each solutioniXi is a position of firefly i at the current iteration. When the fitness function of solution i is higher than that of another solution j, the distance between the firefly i and j is obtained by using the following model.
(18)rij=Xi−Xj2.Then, the updated distance is employed to be substituted into another (19) to calculate a new attractiveness. Then, the new position for the ith considered firefly can be determined corresponding to the generation of a new solution of the ith solution. The procedure of generating a new solution is carried out as (20).
(19)β=β0e−γrij2,(20)Xijnew=Xi+β.rand.ΔXij+rand,where rand is a random number of solution i and β0 is the attractiveness at zero distance and normally set to 1. Xj is a solution having lower fitness function than Xi; and ΔXij is an updated step size calculated by the following model.
(21)ΔXij=Xj−Xi.Equations (18), (19), and (20) continues to be determined for the ith solution until there are no more solutions that have lower fitness function.In summary, for each solutioni, we can have either one or higher than one new solution or no new solution depending on the fitness comparison between solution i and other solutions among the current population. The statement can be described according to the following term.
(22)Xinew=Xi,ifXiisXGbest,XiGbestnew,ifXjisXGbest,XijnewwithFTbest,otherwise.For the first term in (22), if the considered solution i is also the global best solution, there will be no any new solution generated. For the second case, there will be only one new solution, XiGbestnew is generated if the considered solution is the second best solution, and Xj is the global best solution XGbest among the population. For other cases, it means Xi is the third best solution or worse than the third best solution and even it is the worst solution, there will be from two new solutions to Npop−1 new solutions Xijnew. In this regard, the set of new solutions of solution i will be evaluated by fitness function value comparison and the best one with the lowest fitness FTGbest is retained while others are discarded.The whole description of FA is shown in detail in flowchart of Figure4.Figure 4
The flowchart of implementing FA for a general optimization problem.


### 3.2. The Proposed Improved Firefly Algorithm

FA has been applied for solving different optimization problems such as OOTGU problem and benchmark optimization functions, and its results have been better than those from other popular methods like PSO, GA, and DE. But the method has not met the demand of dealing with large-scale systems and complicated constraints of OOTGU problem, and many studies have worked on improving it [49–51, 53]. In fact, FA still suffers from several drawbacks such as premature convergence to local optimum or convergence to near global optimum with high number of iterations [54]. In order to solve the difficulties of using FA, we suggest three improvements of the attractiveness β (namely, that is the radius), updated step size ∆Xij, and random number (rand) in (20). The detail of the three improvements is as follows.


#### 3.2.1. The First Improvement

Instead of using the distance between the considered solutioni and another better solution to determine the radius, the best solution XGbest is recommended to be used for calculating the radius as shown in the expression.
(23)ribest=Xi−XGbest2.We useXGbest to replace Xj in (18) because the change can produce a more effective scaling factor, which is the attractiveness β. The attractiveness β is a function with respect to the radius as seen in (18). The usefulness of the change is demonstrated in the numerical results where IFA1 (FA with the first improvement) is better than FA in terms of optimal solution reflected via the best solution comparison and the stable solution search ability reflected via the standard deviation comparison.


#### 3.2.2. The Second Proposed Improvement

In the second improvement, we propose a new technique to generate updated step size and obtain lower solution fitness than those of FA. It is clear that the manner of producing an updated step size by using (21) is similar to mutation operation of differential evolution algorithm (DEA) in which β acts as mutation factor ranging from 0 to 2. Some studies before [55, 56] have pointed out the disadvantages of DEA such as low convergence to global optimum or trapping into local optimum easily. Consequently, the proposed improvement aims to tackle the constrictions of FA by using the following expressions.
(24)ΔX1ij=XGbest−Xworst,ΔX2ij=Xj−Xi+Xr1−Xr2.It is clear that the first updated step sizeΔX1ij is less than the second updated step size ΔX2ij. Thus, the use purpose of each one is different, in which the former narrows the search zone nearby old solutions while the latter can expand exploitation to avoid falling into the same solutions. Derived from the difference, a condition of using either ΔX1ij or ΔX2ij is suggested based on the fitness ratios as follows:
(25)ΔXij=ΔX1ij,ifFRi>FRNpop,ΔX2ij,otherwise.FRi and FRNpop are not fixed at constant values, but they are not control parameters and their values can be exactly calculated by
(26)FRi=FTi−FTGbestFTGbest,FRNpop=FTNpop−FTGbestFTGbest.It is clear that there are two possibilities for comparison betweenFRi and FRNpop. FRi is either less or higher than FRNpop. Since the first assumption occurs, it means that solution i is better than the average solution of the whole population and if the second assumption happens, solution i is worse than the mentioned average solution. In the case that solution i is better than the mentioned average solution, solution i is very close to current good solutions, which are better than the mentioned average solution. So if we use the first updated step size ΔX1ij, the change from current solution i is small and the new solution can be very close to or the same as the current good solutions. Clearly, the second updated step size ΔX2ij is much larger than ΔX1ij and it can push the obtained new solution far away from the current good solutions, avoiding jumping to local optimum zone and premature convergence. On the contrary, since solution i is worse than the mentioned average solution, it is not close to the current good solutions and the search around it with a not large step size ΔX1ij is better. The impact of the proposed improvement is significant for obtaining good result of the proposed method via the observation of comparison between FA and IFA2 (FA with the second improvement) in numerical results.


#### 3.2.3. The Third Proposed Improvement

In the third proposed improvement, we suggest using normal distribution (randn) instead of using uniform distribution (rand) like the conventional FA. The suggestion aim is to diversify the search zone as used in Levy flights of Cuckoo search algorithm [46]. As a result, the equation of producing new solutions is as follows:
(27)Xijnew=Xi+β.randn.ΔXij+randn.The effectiveness of the proposed improvement is demonstrated by comparing IFA3 (FA with the third modification) with FA in numerical results.The implementation of the proposed IFA for a general optimization problem is similar to the flowchart shown in Figure4 of FA; however, the difference between the two considered algorithms is the way to produce new solutions Xijnew (shown in Step 5). The whole search procedure of the proposed IFA is described in Algorithm 1.Algorithm 1: The pseudo code of the proposed IFA method.
Set control parameters to values
Randomly initialize populationNpop
Calculate fitness of the whole populationFTii=1,…,Npop
Determine the best solutionXGbest and the worst solution XWorst
SetIter=0;
WhileIter<NIter
SetIter=Iter+1;
CalculateFTNpop, FRi and FRNpop
Fori=1:Npop
Forj=1:Npop
IfFTi>FTj
DetermineΔX1ij and ΔX2ij
ΔX1ij=XGbest−XWorst
ΔX2ij=Xj−Xi+Xr1−Xr2
IfFRi>FRNpop
ΔXij=ΔX1ij
Else
ΔXij=ΔX2ij
End
DetermineriBest
riBest=Xi−XGBest2
Determine attractivenessβ and new solution Xijnew
β=β0e−γrij2
Xijnew=Xi+β.randn.ΔXij+randn
Calculate fitness function ofXijnew
End% if
End% for
Compare each old solution and each new solution to keep better one.
Compare these fitness function to choose the best fitness function
SetXijnew with the lowest fitness to Xinew
End% for
Determine the best solutionXGbest and the worst solution XWorst
End% while


## 3.1. The Conventional Firefly Algorithm

Firefly algorithm, a metaheuristic algorithm, was first built by Yang in 2008 for solving optimization problems [45]. The configuration of FA was developed by the three following ideas:
(i)
Each firefly attracts another one by its brightness.
(ii)
The fireflies with higher brightness have higher attractiveness level to other ones.
(iii)
The fireflies with lower brightness level move to other ones with higher brightness level.The three behaviors of nature fireflies have inspired Yang in building an optimization algorithm called firefly algorithm. There is a mutual connection between the behaviors of fireflies and the construction of FA. In fact, each firefly corresponding to each optimal solution will own its brightness corresponding to the fitness function of the optimal solutions. The action that the fireflies with darker brightness will look for and get to other fireflies producing higher brightness level is similar to the newly produced solutions based on old solutions with a better fitness function. Consequently, in FA, each old solution can be newly produced several times depending on the comparison of its brightness with other ones. As a result, only one new solution of each old solution is kept based on the comparison of fitness function.Suppose that each solutioniXi is a position of firefly i at the current iteration. When the fitness function of solution i is higher than that of another solution j, the distance between the firefly i and j is obtained by using the following model.
(18)rij=Xi−Xj2.Then, the updated distance is employed to be substituted into another (19) to calculate a new attractiveness. Then, the new position for the ith considered firefly can be determined corresponding to the generation of a new solution of the ith solution. The procedure of generating a new solution is carried out as (20).
(19)β=β0e−γrij2,(20)Xijnew=Xi+β.rand.ΔXij+rand,where rand is a random number of solution i and β0 is the attractiveness at zero distance and normally set to 1. Xj is a solution having lower fitness function than Xi; and ΔXij is an updated step size calculated by the following model.
(21)ΔXij=Xj−Xi.Equations (18), (19), and (20) continues to be determined for the ith solution until there are no more solutions that have lower fitness function.In summary, for each solutioni, we can have either one or higher than one new solution or no new solution depending on the fitness comparison between solution i and other solutions among the current population. The statement can be described according to the following term.
(22)Xinew=Xi,ifXiisXGbest,XiGbestnew,ifXjisXGbest,XijnewwithFTbest,otherwise.For the first term in (22), if the considered solution i is also the global best solution, there will be no any new solution generated. For the second case, there will be only one new solution, XiGbestnew is generated if the considered solution is the second best solution, and Xj is the global best solution XGbest among the population. For other cases, it means Xi is the third best solution or worse than the third best solution and even it is the worst solution, there will be from two new solutions to Npop−1 new solutions Xijnew. In this regard, the set of new solutions of solution i will be evaluated by fitness function value comparison and the best one with the lowest fitness FTGbest is retained while others are discarded.The whole description of FA is shown in detail in flowchart of Figure4.Figure 4
The flowchart of implementing FA for a general optimization problem.


## 3.2. The Proposed Improved Firefly Algorithm

FA has been applied for solving different optimization problems such as OOTGU problem and benchmark optimization functions, and its results have been better than those from other popular methods like PSO, GA, and DE. But the method has not met the demand of dealing with large-scale systems and complicated constraints of OOTGU problem, and many studies have worked on improving it [49–51, 53]. In fact, FA still suffers from several drawbacks such as premature convergence to local optimum or convergence to near global optimum with high number of iterations [54]. In order to solve the difficulties of using FA, we suggest three improvements of the attractiveness β (namely, that is the radius), updated step size ∆Xij, and random number (rand) in (20). The detail of the three improvements is as follows.


### 3.2.1. The First Improvement

Instead of using the distance between the considered solutioni and another better solution to determine the radius, the best solution XGbest is recommended to be used for calculating the radius as shown in the expression.
(23)ribest=Xi−XGbest2.We useXGbest to replace Xj in (18) because the change can produce a more effective scaling factor, which is the attractiveness β. The attractiveness β is a function with respect to the radius as seen in (18). The usefulness of the change is demonstrated in the numerical results where IFA1 (FA with the first improvement) is better than FA in terms of optimal solution reflected via the best solution comparison and the stable solution search ability reflected via the standard deviation comparison.


### 3.2.2. The Second Proposed Improvement

In the second improvement, we propose a new technique to generate updated step size and obtain lower solution fitness than those of FA. It is clear that the manner of producing an updated step size by using (21) is similar to mutation operation of differential evolution algorithm (DEA) in which β acts as mutation factor ranging from 0 to 2. Some studies before [55, 56] have pointed out the disadvantages of DEA such as low convergence to global optimum or trapping into local optimum easily. Consequently, the proposed improvement aims to tackle the constrictions of FA by using the following expressions.
(24)ΔX1ij=XGbest−Xworst,ΔX2ij=Xj−Xi+Xr1−Xr2.It is clear that the first updated step sizeΔX1ij is less than the second updated step size ΔX2ij. Thus, the use purpose of each one is different, in which the former narrows the search zone nearby old solutions while the latter can expand exploitation to avoid falling into the same solutions. Derived from the difference, a condition of using either ΔX1ij or ΔX2ij is suggested based on the fitness ratios as follows:
(25)ΔXij=ΔX1ij,ifFRi>FRNpop,ΔX2ij,otherwise.FRi and FRNpop are not fixed at constant values, but they are not control parameters and their values can be exactly calculated by
(26)FRi=FTi−FTGbestFTGbest,FRNpop=FTNpop−FTGbestFTGbest.It is clear that there are two possibilities for comparison betweenFRi and FRNpop. FRi is either less or higher than FRNpop. Since the first assumption occurs, it means that solution i is better than the average solution of the whole population and if the second assumption happens, solution i is worse than the mentioned average solution. In the case that solution i is better than the mentioned average solution, solution i is very close to current good solutions, which are better than the mentioned average solution. So if we use the first updated step size ΔX1ij, the change from current solution i is small and the new solution can be very close to or the same as the current good solutions. Clearly, the second updated step size ΔX2ij is much larger than ΔX1ij and it can push the obtained new solution far away from the current good solutions, avoiding jumping to local optimum zone and premature convergence. On the contrary, since solution i is worse than the mentioned average solution, it is not close to the current good solutions and the search around it with a not large step size ΔX1ij is better. The impact of the proposed improvement is significant for obtaining good result of the proposed method via the observation of comparison between FA and IFA2 (FA with the second improvement) in numerical results.


### 3.2.3. The Third Proposed Improvement

In the third proposed improvement, we suggest using normal distribution (randn) instead of using uniform distribution (rand) like the conventional FA. The suggestion aim is to diversify the search zone as used in Levy flights of Cuckoo search algorithm [46]. As a result, the equation of producing new solutions is as follows:
(27)Xijnew=Xi+β.randn.ΔXij+randn.The effectiveness of the proposed improvement is demonstrated by comparing IFA3 (FA with the third modification) with FA in numerical results.The implementation of the proposed IFA for a general optimization problem is similar to the flowchart shown in Figure4 of FA; however, the difference between the two considered algorithms is the way to produce new solutions Xijnew (shown in Step 5). The whole search procedure of the proposed IFA is described in Algorithm 1.Algorithm 1: The pseudo code of the proposed IFA method.
Set control parameters to values
Randomly initialize populationNpop
Calculate fitness of the whole populationFTii=1,…,Npop
Determine the best solutionXGbest and the worst solution XWorst
SetIter=0;
WhileIter<NIter
SetIter=Iter+1;
CalculateFTNpop, FRi and FRNpop
Fori=1:Npop
Forj=1:Npop
IfFTi>FTj
DetermineΔX1ij and ΔX2ij
ΔX1ij=XGbest−XWorst
ΔX2ij=Xj−Xi+Xr1−Xr2
IfFRi>FRNpop
ΔXij=ΔX1ij
Else
ΔXij=ΔX2ij
End
DetermineriBest
riBest=Xi−XGBest2
Determine attractivenessβ and new solution Xijnew
β=β0e−γrij2
Xijnew=Xi+β.randn.ΔXij+randn
Calculate fitness function ofXijnew
End% if
End% for
Compare each old solution and each new solution to keep better one.
Compare these fitness function to choose the best fitness function
SetXijnew with the lowest fitness to Xinew
End% for
Determine the best solutionXGbest and the worst solution XWorst
End% while


## 3.2.1. The First Improvement

Instead of using the distance between the considered solutioni and another better solution to determine the radius, the best solution XGbest is recommended to be used for calculating the radius as shown in the expression.
(23)ribest=Xi−XGbest2.We useXGbest to replace Xj in (18) because the change can produce a more effective scaling factor, which is the attractiveness β. The attractiveness β is a function with respect to the radius as seen in (18). The usefulness of the change is demonstrated in the numerical results where IFA1 (FA with the first improvement) is better than FA in terms of optimal solution reflected via the best solution comparison and the stable solution search ability reflected via the standard deviation comparison.


## 3.2.2. The Second Proposed Improvement

In the second improvement, we propose a new technique to generate updated step size and obtain lower solution fitness than those of FA. It is clear that the manner of producing an updated step size by using (21) is similar to mutation operation of differential evolution algorithm (DEA) in which β acts as mutation factor ranging from 0 to 2. Some studies before [55, 56] have pointed out the disadvantages of DEA such as low convergence to global optimum or trapping into local optimum easily. Consequently, the proposed improvement aims to tackle the constrictions of FA by using the following expressions.
(24)ΔX1ij=XGbest−Xworst,ΔX2ij=Xj−Xi+Xr1−Xr2.It is clear that the first updated step sizeΔX1ij is less than the second updated step size ΔX2ij. Thus, the use purpose of each one is different, in which the former narrows the search zone nearby old solutions while the latter can expand exploitation to avoid falling into the same solutions. Derived from the difference, a condition of using either ΔX1ij or ΔX2ij is suggested based on the fitness ratios as follows:
(25)ΔXij=ΔX1ij,ifFRi>FRNpop,ΔX2ij,otherwise.FRi and FRNpop are not fixed at constant values, but they are not control parameters and their values can be exactly calculated by
(26)FRi=FTi−FTGbestFTGbest,FRNpop=FTNpop−FTGbestFTGbest.It is clear that there are two possibilities for comparison betweenFRi and FRNpop. FRi is either less or higher than FRNpop. Since the first assumption occurs, it means that solution i is better than the average solution of the whole population and if the second assumption happens, solution i is worse than the mentioned average solution. In the case that solution i is better than the mentioned average solution, solution i is very close to current good solutions, which are better than the mentioned average solution. So if we use the first updated step size ΔX1ij, the change from current solution i is small and the new solution can be very close to or the same as the current good solutions. Clearly, the second updated step size ΔX2ij is much larger than ΔX1ij and it can push the obtained new solution far away from the current good solutions, avoiding jumping to local optimum zone and premature convergence. On the contrary, since solution i is worse than the mentioned average solution, it is not close to the current good solutions and the search around it with a not large step size ΔX1ij is better. The impact of the proposed improvement is significant for obtaining good result of the proposed method via the observation of comparison between FA and IFA2 (FA with the second improvement) in numerical results.


## 3.2.3. The Third Proposed Improvement

In the third proposed improvement, we suggest using normal distribution (randn) instead of using uniform distribution (rand) like the conventional FA. The suggestion aim is to diversify the search zone as used in Levy flights of Cuckoo search algorithm [46]. As a result, the equation of producing new solutions is as follows:
(27)Xijnew=Xi+β.randn.ΔXij+randn.The effectiveness of the proposed improvement is demonstrated by comparing IFA3 (FA with the third modification) with FA in numerical results.The implementation of the proposed IFA for a general optimization problem is similar to the flowchart shown in Figure4 of FA; however, the difference between the two considered algorithms is the way to produce new solutions Xijnew (shown in Step 5). The whole search procedure of the proposed IFA is described in Algorithm 1.Algorithm 1: The pseudo code of the proposed IFA method.
Set control parameters to values
Randomly initialize populationNpop
Calculate fitness of the whole populationFTii=1,…,Npop
Determine the best solutionXGbest and the worst solution XWorst
SetIter=0;
WhileIter<NIter
SetIter=Iter+1;
CalculateFTNpop, FRi and FRNpop
Fori=1:Npop
Forj=1:Npop
IfFTi>FTj
DetermineΔX1ij and ΔX2ij
ΔX1ij=XGbest−XWorst
ΔX2ij=Xj−Xi+Xr1−Xr2
IfFRi>FRNpop
ΔXij=ΔX1ij
Else
ΔXij=ΔX2ij
End
DetermineriBest
riBest=Xi−XGBest2
Determine attractivenessβ and new solution Xijnew
β=β0e−γrij2
Xijnew=Xi+β.randn.ΔXij+randn
Calculate fitness function ofXijnew
End% if
End% for
Compare each old solution and each new solution to keep better one.
Compare these fitness function to choose the best fitness function
SetXijnew with the lowest fitness to Xinew
End% for
Determine the best solutionXGbest and the worst solution XWorst
End% while


## 4. The Implementation of IFA for OOTGU Problem

### 4.1. Dealing with Load Demand-Supply Balance and Other Constraints

In order to deal all constraints, decision variables and dependent variables are chosen thoroughly. For the case that constraints of real transmission power network are not taken into account, power output of the first (N−1) thermal units is included in the position of each firefly as shown in
(28)Xi=P1,i,P2,i,…,PN−1,i,i=1,…,Npop,where
(29)Xmin≤Xi≤Xmax,(30)Xmin=P1,min,P2,min,…,PN−1,min,(31)Xmax=P1,max,P2,max,…,PN−1,max.As a result, the load demand-supply balance constraint can be dealt successfully by using the dependent variablePN,i obtained by the below model [15].
(32)PN,i=PLD+PTL−∑s=1N−1Ps,i.For complicated case with the consideration of all transmission power network constraints, active power output of all generators excluding at slack bus, voltage of all generators, tap setting of all transformers, and reactive power output of all capacitors are chosen as decision variables. Thus,Xi, Xmin, and Xmax are not the same as (28), (30), and (31) and can be obtained as follows:
(33)Xi=P2,i,…,PN,i,V1,i,…,VN,i,T1,i,…,TNt,i,Qc1,i,…,QcNc,i,Xmin=P2,min,…,PN,min,V1,min,…,VN,min,T1,min,…,TNt,min,Qc1,min,…,QcNc,min,Xmax=P2,max,…,PN,max,V1,max,…,VN,max,T1,max,…,TNt,max,Qc1,max,…,QcNc,max.Note thatP1,i is not included in the set of decision variable because it is supposed to be the active power of generator located at slack bus. Other dependent variables are obtained by running optimal power flow program, Mathpower.


### 4.2. Penalizing the Violation of Dependent Variables

There is a possibility thatPN,i can violate its limitation, lower than the lowest generation or higher than the highest generation. Therefore, the violation must be controlled and considered in the quality evaluation of solutions. The task is done by calculating the penalty term as indicated in the following equation.
(34)Penalty1i=PN,i−PN,max,ifPN,i>PN,max,PN,min−PN,i,ifPN,i<PN,min,0,ifPN,min≤PN,i≤PN,max,(35)Penalty2i=max0,PSR−∑s=1NPSRs2.Similarly, when the real transmission power network constraints are running, the third penalty term calculation is carried out by the following model.(36)Penalty3i=P1,i−P1,ilim2+∑s=1NQs,i−Qs,ilim2+∑j=1NloadVLj,i−VLj,ilim2+∑j=1NlineMVAj,i−MVAj,ilim2,where the limits related to dependent variables are determined by
(37)P1,ilim=P1,max,ifP1,i>P1,max,P1,min,ifP1,i<P1,min,P1,i,else,Qs,ilim=Qs,max,ifQs,i>Qs,max,Qs,min,ifQs,i<Qs,min,Qs,i,else,VLj,ilim=VLj,max,ifVLj,i>VLj,max,VLj,min,ifVLj,i<VLj,min,VLj,i,else,MVAj,ilim=MVAjmax,ifMVAj,i>MVAjmax,MVAj,i,otherwise.


### 4.3. Fitness Function

The fitness function of all solution should be determined to arrange the effectiveness of all the solutions. Consequently, the fitness function considering objective function and the sum of penalty terms is shown in the following model.(38)FTi=∑s=1NFsPs+PF×penalty1i+penalty2i+penalty3i2,where PF is penalty factor used to amplify the violation of the dependent variable [16].


### 4.4. The Procedure of IFA for OOTGU Problem

The whole search process of the proposed method when solving the considered problem is described in detail as follows.Step 1.
Initialize a population withNpop solutions Xi,i=1,…,Npop satisfying the condition of (29).Step 2.
Calculate dependent variables for solutioni by using (32) or running Mathpower optimal power flow program.Step 3.
Determine penalty terms by using (34), (35), and (36).Step 4.
Find the fuel cost coefficients corresponding to power output of TGU in case of multifuel options.Step 5.
Calculate fitness functionFTi by using (38).Step 6.
Assign the solution with lowest fitness into the best oneXGbest and start the first iteration (Iter=1).Step 7.
Consider solutioni=1.Step 8.
Produce new solutionsXijnew for solution i by using three proposed modifications.Step 9.
Check limits forXijnew and perform Steps (2)–(4) above.Step 10.
Calculate fitness function for all new solutionsXijnew and assign the solution with the lowest fitness into Xinew.Step 11.
Ifi=Npop, go to the next step. Otherwise, set i=i+1 and back to Step 8.Step 12.
Compare each new solutionXinew to each old solution Xi to keep better one.Step 13.
Determine the best solution.Step 14.
If Iter is equal to the highest iteration numberNIter, stop the iterative algorithm. Otherwise, set Iter=Iter+1 and back to Step 7.


## 4.1. Dealing with Load Demand-Supply Balance and Other Constraints

In order to deal all constraints, decision variables and dependent variables are chosen thoroughly. For the case that constraints of real transmission power network are not taken into account, power output of the first (N−1) thermal units is included in the position of each firefly as shown in
(28)Xi=P1,i,P2,i,…,PN−1,i,i=1,…,Npop,where
(29)Xmin≤Xi≤Xmax,(30)Xmin=P1,min,P2,min,…,PN−1,min,(31)Xmax=P1,max,P2,max,…,PN−1,max.As a result, the load demand-supply balance constraint can be dealt successfully by using the dependent variablePN,i obtained by the below model [15].
(32)PN,i=PLD+PTL−∑s=1N−1Ps,i.For complicated case with the consideration of all transmission power network constraints, active power output of all generators excluding at slack bus, voltage of all generators, tap setting of all transformers, and reactive power output of all capacitors are chosen as decision variables. Thus,Xi, Xmin, and Xmax are not the same as (28), (30), and (31) and can be obtained as follows:
(33)Xi=P2,i,…,PN,i,V1,i,…,VN,i,T1,i,…,TNt,i,Qc1,i,…,QcNc,i,Xmin=P2,min,…,PN,min,V1,min,…,VN,min,T1,min,…,TNt,min,Qc1,min,…,QcNc,min,Xmax=P2,max,…,PN,max,V1,max,…,VN,max,T1,max,…,TNt,max,Qc1,max,…,QcNc,max.Note thatP1,i is not included in the set of decision variable because it is supposed to be the active power of generator located at slack bus. Other dependent variables are obtained by running optimal power flow program, Mathpower.


## 4.2. Penalizing the Violation of Dependent Variables

There is a possibility thatPN,i can violate its limitation, lower than the lowest generation or higher than the highest generation. Therefore, the violation must be controlled and considered in the quality evaluation of solutions. The task is done by calculating the penalty term as indicated in the following equation.
(34)Penalty1i=PN,i−PN,max,ifPN,i>PN,max,PN,min−PN,i,ifPN,i<PN,min,0,ifPN,min≤PN,i≤PN,max,(35)Penalty2i=max0,PSR−∑s=1NPSRs2.Similarly, when the real transmission power network constraints are running, the third penalty term calculation is carried out by the following model.(36)Penalty3i=P1,i−P1,ilim2+∑s=1NQs,i−Qs,ilim2+∑j=1NloadVLj,i−VLj,ilim2+∑j=1NlineMVAj,i−MVAj,ilim2,where the limits related to dependent variables are determined by
(37)P1,ilim=P1,max,ifP1,i>P1,max,P1,min,ifP1,i<P1,min,P1,i,else,Qs,ilim=Qs,max,ifQs,i>Qs,max,Qs,min,ifQs,i<Qs,min,Qs,i,else,VLj,ilim=VLj,max,ifVLj,i>VLj,max,VLj,min,ifVLj,i<VLj,min,VLj,i,else,MVAj,ilim=MVAjmax,ifMVAj,i>MVAjmax,MVAj,i,otherwise.


## 4.3. Fitness Function

The fitness function of all solution should be determined to arrange the effectiveness of all the solutions. Consequently, the fitness function considering objective function and the sum of penalty terms is shown in the following model.(38)FTi=∑s=1NFsPs+PF×penalty1i+penalty2i+penalty3i2,where PF is penalty factor used to amplify the violation of the dependent variable [16].


## 4.4. The Procedure of IFA for OOTGU Problem

The whole search process of the proposed method when solving the considered problem is described in detail as follows.Step 1.
Initialize a population withNpop solutions Xi,i=1,…,Npop satisfying the condition of (29).Step 2.
Calculate dependent variables for solutioni by using (32) or running Mathpower optimal power flow program.Step 3.
Determine penalty terms by using (34), (35), and (36).Step 4.
Find the fuel cost coefficients corresponding to power output of TGU in case of multifuel options.Step 5.
Calculate fitness functionFTi by using (38).Step 6.
Assign the solution with lowest fitness into the best oneXGbest and start the first iteration (Iter=1).Step 7.
Consider solutioni=1.Step 8.
Produce new solutionsXijnew for solution i by using three proposed modifications.Step 9.
Check limits forXijnew and perform Steps (2)–(4) above.Step 10.
Calculate fitness function for all new solutionsXijnew and assign the solution with the lowest fitness into Xinew.Step 11.
Ifi=Npop, go to the next step. Otherwise, set i=i+1 and back to Step 8.Step 12.
Compare each new solutionXinew to each old solution Xi to keep better one.Step 13.
Determine the best solution.Step 14.
If Iter is equal to the highest iteration numberNIter, stop the iterative algorithm. Otherwise, set Iter=Iter+1 and back to Step 7.


## 5. Numerical Results

In the paper, the efficiency of proposed IFA method is verified by testing on six main cases with different characteristics of fuel cost function, thermal unit constraints, and power network constraints. The detail of test systems is described as follows.Case 1.
Six thermal generating unit test system with POZ constraint and power loss constraints [57].Case 2.
Fifteen thermal generating unit test system with POZ constraint and without power loss constraints [22].Case 3.
Twenty thermal generating unit test system with power loss constraints [3].Case 4.
Forty thermal generating unit test system with valve-point loading effects [58].Case 5.
Ten thermal generating units with multifuels.Subcase 5.1.
Neglect valve-point loading effects and power losses [7].
Subcase 5.2.
Consider valve-point loading effects and neglect power losses [20].
Subcase 5.3.
Consider valve-point loading effects and the constraints of POZ, ramp rate limits, and active power spinning reserve [20].Case 6.
The IEEE-118 bus system with all constraints of real transmission power networks [59].The population sizeNpop and the highest iteration number NIter selected for implementation of IFA are reported in Table 1. For each study case, each method is run one hundred independent trials by using MATLAB and computer with 4 GB of Ram and 2.4 GHz processor.Table 1
The selection of population size and the highest iteration number.
Case
Number of units
Npop
NIter
1
6
10
30
2
15
10
100
3
20
10
300
4
40
10
10,000
5
10
10
100
10
10
100
10
10
300
6
54
20
500


### 5.1. Impact Analysis of Each Proposed Improvement

In this section, the contribution of each proposed modification to the performance of the proposed method is investigated by running conventional FA, FA with the first proposed improvement (IFA1), FA with the second proposed improvement (IFA2), FA with the third proposed improvement (IFA3), and the proposed IFA method on five benchmark functions and two power systems. Such benchmark function is given in Table2 [17, 60] while the two power systems with 6 units and 15 units considering POZ constraints are taken from Cases 1 and 2.Table 2
Definition of five applied benchmark functions.
Function
Definition
Domain
f1(Sphere)
f1=∑i=1Dxi2
xi∈−5.12,5.12D=30
f2(Rastrigin)
f2=∑i=1D10+xi2−10cos2πxi
xi∈−5.12,5.12D=30
f3(Ackley)
f3=20+e−20exp−0.21D∑i=1Dxi2−exp1D∑i=1Dcos2πxi
xi∈−30,30D=30
f4(Griewangk)
f4=1+∑i=1Dxi24000−∏i=1Dcosxii
xi∈−600,600D=30
f5(Schwelfel)
f5=∑i=1Dxi+∏iDxi
xi∈−10,10D=30


#### 5.1.1. Testing on Benchmark Functions

In order to implement FA, IFA1, IFA2, IFA3, and the proposed method, the population size and number of iterations are set to 50 and 1000, respectively. In addition, we also employ Wilcoxon signed-rank test [61] and set the significance improvement level to 0.01 for calculating Z value, W value, and p value for better comparison between the proposed method and others. As a result, the obtained results in terms of minimum fitness, standard deviation of 50 trial runs, Z value, W value, and p value are reported in Table 3. Among the result numbers, the comparison of the best cost is to reflect the best optimal solution and comparison of the standard deviation cost is to reflect the stabilization of search ability. As pointed in [62], p value can be either equal to or smaller than 0.01 or higher than 0.01. In case that p value of a compared method is higher than 0.01, the significance improvement level of the proposed method over the compared method is not adopted but the significance level is accepted in case such p value is less than 0.01 or equal to 0.01 [62]. Table 3 shows that all minimum values and standard deviation of the proposed method are smaller those of other ones meanwhile Z value, W value, and p value are, respectively, equal to −6.154, 0, and 0. p value = 0 indicates that the proposed method reaches the significance improvement level over other FA variants while W value = 0 determines the accuracy of Z value and p value. For other comparisons among FA, IFA1, IFA2, and IFA3, the minimum fitness and standard deviation reports can indicate that IFA2 is the best modification; IFA3 is the second best modification while IFA1 is only better than conventional FA. The further investigation of performance of each FA methods can be seen clearly in the following section since two power systems with 6 units and 15 units considering POZ constraints are employed for testing.Table 3
Result comparison for five benchmark functions.
Function
Method
Min.
Std. dev.
Z value
W value
P value
f1(Sphere)
FA
47,343
5474.1
−6.154
0
0
IFA1
2381.5
488.8
−6.154
0
0
IFA2
21.77
8.72
−6.154
0
0
IFA3
154.198
24.3
−6.154
0
0
IFA
0.6537
0.1865
f2(Rastrigin)
FA
98.6076
21.2183
−6.154
0
0
IFA1
60.57
125.05
−6.154
0
0
IFA2
6.42
58.08
−6.154
0
0
IFA3
48.58
32.72
−6.154
0
0
IFA
0.17654
11.1474
f3(Ackly)
FA
19.4832
0.2661
−6.154
0
0
IFA1
11.521
1.1903
−6.154
0
0
IFA2
2.7034
0.1023
−6.154
0
0
IFA3
5.577
2.1884
−6.154
0
0
IFA
0.69
0.0481
f4(Griewank)
FA
436.1885
55.157
−6.154
0
0
IFA1
16.67
42.97
−6.154
0
0
IFA2
0.009611
0.165
−6.154
0
0
IFA3
0.33162
5.993
−6.154
0
0
IFA
0.00020344
5.1E − 05
f5(Schwelfel)
FA
49,347
7517.8
−6.154
0
0
IFA1
249.4881
76.65
−6.154
0
0
IFA2
16.19
15.59
−6.154
0
0
IFA3
36.72
36.12
−6.154
0
0
IFA
4.0252
0.292


#### 5.1.2. Testing on Power Systems with POZ Constraints

In this section, the performance of FA, IFA1, IFA2, IFA3, and the proposed IFA method is illustrated by running on power systems of Study Cases1 and 2. Furthermore, four other popular methods consisting of particle swarm optimization (PSO) [63], differential evolution (DE) [64], flower pollination algorithm (FPA) [65], and crow search algorithm (CSA) [66] are also implemented for solving two test systems. For running the methods, we set the population and the number of iteration to those values of the proposed method shown in Table 1. In addition, we also tune other control parameters of these methods to their predetermined ranges. For the implementation of PSO, acceleration factors c1 and c2 are set to the range of [0.2, 2.0] with a step size of 0.2. For the implementation of DE, mutation factor and crossover factor are, respectively, set to range of [0.2, 1.2] and [0.2, 1.0] with a step size of 0.2. For the implementation of FPA, probability of using either local search or global search is set to the range of [0.1, 1.0] with a step size of 0.1. For the implementation of CSA, awareness probability and flight length are, respectively, set to the range of [0.2, 1.0] and [0.2, 2.0] with a step size of 0.2. The obtained results from the proposed method and these compared methods such as minimum cost and standard deviation and the results from Wilcoxon signed-rank test such as Z value, W value, and p value are reported in Tables 4 and 5 for Study Cases 1 and 2.Table 4
Result comparisons for Case1.
Method
Best cost ($)
Std. dev. ($)
Z value
W value
p value
FA
15445.215
26.815
−6.154
0
0
IFA1
15444.918
18.01
−6.154
0
0
IFA2
15443.076
1.371
−6.154
0
0
IFA3
15443.257
17.163
−6.154
0
0
PSO
15444.447
16.5186
−6.154
0
0
FPA
15443.55
2.9594
−6.154
0
0
DE
15443.5
26.8126
−6.154
0
0
CSA
15444.75
6.5271
−6.154
0
0
IFA
15443.075
0.050Table 5
Result comparisons for Case2.
Method
Best cost ($)
Std. dev. ($)
Z value
W value
p value
FA
32885.8432
104.534
−6.154
0
0
IFA1
32781.003
88.441
−6.154
0
0
IFA2
32544.9906
8.3011
−6.154
0
0
IFA3
32712.6581
91.8623
−6.154
0
0
PSO
32571.62
94.63377
−6.154
0
0
FPA
32545.1
12.26141
−6.1057
5
0
DE
32635.04
114.94
−6.154
0
0
CSA
32553.03
61.97552
−6.154
0
0
IFA
32544.9704
0.1438For Case1, the proposed method has obtained the best results with less best cost than FA, IFA1, IFA2, and IFA3, respectively, by $2.14, $1.843, $0.001, and $0.182. Similarly, the standard deviation cost of IFA is, respectively, less than that of FA, IFA1, IFA2, and IFA3 by $26.765, $18.05, $1.322, and $17.113. Furthermore, Z value, W value, and p value also indicate the significance improvement of the proposed method over FA variants. Comparison among IFA1, IFA2, and IFA3 can see that such three improved versions has better performance than FA; meanwhile, IFA2 is the best improvement and IFA1 is the worst improvement. Analysis of Case 2 also gives the same evaluation for FA methods.iFor comparison with PSO, FPA, DE, and CSA, all the costs yielded by the proposed method are smaller than those from PSO, FPA, DE, and CSA. In particular, the standard deviation of the proposed method is approximately zero while this value of other methods is much higher than zero. Furthermore,Z value, W value, and p value calculated by setting the significance improvement level to 0.01 also indicate that the proposed method can improve working performance significantly when compared to PSO, FPA, DE, and CSA. Among IFA1, IFA2, and IFA3, only IFA2 outperforms PSO, FPA, DE, and CSA since the best cost and standard deviation cost of IFA2 are smaller than those of PSO, FPA, DE, and CSA. Consequently, it can conclude that the proposed method has high performance for searching optimal solutions while the second modification has the most impact on results.Besides, the more accurate comparison of the best cost improvement and the standard deviation improvement can be quantitative by using the following equation.(39)IFAresultimprovement%=resultofcomparedmethod‐resultofIFAresultofcomparedmethod×100%.In (39), result of IFA or result of compared method can be the best cost, mean cost, worst cost, or standard deviation cost dependent on decision of evaluation but the two most important evaluations must be used including the best cost improvement and the standard deviation improvement. As result is the best cost, IFA result improvement is the improvement of optimal solution quality of IFA over the compared method. In addition, as the result is standard deviation, the IFA result improvement is the improvement of stabilization of IFA over the compared method. As a result, the improvement of IFA over others in terms of the best optimal solution quality and the stabilization is reported in Table 6 for Cases 1 and 2. As observed from Table 6, the best cost improvement of IFA over FA and other methods is not high but the improvement of standard deviation cost is approximately 100% for the two cases.Table 6
Improvement of IFA over FA for Cases1 and 2.
Case
Compared method
Best cost improvement (%)
Std. dev. improvement (%)
1
FA
0.014
99.81
PSO
0.009
99.70
FPA
0.003
98.31
DE
0.003
99.81
CSA
0.011
99.23
2
FA
1.037
99.86
PSO
0.082
99.85
FPA
0.000
98.83
DE
0.276
99.87
CSA
0.025
99.77The best cost of fifty independent trials for Cases1 and 2 obtained by the proposed IFA and FA is, respectively, plotted in Figures 5 and 6. The figures can evaluate the standout of IFA over FA clearly once all minimum costs of IFA are nearly on a flat line and much lower than those of FA, especially for Case 2 while fifty cost values of FA have high fluctuations for the two cases.Figure 5
Fitness-iteration characteristics obtained FA and IFA for Study Case1.Figure 6
Fitness-iteration characteristics obtained by FA and IFA for Study Case2.


### 5.2. Comparison and Discussion

In order to investigate the further performance of IFA, we carry out the comparisons of results obtained by IFA and many optimization algorithms. In addition to the comparison of the best costs, another comparison criterion is also considered to be the number of fitness evaluationsNFES, which is shown in the following equation.
(40)NFES=ω×Npop×NIter,where ω is the number of generations in each iteration. For some optimization algorithms possessing two new solution generations such as CSA [15] and ICSA [16], ω is 2 while for others possessing one new solution generation like PSO [11] and DE [11], ω is only 1. For the proposed IFA, there is only one new solution generation in each iteration; thus, ω is also equal to 1. As a result, the value of NFES is added into each table for the comparison and the indication is that a method with lower NFES is considered to be more efficient if its best cost is also lower or equal.


#### 5.2.1. Comparisons of Test Case1

Table7 reports comparisons of the best cost, mean cost, worst cost, and NFES from IFA and other methods for Case 1, and especially the best cost improvement of IFA compared to each one is also reported at the final column. It notes that if the best cost improvement (BCI) in % has plus sign (+), the result from IFA is better than the compared method, and if BCI has minus sign (−), the compared one has better optimal solution than IFA. Besides, if BCI equals zero, the proposed method and compared methods have the same best cost and the same quality of optimal solution. The observation of BCI can lead to an evaluation that the proposed method is slightly better than most methods (the improvement is from 0.0002% to +0.1030%). KHA4 [23], OBKHA [24], and EMA [30] have the same cost as the proposed method, but α-MFA [53] has slightly lower cost than IFA. However, the recalculated cost of α-MFA [53] is $15445.46, which is higher than the cost of the proposed method. The comparison of NFES can indicate that the proposed method is at least 4 times faster than other ones since NFES of IFA is 300 but that of other ones is from 1200 (α-MFA [53]) to 50,000 (EMA [30]). In spite of using lower value of NFES, the mean cost and the highest cost of the proposed method are also much less or approximately equal to those from other ones. As a result, it can result in a conclusion that the proposed method is much effective than other methods when solving test Case 1 with quadratic fuel cost function, POZ constraint, and power loss constraint.Table 7
Result comparisons for Case1.
Method
Best cost ($)
Mean cost ($)
Worst cost ($)
NFES
BCI (%)
DE [11]
15,449.77
15,449.87
15,449.78
3600
+0.0433
GA [11]
15,459
15,524
15,469
3600
+0.1030
PSO [11]
15,450
15,492
15,454
3600
+0.0448
KHA1 [23]
15,450.75
15,452.82
15,455.50
10,000
+0.0497
KHA2 [23]
15,448.21
15,450.83
15,453.40
10,000
+0.0333
KHA3 [23]
15,445.36
15,447.22
15,449.60
10,000
+0.0148
KHA4 [23]
15,443.08
15,443.19
15,443.30
10,000
0.0000
OBKHA [24]
15,443.08
15,443.92
15,443.33
—
0.0000
BBO [26]
15,443.10
—
—
10,000
+0.0002
EMA [30]
15,443.07
15,443.08
—
50,000
0.0000
AIS [34]
15,448
—
—
4000
+0.0319
NAPSO [18]
15,443.77
15,443.77
15,443.77
3000
+0.0045
NAPSO1 [18]
15,449.14
15,451.04
15,452.92
3000
+0.0393
NAPSO 2 [18]
15,444.77
15,446.77
15,449.77
3000
+0.0110
NAPSO 3 [18]
15,449.08
15,450.33
15,452.78
3000
+0.0389
α-MFA [53]
15442.9
—
—
1200
−0.0011
IFA
15,443.075
15,443.12
15,443.52
300


#### 5.2.2. Comparisons of Test Case2

This section compares the performance of the proposed method with other ones when dealing with a larger scale system compared to Case1. The obtained results from the proposed method and other compared methods are given in Table 8. It is obviously seen that ICSA [16] has the best optimal solution; however, the confirmation of a valid solution is not done for ICSA because there was no optimal solution reported for the case in [16]. In comparison with mean cost and maximum cost of other ones, the proposed method has higher values because we run the proposed method with NFES of 1000 while others were run by a very high value. In case that we increase NFES of the proposed method to 6000 and higher than 6000, the result is that the mean cost and the highest cost are approximately equal to the best cost. The illustration aims to demonstrate the fast convergence of the proposed method compared to others, and it is clear that the proposed method is better than this compared method for Case 2 with 15 units considering quadratic fuel cost function and POZ constraint.Table 8
Result comparisons for Case2.
Method
Best cost ($)
Mean cost ($)
Worst cost ($)
NFES
CSA [15]
32,544.9704
32,545.006
32,546.6734
8000
ICSA [16]
32,542.5593
32,543.16884
32,546.6574
6000
IGA [22]
—
32,04.736
—
250,000
IGAMU [22]
—
32,544.990
—
50,000
IFA
32544.9704
32,545.2241
32545.54
1000


#### 5.2.3. Comparisons of Test Case3

In the section, a larger system with 20 units considering quadratic fuel cost function and power loss constraint is employed for comparison. The information of the best cost andNFES is shown in Table 9 for comparison with deterministic methods and metaheuristic methods. The best can see that the proposed method can obtain less cost than FA, Lamda, and HM methods in [3] and BBOA in [26] but slightly higher than cost of CSA and ICSA in [16]. The deviation implies that the proposed method is more effective than lambda, HM, and BBOA, but it is not less effective than methods in [16]. Besides, the proposed method has used much less NFES than other ones including FA. As we use a higher number of NFES equal to 10,000, the best cost obtained by the proposed method is also as good as that of methods in [16]. The numerical comparison can show a good performance of the proposed method when dealing with a system with 20 units considering quadratic fuel cost function and power loss constraint.Table 9
Result comparisons for Case3.
Method
Cost ($)
NFES
Lambda [3]
624656.639
—
HM [3]
62456.6341
—
CSA [16]
62456.633
10,000
ICSA [16]
62456.633
10,000
BBOA [26]
62456.7926
20,000
FA
62458.881
10,000
IFA
62456.634
3000


#### 5.2.4. Comparisons of Test Case4

Test Case4 is more complicated than above cases due to the consideration of 40 units and nonconvex fuel cost functions. The result comparisons of the case are shown in Table 10. The values of BCI can see that the proposed method has better cost than methods with plus sign “+” of BCI such as BBO, OLCSO, KHA1, KHA2, KHA3, FA [47], and CFA [51]. In comparison with remaining methods, the proposed method has lower quality of solution with higher cost excluding FPA, which did not report optimal solution and had the surprised lowest value. The comparison of NFES also goes against the proposed method in the evaluation of performance. In fact, the proposed method uses a higher number of NFES compared to CCDE, FA, and KHA methods but it uses smaller NFES than the rest of methods including CFA, which is an improved version of FA. The analysis of best cost and NFES can point out that the application of the proposed method is not as effective as the applications with quadratic fuel cost function in Cases 1, 2, and 3 above.Table 10
Result comparisons for Case4.
Method
Best cost ($)
Mean cost ($)
Worst cost ($)
NFES
BCI (%)
CCDE [13]
121,412.58
121,413.56
121,416.53
70,000
−0.002
CSA [16]
121,412.54
121,520.41
121810.25
160,000
−0.002
BBO [26]
121,426.95
121,508.03
121,688.66
200,000
+0.010
FPA [39]
121,074.50
121,095.70
121,196.30
2000
−0.281
OLCSO [44]
121,415.82
121,504.05
121,460.78
400,000
+0.001
FA [47]
121,415.05
121,416.57
121,424.56
50,000
+0.000
KHA1 [23]
121,460.42
121,468.98
121,477.46
10,000
+0.038
KHA2 [23]
121,448.36
121,453.68
121,461.39
10,000
+0.028
KHA3 [23]
121,423.46
121,428.23
121,433.56
10,000
+0.007
KHA4 [23]
121,412.60
121,413.15
121,415.00
10,000
−0.002
CFA [51]
121,415.02
121,416.23
121,420.67
336,000
+0.00035
IFA
121,414.6
121,549.038
121,787.5
100,000
0.000


#### 5.2.5. Comparisons of Test Case5

This study case considers different fuel cost functions with nonsmooth forms such as the sum of several quadratic functions for Subcase 5.1, the sum of several nonconvex functions for Subcases 5.2 and 5.3, and different complicated constraints such as POZ, ramp rate limits, active power spinning reserve, and power losses for Subcase 5.3. Subcase 5.1 is the simplest among the three subcases considering four load values, 2400, 2500, 2600, and 2700 MW. Table11 reports the comparison of the best cost and the value of NFES from methods. It can see that the proposed method can yield optimal solutions as good as others excluding FA with the worst results and the proposed method’s NFES is also the lowest value excluding comparison to ALHN and EALHN, which are not metaheuristic methods and face to the difficulties of applying for nonconvex fuel cost function. The comparison can imply that the proposed method is very promising for dealing with ten-unit system considering multifuel options where nonconvex fuel cost function or valve-point loading effects are not taken into account.Table 11
The comparison of the best cost (in $/h) for Subcase 5.1.
Method
PD (MW)
NFES
2400
2500
2600
2700
ALHN [8]
481.723
526.239
574.381
623.809
—
EALHN [9]
481.723
526.239
574.381
623.809
—
DE [11]
481.723
526.239
574.381
623.809
12,000
HDP-ICDE [12]
481.7226
526.2388
574.3808
623.809
8000
HRCGA [21]
481.7226
526.2388
574.3808
623.8092
30,000
AIS [34]
481.723
526.24
574.381
623.8092
4000
CSA [15]
—
—
—
623.8092
6000
CFA [51]
—
—
—
623.8339
156,000
FA
505.2337
580.4417
639.287
679.9525
2000
IFA
481.7226
526.2388
574.3808
623.809
1000Comparisons for Subcases 5.2 reported in Table12 can see that KHA methods in [23] are the best ones with the surprised lowest BCI compared to the proposed method; however, the validation confirmation of solutions cannot be done due to the lack of optimal solution report in [23]. Compared to popular methods such as PSO, GA, TSA (Tabu search algorithm), and FA, the proposed method is better in terms of lower cost and lower NFES. On the contrary, the proposed method is less effective than remaining methods in [13, 15, 16, 19]. As we set NFES to higher values such as 3000, 7000, and 10,000, the results obtained are also worse than these methods. Thus, we conclude that the proposed method is not much effective for the subcase considering valve-point loading effects.Table 12
The comparison of costs (in $/h) for Subcase 5.2.
Method
Best cost ($)
Mean cost ($)
Worst cost ($)
NFES
BCI (%)
ICSA [16]
623.8684
623.9495
626.3666
12,000
−0.0013
CCDE [13]
623.8288
623.8574
623.8904
7000
−0.0077
CSA [15]
623.8684
623.9495
626.3666
10,000
−0.0013
DSPSO-TSA [19]
623.8375
623.8625
623.9001
3000
−0.0063
PSO [19]
624.3045
624.5054
625.9252
3000
0.0685
GA [19]
624.505
624.7419
624.8169
3000
0.1006
TSA [19]
624.3078
635.0623
624.8285
3000
0.0690
KHA1 [23]
611.3276
613.0895
614.8293
10,000
−2.0528
KHA2 [23]
609.0768
610.3271
611.2105
10,000
−2.4299
KHA3 [23]
607.5437
608.1164.
608.5431
10,000
−2.6884
KHA4 [23]
605.7582
605.8043
605.9426
10,000
−2.9911
FA
664.5306
675.5344
679.426
2000
6.1177
IFA
623.8768
625.2704
629.2765
1000The same conclusion is also seen for Subcase 5.3 reported in Table13 since the proposed method can improve the results 0.026% compared to RCGA and 7.284% compared to FA, but the results are 0.022% worse than that of NRCGA. It is clear that the proposed method cannot be better than further improved method like NRCGA. Once again, the proposed method is not effective for system considering VPLE.Table 13
The comparison of costs (in $/h) for Subcase 5.3.
Method
Best cost ($)
Mean cost ($)
Worst cost ($)
NFES
BCI (%)
RCGA [20]
624.6605
625.9201
628.9253
3000
0.026
NRCGA [20]
624.355
624.5792
624.7541
3000
−0.022
FA
673.5544
685.2872
699.2855
6000
7.284
IFA
624.49507
625.2647
629.3951
3000


#### 5.2.6. Comparisons of Test Case6

An IEEE-118 bus system is employed for the test case. The system is similar to a real system with the presence of nearly all electric components such as transformers, capacitors, transmission lines with impedance and admittance, active and reactive power of generators, voltage of load buses and generator buses, the transmission capacity of conductors, and tap setting of transformers. The system has 64 load buses, 186 transmission lines, 54 thermal generating units with quadratic fuel cost function, 14 capacitors, and 9 transformers. The results of minimum cost, mean cost, worst cost andNFES, and the best cost improvement of the proposed method compared to others are reported in Table 14. It can be seen from BCI that MSA [43] is the best method with the best cost improvement of 0.93% compared to the proposed method. On the other hand, the proposed method can improve the best cost from 1.7% (compared to COA) to 10.08% (compared to PSO); meanwhile, its NFES is equal to NFES of PSO and PG-CF-PSO, and much lower than that of all remaining methods. In summary, the proposed method is a promising method for the system with 54 units, quadratic fuel cost function, and all constraints in transmission power networks.Table 14
Result comparisons for Case6.
Method
Best cost ($)
Mean cost ($)
Worst cost ($)
NFES
BCI (%)
PSO [17]
145,520.01
158,596.173
184,686.825
10,000
10.08
PG-CF-PSO [17]
139,604.13
152,204.261
170,022.973
10,000
6.27
COA [39]
133,110.43
153,110.432
138,260.403
—
1.70
ICBOA [38]
135,121.57
—
—
450,000
3.16
CBO [38]
135,073.00
—
—
450,000
3.13
ECBO [38]
135,172.27
—
—
450,000
3.20
DE [38]
142,751.118
—
—
225,000
8.34
ABC [38]
135,145.189
—
—
450,000
3.18
BBOA [38]
135,272.196
—
—
450,000
3.27
MSA [43]
129,640.719
—
—
—
−0.93
FA
136,472.670
139,845.99
146,821.82
10,000
4.12
IFA
130,847.3503
135,391.09
138,816.24
10,000The optimal solutions obtained by the proposed method for Cases1, 2, 4, and 6 and Subcase 5.3 are given in Tables 15–17.Table 15
The optimal solution obtained by the proposed method for Case1, Case 2, and Subcase 5.3.
Case 1
Case 2
Subcase 5.3
Unit
Ps (MW)
Unit
Ps (MW)
Unit
Ps (MW)
1
447.3425
`1
450
1
221.3859
2
173.2211
2
450
2
213.1943
3
263.4287
3
130
3
285.3037
4
138.9756
4
130
4
239.1047
5
165.4668
5
335
5
260
6
87.0112
6
455
6
240.7162
7
465
7
293.1166
8
60
8
239.2326
9
25
9
438.3883
10
20
10
269.5577
11
20
12
55
13
25
14
15
15
15Table 16
The optimal solution obtained by the proposed method for Case4.
Unit
Ps (MW)
Unit
Ps (MW)
Unit
Ps (MW)
Unit
Ps (MW)
1
110.7998
11
168.7998
21
523.2794
31
190
2
110.7998
12
168.7998
22
523.2794
32
190
3
97.3999
13
214.7598
23
523.2794
33
190
4
179.7331
14
394.2794
24
523.2794
34
164.7998
5
92.7582
15
394.2794
25
523.2794
35
164.7998
6
140
16
304.5196
26
523.2794
36
164.7998
7
259.5997
17
489.2794
27
10
37
110
8
284.5997
18
489.2794
28
10
38
110
9
284.5997
19
511.2794
29
10
39
110
10
130
20
511.2794
30
87.7999
40
511.2794Table 17
Optimal solution obtained by the proposed method for Case6.
Variable
Value
Variable
Value
Variable
Value
P1 (MW)
0.8976
P103 (MW)
34.7545
V77 (pu)
1.0133
P4 (MW)
1.7628
P104 (MW)
34.628
V80 (pu)
1.0264
P6 (MW)
0.1324
P105 (MW)
19.5482
V85 (pu)
1.0033
P8 (MW)
12.2424
P107 (MW)
56.655
V87 (pu)
1.0053
P10 (MW)
395.304
P110 (MW)
6.0809
V89 (pu)
0.993
P12 (MW)
76.2481
P111 (MW)
29.5586
V90 (pu)
0.9881
P15 (MW)
10.092
P112 (MW)
15.4431
V91 (pu)
0.9933
P18 (MW)
40.2488
P113 (MW)
0.1302
V92 (pu)
0.9868
P19 (MW)
29.3122
P116 (MW)
30.0209
V99 (pu)
1.0214
P24 (MW)
0.6953
V1 (pu)
1.0489
V100 (pu)
1.0158
P25 (MW)
188.9729
V4 (pu)
1.0345
V103 (pu)
1.0278
P26 (MW)
269.9192
V6 (pu)
1.0035
V104 (pu)
1.0094
P27 (MW)
21.4182
V8 (pu)
0.9753
V105 (pu)
1.0106
P31 (MW)
7.1482
V10 (pu)
1.0454
V107 (pu)
1.0287
P32 (MW)
6.0145
V12 (pu)
1.0275
V110 (pu)
0.996
P34 (MW)
2.3381
V15 (pu)
1.06
V111 (pu)
1.0191
P36 (MW)
0.892
V18 (pu)
1.0618
V112 (pu)
0.9841
P40 (MW)
36.72
V19 (pu)
1.0581
V113 (pu)
1.0346
P42 (MW)
28.3862
V24 (pu)
1.0388
V116 (pu)
1.0382
P46 (MW)
20.5522
V25 (pu)
1.0586
Qc5 (MVAr)
−14.3
P49 (MW)
193.6641
V26 (pu)
1.0993
Qc34 (MVAr)
9.8
P54 (MW)
50.635
V27 (pu)
1.0322
Qc37 (MVAr)
−15.2
P55 (MW)
66.2727
V31 (pu)
1.035
Qc44 (MVAr)
3.9
P56 (MW)
70.6712
V32 (pu)
1.0407
Qc45 (MVAr)
5.5
P59 (MW)
131.0352
V34 (pu)
1.0137
Qc46 (MVAr)
4.7
P61 (MW)
151.3198
V36 (pu)
1.0087
Qc48 (MVAr)
8.5
P62 (MW)
0.9109
V40 (pu)
1.0247
Qc74 (MVAr)
11.3
P65 (MW)
341.1965
V42 (pu)
1.0525
Qc79 (MVAr)
11.8
P66 (MW)
295.7935
V46 (pu)
1.028
Qc82 (MVAr)
0
P70 (MW)
22.3211
V49 (pu)
1.0495
Qc83 (MVAr)
6.3
P72 (MW)
0.0596
V54 (pu)
1.0207
Qc105 (MVAr)
14.8
P73 (MW)
0.2632
V55 (pu)
1.0253
Qc107 (MVAr)
5.3
P74 (MW)
22.6539
V56 (pu)
1.0209
Qc110 (MVAr)
2.6
P76 (MW)
21.9574
V59 (pu)
1.0345
T8 (pu)
0.98
P77 (MW)
47.6506
V61 (pu)
1.0311
T32 (pu)
0.99
P80 (MW)
389.276
V62 (pu)
1.0269
T36 (pu)
1.03
P85 (MW)
0.0008
V65 (pu)
1.0395
T51 (pu)
1.1
P87 (MW)
2.0455
V66 (pu)
1.0644
T93 (pu)
0.91
P89 (MW)
442.1625
V69 (pu)
1.0447
T95 (pu)
1.06
P90 (MW)
41.6673
V70 (pu)
1.0152
T102 (pu)
0.97
P91 (MW)
17.3391
V72 (pu)
1.0264
T107 (pu)
0.9
P92 (MW)
0.0197
V73 (pu)
1.0137
T127 (pu)
1.03
P99 (MW)
0.1407
V74 (pu)
1.0168
V77 (pu)
1.0133
P100 (MW)
212.6985
V76 (pu)
1.0206


## 5.1. Impact Analysis of Each Proposed Improvement

In this section, the contribution of each proposed modification to the performance of the proposed method is investigated by running conventional FA, FA with the first proposed improvement (IFA1), FA with the second proposed improvement (IFA2), FA with the third proposed improvement (IFA3), and the proposed IFA method on five benchmark functions and two power systems. Such benchmark function is given in Table2 [17, 60] while the two power systems with 6 units and 15 units considering POZ constraints are taken from Cases 1 and 2.Table 2
Definition of five applied benchmark functions.
Function
Definition
Domain
f1(Sphere)
f1=∑i=1Dxi2
xi∈−5.12,5.12D=30
f2(Rastrigin)
f2=∑i=1D10+xi2−10cos2πxi
xi∈−5.12,5.12D=30
f3(Ackley)
f3=20+e−20exp−0.21D∑i=1Dxi2−exp1D∑i=1Dcos2πxi
xi∈−30,30D=30
f4(Griewangk)
f4=1+∑i=1Dxi24000−∏i=1Dcosxii
xi∈−600,600D=30
f5(Schwelfel)
f5=∑i=1Dxi+∏iDxi
xi∈−10,10D=30


### 5.1.1. Testing on Benchmark Functions

In order to implement FA, IFA1, IFA2, IFA3, and the proposed method, the population size and number of iterations are set to 50 and 1000, respectively. In addition, we also employ Wilcoxon signed-rank test [61] and set the significance improvement level to 0.01 for calculating Z value, W value, and p value for better comparison between the proposed method and others. As a result, the obtained results in terms of minimum fitness, standard deviation of 50 trial runs, Z value, W value, and p value are reported in Table 3. Among the result numbers, the comparison of the best cost is to reflect the best optimal solution and comparison of the standard deviation cost is to reflect the stabilization of search ability. As pointed in [62], p value can be either equal to or smaller than 0.01 or higher than 0.01. In case that p value of a compared method is higher than 0.01, the significance improvement level of the proposed method over the compared method is not adopted but the significance level is accepted in case such p value is less than 0.01 or equal to 0.01 [62]. Table 3 shows that all minimum values and standard deviation of the proposed method are smaller those of other ones meanwhile Z value, W value, and p value are, respectively, equal to −6.154, 0, and 0. p value = 0 indicates that the proposed method reaches the significance improvement level over other FA variants while W value = 0 determines the accuracy of Z value and p value. For other comparisons among FA, IFA1, IFA2, and IFA3, the minimum fitness and standard deviation reports can indicate that IFA2 is the best modification; IFA3 is the second best modification while IFA1 is only better than conventional FA. The further investigation of performance of each FA methods can be seen clearly in the following section since two power systems with 6 units and 15 units considering POZ constraints are employed for testing.Table 3
Result comparison for five benchmark functions.
Function
Method
Min.
Std. dev.
Z value
W value
P value
f1(Sphere)
FA
47,343
5474.1
−6.154
0
0
IFA1
2381.5
488.8
−6.154
0
0
IFA2
21.77
8.72
−6.154
0
0
IFA3
154.198
24.3
−6.154
0
0
IFA
0.6537
0.1865
f2(Rastrigin)
FA
98.6076
21.2183
−6.154
0
0
IFA1
60.57
125.05
−6.154
0
0
IFA2
6.42
58.08
−6.154
0
0
IFA3
48.58
32.72
−6.154
0
0
IFA
0.17654
11.1474
f3(Ackly)
FA
19.4832
0.2661
−6.154
0
0
IFA1
11.521
1.1903
−6.154
0
0
IFA2
2.7034
0.1023
−6.154
0
0
IFA3
5.577
2.1884
−6.154
0
0
IFA
0.69
0.0481
f4(Griewank)
FA
436.1885
55.157
−6.154
0
0
IFA1
16.67
42.97
−6.154
0
0
IFA2
0.009611
0.165
−6.154
0
0
IFA3
0.33162
5.993
−6.154
0
0
IFA
0.00020344
5.1E − 05
f5(Schwelfel)
FA
49,347
7517.8
−6.154
0
0
IFA1
249.4881
76.65
−6.154
0
0
IFA2
16.19
15.59
−6.154
0
0
IFA3
36.72
36.12
−6.154
0
0
IFA
4.0252
0.292


### 5.1.2. Testing on Power Systems with POZ Constraints

In this section, the performance of FA, IFA1, IFA2, IFA3, and the proposed IFA method is illustrated by running on power systems of Study Cases1 and 2. Furthermore, four other popular methods consisting of particle swarm optimization (PSO) [63], differential evolution (DE) [64], flower pollination algorithm (FPA) [65], and crow search algorithm (CSA) [66] are also implemented for solving two test systems. For running the methods, we set the population and the number of iteration to those values of the proposed method shown in Table 1. In addition, we also tune other control parameters of these methods to their predetermined ranges. For the implementation of PSO, acceleration factors c1 and c2 are set to the range of [0.2, 2.0] with a step size of 0.2. For the implementation of DE, mutation factor and crossover factor are, respectively, set to range of [0.2, 1.2] and [0.2, 1.0] with a step size of 0.2. For the implementation of FPA, probability of using either local search or global search is set to the range of [0.1, 1.0] with a step size of 0.1. For the implementation of CSA, awareness probability and flight length are, respectively, set to the range of [0.2, 1.0] and [0.2, 2.0] with a step size of 0.2. The obtained results from the proposed method and these compared methods such as minimum cost and standard deviation and the results from Wilcoxon signed-rank test such as Z value, W value, and p value are reported in Tables 4 and 5 for Study Cases 1 and 2.Table 4
Result comparisons for Case1.
Method
Best cost ($)
Std. dev. ($)
Z value
W value
p value
FA
15445.215
26.815
−6.154
0
0
IFA1
15444.918
18.01
−6.154
0
0
IFA2
15443.076
1.371
−6.154
0
0
IFA3
15443.257
17.163
−6.154
0
0
PSO
15444.447
16.5186
−6.154
0
0
FPA
15443.55
2.9594
−6.154
0
0
DE
15443.5
26.8126
−6.154
0
0
CSA
15444.75
6.5271
−6.154
0
0
IFA
15443.075
0.050Table 5
Result comparisons for Case2.
Method
Best cost ($)
Std. dev. ($)
Z value
W value
p value
FA
32885.8432
104.534
−6.154
0
0
IFA1
32781.003
88.441
−6.154
0
0
IFA2
32544.9906
8.3011
−6.154
0
0
IFA3
32712.6581
91.8623
−6.154
0
0
PSO
32571.62
94.63377
−6.154
0
0
FPA
32545.1
12.26141
−6.1057
5
0
DE
32635.04
114.94
−6.154
0
0
CSA
32553.03
61.97552
−6.154
0
0
IFA
32544.9704
0.1438For Case1, the proposed method has obtained the best results with less best cost than FA, IFA1, IFA2, and IFA3, respectively, by $2.14, $1.843, $0.001, and $0.182. Similarly, the standard deviation cost of IFA is, respectively, less than that of FA, IFA1, IFA2, and IFA3 by $26.765, $18.05, $1.322, and $17.113. Furthermore, Z value, W value, and p value also indicate the significance improvement of the proposed method over FA variants. Comparison among IFA1, IFA2, and IFA3 can see that such three improved versions has better performance than FA; meanwhile, IFA2 is the best improvement and IFA1 is the worst improvement. Analysis of Case 2 also gives the same evaluation for FA methods.iFor comparison with PSO, FPA, DE, and CSA, all the costs yielded by the proposed method are smaller than those from PSO, FPA, DE, and CSA. In particular, the standard deviation of the proposed method is approximately zero while this value of other methods is much higher than zero. Furthermore,Z value, W value, and p value calculated by setting the significance improvement level to 0.01 also indicate that the proposed method can improve working performance significantly when compared to PSO, FPA, DE, and CSA. Among IFA1, IFA2, and IFA3, only IFA2 outperforms PSO, FPA, DE, and CSA since the best cost and standard deviation cost of IFA2 are smaller than those of PSO, FPA, DE, and CSA. Consequently, it can conclude that the proposed method has high performance for searching optimal solutions while the second modification has the most impact on results.Besides, the more accurate comparison of the best cost improvement and the standard deviation improvement can be quantitative by using the following equation.(39)IFAresultimprovement%=resultofcomparedmethod‐resultofIFAresultofcomparedmethod×100%.In (39), result of IFA or result of compared method can be the best cost, mean cost, worst cost, or standard deviation cost dependent on decision of evaluation but the two most important evaluations must be used including the best cost improvement and the standard deviation improvement. As result is the best cost, IFA result improvement is the improvement of optimal solution quality of IFA over the compared method. In addition, as the result is standard deviation, the IFA result improvement is the improvement of stabilization of IFA over the compared method. As a result, the improvement of IFA over others in terms of the best optimal solution quality and the stabilization is reported in Table 6 for Cases 1 and 2. As observed from Table 6, the best cost improvement of IFA over FA and other methods is not high but the improvement of standard deviation cost is approximately 100% for the two cases.Table 6
Improvement of IFA over FA for Cases1 and 2.
Case
Compared method
Best cost improvement (%)
Std. dev. improvement (%)
1
FA
0.014
99.81
PSO
0.009
99.70
FPA
0.003
98.31
DE
0.003
99.81
CSA
0.011
99.23
2
FA
1.037
99.86
PSO
0.082
99.85
FPA
0.000
98.83
DE
0.276
99.87
CSA
0.025
99.77The best cost of fifty independent trials for Cases1 and 2 obtained by the proposed IFA and FA is, respectively, plotted in Figures 5 and 6. The figures can evaluate the standout of IFA over FA clearly once all minimum costs of IFA are nearly on a flat line and much lower than those of FA, especially for Case 2 while fifty cost values of FA have high fluctuations for the two cases.Figure 5
Fitness-iteration characteristics obtained FA and IFA for Study Case1.Figure 6
Fitness-iteration characteristics obtained by FA and IFA for Study Case2.


## 5.1.1. Testing on Benchmark Functions

In order to implement FA, IFA1, IFA2, IFA3, and the proposed method, the population size and number of iterations are set to 50 and 1000, respectively. In addition, we also employ Wilcoxon signed-rank test [61] and set the significance improvement level to 0.01 for calculating Z value, W value, and p value for better comparison between the proposed method and others. As a result, the obtained results in terms of minimum fitness, standard deviation of 50 trial runs, Z value, W value, and p value are reported in Table 3. Among the result numbers, the comparison of the best cost is to reflect the best optimal solution and comparison of the standard deviation cost is to reflect the stabilization of search ability. As pointed in [62], p value can be either equal to or smaller than 0.01 or higher than 0.01. In case that p value of a compared method is higher than 0.01, the significance improvement level of the proposed method over the compared method is not adopted but the significance level is accepted in case such p value is less than 0.01 or equal to 0.01 [62]. Table 3 shows that all minimum values and standard deviation of the proposed method are smaller those of other ones meanwhile Z value, W value, and p value are, respectively, equal to −6.154, 0, and 0. p value = 0 indicates that the proposed method reaches the significance improvement level over other FA variants while W value = 0 determines the accuracy of Z value and p value. For other comparisons among FA, IFA1, IFA2, and IFA3, the minimum fitness and standard deviation reports can indicate that IFA2 is the best modification; IFA3 is the second best modification while IFA1 is only better than conventional FA. The further investigation of performance of each FA methods can be seen clearly in the following section since two power systems with 6 units and 15 units considering POZ constraints are employed for testing.Table 3
Result comparison for five benchmark functions.
Function
Method
Min.
Std. dev.
Z value
W value
P value
f1(Sphere)
FA
47,343
5474.1
−6.154
0
0
IFA1
2381.5
488.8
−6.154
0
0
IFA2
21.77
8.72
−6.154
0
0
IFA3
154.198
24.3
−6.154
0
0
IFA
0.6537
0.1865
f2(Rastrigin)
FA
98.6076
21.2183
−6.154
0
0
IFA1
60.57
125.05
−6.154
0
0
IFA2
6.42
58.08
−6.154
0
0
IFA3
48.58
32.72
−6.154
0
0
IFA
0.17654
11.1474
f3(Ackly)
FA
19.4832
0.2661
−6.154
0
0
IFA1
11.521
1.1903
−6.154
0
0
IFA2
2.7034
0.1023
−6.154
0
0
IFA3
5.577
2.1884
−6.154
0
0
IFA
0.69
0.0481
f4(Griewank)
FA
436.1885
55.157
−6.154
0
0
IFA1
16.67
42.97
−6.154
0
0
IFA2
0.009611
0.165
−6.154
0
0
IFA3
0.33162
5.993
−6.154
0
0
IFA
0.00020344
5.1E − 05
f5(Schwelfel)
FA
49,347
7517.8
−6.154
0
0
IFA1
249.4881
76.65
−6.154
0
0
IFA2
16.19
15.59
−6.154
0
0
IFA3
36.72
36.12
−6.154
0
0
IFA
4.0252
0.292


## 5.1.2. Testing on Power Systems with POZ Constraints

In this section, the performance of FA, IFA1, IFA2, IFA3, and the proposed IFA method is illustrated by running on power systems of Study Cases1 and 2. Furthermore, four other popular methods consisting of particle swarm optimization (PSO) [63], differential evolution (DE) [64], flower pollination algorithm (FPA) [65], and crow search algorithm (CSA) [66] are also implemented for solving two test systems. For running the methods, we set the population and the number of iteration to those values of the proposed method shown in Table 1. In addition, we also tune other control parameters of these methods to their predetermined ranges. For the implementation of PSO, acceleration factors c1 and c2 are set to the range of [0.2, 2.0] with a step size of 0.2. For the implementation of DE, mutation factor and crossover factor are, respectively, set to range of [0.2, 1.2] and [0.2, 1.0] with a step size of 0.2. For the implementation of FPA, probability of using either local search or global search is set to the range of [0.1, 1.0] with a step size of 0.1. For the implementation of CSA, awareness probability and flight length are, respectively, set to the range of [0.2, 1.0] and [0.2, 2.0] with a step size of 0.2. The obtained results from the proposed method and these compared methods such as minimum cost and standard deviation and the results from Wilcoxon signed-rank test such as Z value, W value, and p value are reported in Tables 4 and 5 for Study Cases 1 and 2.Table 4
Result comparisons for Case1.
Method
Best cost ($)
Std. dev. ($)
Z value
W value
p value
FA
15445.215
26.815
−6.154
0
0
IFA1
15444.918
18.01
−6.154
0
0
IFA2
15443.076
1.371
−6.154
0
0
IFA3
15443.257
17.163
−6.154
0
0
PSO
15444.447
16.5186
−6.154
0
0
FPA
15443.55
2.9594
−6.154
0
0
DE
15443.5
26.8126
−6.154
0
0
CSA
15444.75
6.5271
−6.154
0
0
IFA
15443.075
0.050Table 5
Result comparisons for Case2.
Method
Best cost ($)
Std. dev. ($)
Z value
W value
p value
FA
32885.8432
104.534
−6.154
0
0
IFA1
32781.003
88.441
−6.154
0
0
IFA2
32544.9906
8.3011
−6.154
0
0
IFA3
32712.6581
91.8623
−6.154
0
0
PSO
32571.62
94.63377
−6.154
0
0
FPA
32545.1
12.26141
−6.1057
5
0
DE
32635.04
114.94
−6.154
0
0
CSA
32553.03
61.97552
−6.154
0
0
IFA
32544.9704
0.1438For Case1, the proposed method has obtained the best results with less best cost than FA, IFA1, IFA2, and IFA3, respectively, by $2.14, $1.843, $0.001, and $0.182. Similarly, the standard deviation cost of IFA is, respectively, less than that of FA, IFA1, IFA2, and IFA3 by $26.765, $18.05, $1.322, and $17.113. Furthermore, Z value, W value, and p value also indicate the significance improvement of the proposed method over FA variants. Comparison among IFA1, IFA2, and IFA3 can see that such three improved versions has better performance than FA; meanwhile, IFA2 is the best improvement and IFA1 is the worst improvement. Analysis of Case 2 also gives the same evaluation for FA methods.iFor comparison with PSO, FPA, DE, and CSA, all the costs yielded by the proposed method are smaller than those from PSO, FPA, DE, and CSA. In particular, the standard deviation of the proposed method is approximately zero while this value of other methods is much higher than zero. Furthermore,Z value, W value, and p value calculated by setting the significance improvement level to 0.01 also indicate that the proposed method can improve working performance significantly when compared to PSO, FPA, DE, and CSA. Among IFA1, IFA2, and IFA3, only IFA2 outperforms PSO, FPA, DE, and CSA since the best cost and standard deviation cost of IFA2 are smaller than those of PSO, FPA, DE, and CSA. Consequently, it can conclude that the proposed method has high performance for searching optimal solutions while the second modification has the most impact on results.Besides, the more accurate comparison of the best cost improvement and the standard deviation improvement can be quantitative by using the following equation.(39)IFAresultimprovement%=resultofcomparedmethod‐resultofIFAresultofcomparedmethod×100%.In (39), result of IFA or result of compared method can be the best cost, mean cost, worst cost, or standard deviation cost dependent on decision of evaluation but the two most important evaluations must be used including the best cost improvement and the standard deviation improvement. As result is the best cost, IFA result improvement is the improvement of optimal solution quality of IFA over the compared method. In addition, as the result is standard deviation, the IFA result improvement is the improvement of stabilization of IFA over the compared method. As a result, the improvement of IFA over others in terms of the best optimal solution quality and the stabilization is reported in Table 6 for Cases 1 and 2. As observed from Table 6, the best cost improvement of IFA over FA and other methods is not high but the improvement of standard deviation cost is approximately 100% for the two cases.Table 6
Improvement of IFA over FA for Cases1 and 2.
Case
Compared method
Best cost improvement (%)
Std. dev. improvement (%)
1
FA
0.014
99.81
PSO
0.009
99.70
FPA
0.003
98.31
DE
0.003
99.81
CSA
0.011
99.23
2
FA
1.037
99.86
PSO
0.082
99.85
FPA
0.000
98.83
DE
0.276
99.87
CSA
0.025
99.77The best cost of fifty independent trials for Cases1 and 2 obtained by the proposed IFA and FA is, respectively, plotted in Figures 5 and 6. The figures can evaluate the standout of IFA over FA clearly once all minimum costs of IFA are nearly on a flat line and much lower than those of FA, especially for Case 2 while fifty cost values of FA have high fluctuations for the two cases.Figure 5
Fitness-iteration characteristics obtained FA and IFA for Study Case1.Figure 6
Fitness-iteration characteristics obtained by FA and IFA for Study Case2.


## 5.2. Comparison and Discussion

In order to investigate the further performance of IFA, we carry out the comparisons of results obtained by IFA and many optimization algorithms. In addition to the comparison of the best costs, another comparison criterion is also considered to be the number of fitness evaluationsNFES, which is shown in the following equation.
(40)NFES=ω×Npop×NIter,where ω is the number of generations in each iteration. For some optimization algorithms possessing two new solution generations such as CSA [15] and ICSA [16], ω is 2 while for others possessing one new solution generation like PSO [11] and DE [11], ω is only 1. For the proposed IFA, there is only one new solution generation in each iteration; thus, ω is also equal to 1. As a result, the value of NFES is added into each table for the comparison and the indication is that a method with lower NFES is considered to be more efficient if its best cost is also lower or equal.


### 5.2.1. Comparisons of Test Case1

Table7 reports comparisons of the best cost, mean cost, worst cost, and NFES from IFA and other methods for Case 1, and especially the best cost improvement of IFA compared to each one is also reported at the final column. It notes that if the best cost improvement (BCI) in % has plus sign (+), the result from IFA is better than the compared method, and if BCI has minus sign (−), the compared one has better optimal solution than IFA. Besides, if BCI equals zero, the proposed method and compared methods have the same best cost and the same quality of optimal solution. The observation of BCI can lead to an evaluation that the proposed method is slightly better than most methods (the improvement is from 0.0002% to +0.1030%). KHA4 [23], OBKHA [24], and EMA [30] have the same cost as the proposed method, but α-MFA [53] has slightly lower cost than IFA. However, the recalculated cost of α-MFA [53] is $15445.46, which is higher than the cost of the proposed method. The comparison of NFES can indicate that the proposed method is at least 4 times faster than other ones since NFES of IFA is 300 but that of other ones is from 1200 (α-MFA [53]) to 50,000 (EMA [30]). In spite of using lower value of NFES, the mean cost and the highest cost of the proposed method are also much less or approximately equal to those from other ones. As a result, it can result in a conclusion that the proposed method is much effective than other methods when solving test Case 1 with quadratic fuel cost function, POZ constraint, and power loss constraint.Table 7
Result comparisons for Case1.
Method
Best cost ($)
Mean cost ($)
Worst cost ($)
NFES
BCI (%)
DE [11]
15,449.77
15,449.87
15,449.78
3600
+0.0433
GA [11]
15,459
15,524
15,469
3600
+0.1030
PSO [11]
15,450
15,492
15,454
3600
+0.0448
KHA1 [23]
15,450.75
15,452.82
15,455.50
10,000
+0.0497
KHA2 [23]
15,448.21
15,450.83
15,453.40
10,000
+0.0333
KHA3 [23]
15,445.36
15,447.22
15,449.60
10,000
+0.0148
KHA4 [23]
15,443.08
15,443.19
15,443.30
10,000
0.0000
OBKHA [24]
15,443.08
15,443.92
15,443.33
—
0.0000
BBO [26]
15,443.10
—
—
10,000
+0.0002
EMA [30]
15,443.07
15,443.08
—
50,000
0.0000
AIS [34]
15,448
—
—
4000
+0.0319
NAPSO [18]
15,443.77
15,443.77
15,443.77
3000
+0.0045
NAPSO1 [18]
15,449.14
15,451.04
15,452.92
3000
+0.0393
NAPSO 2 [18]
15,444.77
15,446.77
15,449.77
3000
+0.0110
NAPSO 3 [18]
15,449.08
15,450.33
15,452.78
3000
+0.0389
α-MFA [53]
15442.9
—
—
1200
−0.0011
IFA
15,443.075
15,443.12
15,443.52
300


### 5.2.2. Comparisons of Test Case2

This section compares the performance of the proposed method with other ones when dealing with a larger scale system compared to Case1. The obtained results from the proposed method and other compared methods are given in Table 8. It is obviously seen that ICSA [16] has the best optimal solution; however, the confirmation of a valid solution is not done for ICSA because there was no optimal solution reported for the case in [16]. In comparison with mean cost and maximum cost of other ones, the proposed method has higher values because we run the proposed method with NFES of 1000 while others were run by a very high value. In case that we increase NFES of the proposed method to 6000 and higher than 6000, the result is that the mean cost and the highest cost are approximately equal to the best cost. The illustration aims to demonstrate the fast convergence of the proposed method compared to others, and it is clear that the proposed method is better than this compared method for Case 2 with 15 units considering quadratic fuel cost function and POZ constraint.Table 8
Result comparisons for Case2.
Method
Best cost ($)
Mean cost ($)
Worst cost ($)
NFES
CSA [15]
32,544.9704
32,545.006
32,546.6734
8000
ICSA [16]
32,542.5593
32,543.16884
32,546.6574
6000
IGA [22]
—
32,04.736
—
250,000
IGAMU [22]
—
32,544.990
—
50,000
IFA
32544.9704
32,545.2241
32545.54
1000


### 5.2.3. Comparisons of Test Case3

In the section, a larger system with 20 units considering quadratic fuel cost function and power loss constraint is employed for comparison. The information of the best cost andNFES is shown in Table 9 for comparison with deterministic methods and metaheuristic methods. The best can see that the proposed method can obtain less cost than FA, Lamda, and HM methods in [3] and BBOA in [26] but slightly higher than cost of CSA and ICSA in [16]. The deviation implies that the proposed method is more effective than lambda, HM, and BBOA, but it is not less effective than methods in [16]. Besides, the proposed method has used much less NFES than other ones including FA. As we use a higher number of NFES equal to 10,000, the best cost obtained by the proposed method is also as good as that of methods in [16]. The numerical comparison can show a good performance of the proposed method when dealing with a system with 20 units considering quadratic fuel cost function and power loss constraint.Table 9
Result comparisons for Case3.
Method
Cost ($)
NFES
Lambda [3]
624656.639
—
HM [3]
62456.6341
—
CSA [16]
62456.633
10,000
ICSA [16]
62456.633
10,000
BBOA [26]
62456.7926
20,000
FA
62458.881
10,000
IFA
62456.634
3000


### 5.2.4. Comparisons of Test Case4

Test Case4 is more complicated than above cases due to the consideration of 40 units and nonconvex fuel cost functions. The result comparisons of the case are shown in Table 10. The values of BCI can see that the proposed method has better cost than methods with plus sign “+” of BCI such as BBO, OLCSO, KHA1, KHA2, KHA3, FA [47], and CFA [51]. In comparison with remaining methods, the proposed method has lower quality of solution with higher cost excluding FPA, which did not report optimal solution and had the surprised lowest value. The comparison of NFES also goes against the proposed method in the evaluation of performance. In fact, the proposed method uses a higher number of NFES compared to CCDE, FA, and KHA methods but it uses smaller NFES than the rest of methods including CFA, which is an improved version of FA. The analysis of best cost and NFES can point out that the application of the proposed method is not as effective as the applications with quadratic fuel cost function in Cases 1, 2, and 3 above.Table 10
Result comparisons for Case4.
Method
Best cost ($)
Mean cost ($)
Worst cost ($)
NFES
BCI (%)
CCDE [13]
121,412.58
121,413.56
121,416.53
70,000
−0.002
CSA [16]
121,412.54
121,520.41
121810.25
160,000
−0.002
BBO [26]
121,426.95
121,508.03
121,688.66
200,000
+0.010
FPA [39]
121,074.50
121,095.70
121,196.30
2000
−0.281
OLCSO [44]
121,415.82
121,504.05
121,460.78
400,000
+0.001
FA [47]
121,415.05
121,416.57
121,424.56
50,000
+0.000
KHA1 [23]
121,460.42
121,468.98
121,477.46
10,000
+0.038
KHA2 [23]
121,448.36
121,453.68
121,461.39
10,000
+0.028
KHA3 [23]
121,423.46
121,428.23
121,433.56
10,000
+0.007
KHA4 [23]
121,412.60
121,413.15
121,415.00
10,000
−0.002
CFA [51]
121,415.02
121,416.23
121,420.67
336,000
+0.00035
IFA
121,414.6
121,549.038
121,787.5
100,000
0.000


### 5.2.5. Comparisons of Test Case5

This study case considers different fuel cost functions with nonsmooth forms such as the sum of several quadratic functions for Subcase 5.1, the sum of several nonconvex functions for Subcases 5.2 and 5.3, and different complicated constraints such as POZ, ramp rate limits, active power spinning reserve, and power losses for Subcase 5.3. Subcase 5.1 is the simplest among the three subcases considering four load values, 2400, 2500, 2600, and 2700 MW. Table11 reports the comparison of the best cost and the value of NFES from methods. It can see that the proposed method can yield optimal solutions as good as others excluding FA with the worst results and the proposed method’s NFES is also the lowest value excluding comparison to ALHN and EALHN, which are not metaheuristic methods and face to the difficulties of applying for nonconvex fuel cost function. The comparison can imply that the proposed method is very promising for dealing with ten-unit system considering multifuel options where nonconvex fuel cost function or valve-point loading effects are not taken into account.Table 11
The comparison of the best cost (in $/h) for Subcase 5.1.
Method
PD (MW)
NFES
2400
2500
2600
2700
ALHN [8]
481.723
526.239
574.381
623.809
—
EALHN [9]
481.723
526.239
574.381
623.809
—
DE [11]
481.723
526.239
574.381
623.809
12,000
HDP-ICDE [12]
481.7226
526.2388
574.3808
623.809
8000
HRCGA [21]
481.7226
526.2388
574.3808
623.8092
30,000
AIS [34]
481.723
526.24
574.381
623.8092
4000
CSA [15]
—
—
—
623.8092
6000
CFA [51]
—
—
—
623.8339
156,000
FA
505.2337
580.4417
639.287
679.9525
2000
IFA
481.7226
526.2388
574.3808
623.809
1000Comparisons for Subcases 5.2 reported in Table12 can see that KHA methods in [23] are the best ones with the surprised lowest BCI compared to the proposed method; however, the validation confirmation of solutions cannot be done due to the lack of optimal solution report in [23]. Compared to popular methods such as PSO, GA, TSA (Tabu search algorithm), and FA, the proposed method is better in terms of lower cost and lower NFES. On the contrary, the proposed method is less effective than remaining methods in [13, 15, 16, 19]. As we set NFES to higher values such as 3000, 7000, and 10,000, the results obtained are also worse than these methods. Thus, we conclude that the proposed method is not much effective for the subcase considering valve-point loading effects.Table 12
The comparison of costs (in $/h) for Subcase 5.2.
Method
Best cost ($)
Mean cost ($)
Worst cost ($)
NFES
BCI (%)
ICSA [16]
623.8684
623.9495
626.3666
12,000
−0.0013
CCDE [13]
623.8288
623.8574
623.8904
7000
−0.0077
CSA [15]
623.8684
623.9495
626.3666
10,000
−0.0013
DSPSO-TSA [19]
623.8375
623.8625
623.9001
3000
−0.0063
PSO [19]
624.3045
624.5054
625.9252
3000
0.0685
GA [19]
624.505
624.7419
624.8169
3000
0.1006
TSA [19]
624.3078
635.0623
624.8285
3000
0.0690
KHA1 [23]
611.3276
613.0895
614.8293
10,000
−2.0528
KHA2 [23]
609.0768
610.3271
611.2105
10,000
−2.4299
KHA3 [23]
607.5437
608.1164.
608.5431
10,000
−2.6884
KHA4 [23]
605.7582
605.8043
605.9426
10,000
−2.9911
FA
664.5306
675.5344
679.426
2000
6.1177
IFA
623.8768
625.2704
629.2765
1000The same conclusion is also seen for Subcase 5.3 reported in Table13 since the proposed method can improve the results 0.026% compared to RCGA and 7.284% compared to FA, but the results are 0.022% worse than that of NRCGA. It is clear that the proposed method cannot be better than further improved method like NRCGA. Once again, the proposed method is not effective for system considering VPLE.Table 13
The comparison of costs (in $/h) for Subcase 5.3.
Method
Best cost ($)
Mean cost ($)
Worst cost ($)
NFES
BCI (%)
RCGA [20]
624.6605
625.9201
628.9253
3000
0.026
NRCGA [20]
624.355
624.5792
624.7541
3000
−0.022
FA
673.5544
685.2872
699.2855
6000
7.284
IFA
624.49507
625.2647
629.3951
3000


### 5.2.6. Comparisons of Test Case6

An IEEE-118 bus system is employed for the test case. The system is similar to a real system with the presence of nearly all electric components such as transformers, capacitors, transmission lines with impedance and admittance, active and reactive power of generators, voltage of load buses and generator buses, the transmission capacity of conductors, and tap setting of transformers. The system has 64 load buses, 186 transmission lines, 54 thermal generating units with quadratic fuel cost function, 14 capacitors, and 9 transformers. The results of minimum cost, mean cost, worst cost andNFES, and the best cost improvement of the proposed method compared to others are reported in Table 14. It can be seen from BCI that MSA [43] is the best method with the best cost improvement of 0.93% compared to the proposed method. On the other hand, the proposed method can improve the best cost from 1.7% (compared to COA) to 10.08% (compared to PSO); meanwhile, its NFES is equal to NFES of PSO and PG-CF-PSO, and much lower than that of all remaining methods. In summary, the proposed method is a promising method for the system with 54 units, quadratic fuel cost function, and all constraints in transmission power networks.Table 14
Result comparisons for Case6.
Method
Best cost ($)
Mean cost ($)
Worst cost ($)
NFES
BCI (%)
PSO [17]
145,520.01
158,596.173
184,686.825
10,000
10.08
PG-CF-PSO [17]
139,604.13
152,204.261
170,022.973
10,000
6.27
COA [39]
133,110.43
153,110.432
138,260.403
—
1.70
ICBOA [38]
135,121.57
—
—
450,000
3.16
CBO [38]
135,073.00
—
—
450,000
3.13
ECBO [38]
135,172.27
—
—
450,000
3.20
DE [38]
142,751.118
—
—
225,000
8.34
ABC [38]
135,145.189
—
—
450,000
3.18
BBOA [38]
135,272.196
—
—
450,000
3.27
MSA [43]
129,640.719
—
—
—
−0.93
FA
136,472.670
139,845.99
146,821.82
10,000
4.12
IFA
130,847.3503
135,391.09
138,816.24
10,000The optimal solutions obtained by the proposed method for Cases1, 2, 4, and 6 and Subcase 5.3 are given in Tables 15–17.Table 15
The optimal solution obtained by the proposed method for Case1, Case 2, and Subcase 5.3.
Case 1
Case 2
Subcase 5.3
Unit
Ps (MW)
Unit
Ps (MW)
Unit
Ps (MW)
1
447.3425
`1
450
1
221.3859
2
173.2211
2
450
2
213.1943
3
263.4287
3
130
3
285.3037
4
138.9756
4
130
4
239.1047
5
165.4668
5
335
5
260
6
87.0112
6
455
6
240.7162
7
465
7
293.1166
8
60
8
239.2326
9
25
9
438.3883
10
20
10
269.5577
11
20
12
55
13
25
14
15
15
15Table 16
The optimal solution obtained by the proposed method for Case4.
Unit
Ps (MW)
Unit
Ps (MW)
Unit
Ps (MW)
Unit
Ps (MW)
1
110.7998
11
168.7998
21
523.2794
31
190
2
110.7998
12
168.7998
22
523.2794
32
190
3
97.3999
13
214.7598
23
523.2794
33
190
4
179.7331
14
394.2794
24
523.2794
34
164.7998
5
92.7582
15
394.2794
25
523.2794
35
164.7998
6
140
16
304.5196
26
523.2794
36
164.7998
7
259.5997
17
489.2794
27
10
37
110
8
284.5997
18
489.2794
28
10
38
110
9
284.5997
19
511.2794
29
10
39
110
10
130
20
511.2794
30
87.7999
40
511.2794Table 17
Optimal solution obtained by the proposed method for Case6.
Variable
Value
Variable
Value
Variable
Value
P1 (MW)
0.8976
P103 (MW)
34.7545
V77 (pu)
1.0133
P4 (MW)
1.7628
P104 (MW)
34.628
V80 (pu)
1.0264
P6 (MW)
0.1324
P105 (MW)
19.5482
V85 (pu)
1.0033
P8 (MW)
12.2424
P107 (MW)
56.655
V87 (pu)
1.0053
P10 (MW)
395.304
P110 (MW)
6.0809
V89 (pu)
0.993
P12 (MW)
76.2481
P111 (MW)
29.5586
V90 (pu)
0.9881
P15 (MW)
10.092
P112 (MW)
15.4431
V91 (pu)
0.9933
P18 (MW)
40.2488
P113 (MW)
0.1302
V92 (pu)
0.9868
P19 (MW)
29.3122
P116 (MW)
30.0209
V99 (pu)
1.0214
P24 (MW)
0.6953
V1 (pu)
1.0489
V100 (pu)
1.0158
P25 (MW)
188.9729
V4 (pu)
1.0345
V103 (pu)
1.0278
P26 (MW)
269.9192
V6 (pu)
1.0035
V104 (pu)
1.0094
P27 (MW)
21.4182
V8 (pu)
0.9753
V105 (pu)
1.0106
P31 (MW)
7.1482
V10 (pu)
1.0454
V107 (pu)
1.0287
P32 (MW)
6.0145
V12 (pu)
1.0275
V110 (pu)
0.996
P34 (MW)
2.3381
V15 (pu)
1.06
V111 (pu)
1.0191
P36 (MW)
0.892
V18 (pu)
1.0618
V112 (pu)
0.9841
P40 (MW)
36.72
V19 (pu)
1.0581
V113 (pu)
1.0346
P42 (MW)
28.3862
V24 (pu)
1.0388
V116 (pu)
1.0382
P46 (MW)
20.5522
V25 (pu)
1.0586
Qc5 (MVAr)
−14.3
P49 (MW)
193.6641
V26 (pu)
1.0993
Qc34 (MVAr)
9.8
P54 (MW)
50.635
V27 (pu)
1.0322
Qc37 (MVAr)
−15.2
P55 (MW)
66.2727
V31 (pu)
1.035
Qc44 (MVAr)
3.9
P56 (MW)
70.6712
V32 (pu)
1.0407
Qc45 (MVAr)
5.5
P59 (MW)
131.0352
V34 (pu)
1.0137
Qc46 (MVAr)
4.7
P61 (MW)
151.3198
V36 (pu)
1.0087
Qc48 (MVAr)
8.5
P62 (MW)
0.9109
V40 (pu)
1.0247
Qc74 (MVAr)
11.3
P65 (MW)
341.1965
V42 (pu)
1.0525
Qc79 (MVAr)
11.8
P66 (MW)
295.7935
V46 (pu)
1.028
Qc82 (MVAr)
0
P70 (MW)
22.3211
V49 (pu)
1.0495
Qc83 (MVAr)
6.3
P72 (MW)
0.0596
V54 (pu)
1.0207
Qc105 (MVAr)
14.8
P73 (MW)
0.2632
V55 (pu)
1.0253
Qc107 (MVAr)
5.3
P74 (MW)
22.6539
V56 (pu)
1.0209
Qc110 (MVAr)
2.6
P76 (MW)
21.9574
V59 (pu)
1.0345
T8 (pu)
0.98
P77 (MW)
47.6506
V61 (pu)
1.0311
T32 (pu)
0.99
P80 (MW)
389.276
V62 (pu)
1.0269
T36 (pu)
1.03
P85 (MW)
0.0008
V65 (pu)
1.0395
T51 (pu)
1.1
P87 (MW)
2.0455
V66 (pu)
1.0644
T93 (pu)
0.91
P89 (MW)
442.1625
V69 (pu)
1.0447
T95 (pu)
1.06
P90 (MW)
41.6673
V70 (pu)
1.0152
T102 (pu)
0.97
P91 (MW)
17.3391
V72 (pu)
1.0264
T107 (pu)
0.9
P92 (MW)
0.0197
V73 (pu)
1.0137
T127 (pu)
1.03
P99 (MW)
0.1407
V74 (pu)
1.0168
V77 (pu)
1.0133
P100 (MW)
212.6985
V76 (pu)
1.0206


## 5.2.1. Comparisons of Test Case1

Table7 reports comparisons of the best cost, mean cost, worst cost, and NFES from IFA and other methods for Case 1, and especially the best cost improvement of IFA compared to each one is also reported at the final column. It notes that if the best cost improvement (BCI) in % has plus sign (+), the result from IFA is better than the compared method, and if BCI has minus sign (−), the compared one has better optimal solution than IFA. Besides, if BCI equals zero, the proposed method and compared methods have the same best cost and the same quality of optimal solution. The observation of BCI can lead to an evaluation that the proposed method is slightly better than most methods (the improvement is from 0.0002% to +0.1030%). KHA4 [23], OBKHA [24], and EMA [30] have the same cost as the proposed method, but α-MFA [53] has slightly lower cost than IFA. However, the recalculated cost of α-MFA [53] is $15445.46, which is higher than the cost of the proposed method. The comparison of NFES can indicate that the proposed method is at least 4 times faster than other ones since NFES of IFA is 300 but that of other ones is from 1200 (α-MFA [53]) to 50,000 (EMA [30]). In spite of using lower value of NFES, the mean cost and the highest cost of the proposed method are also much less or approximately equal to those from other ones. As a result, it can result in a conclusion that the proposed method is much effective than other methods when solving test Case 1 with quadratic fuel cost function, POZ constraint, and power loss constraint.Table 7
Result comparisons for Case1.
Method
Best cost ($)
Mean cost ($)
Worst cost ($)
NFES
BCI (%)
DE [11]
15,449.77
15,449.87
15,449.78
3600
+0.0433
GA [11]
15,459
15,524
15,469
3600
+0.1030
PSO [11]
15,450
15,492
15,454
3600
+0.0448
KHA1 [23]
15,450.75
15,452.82
15,455.50
10,000
+0.0497
KHA2 [23]
15,448.21
15,450.83
15,453.40
10,000
+0.0333
KHA3 [23]
15,445.36
15,447.22
15,449.60
10,000
+0.0148
KHA4 [23]
15,443.08
15,443.19
15,443.30
10,000
0.0000
OBKHA [24]
15,443.08
15,443.92
15,443.33
—
0.0000
BBO [26]
15,443.10
—
—
10,000
+0.0002
EMA [30]
15,443.07
15,443.08
—
50,000
0.0000
AIS [34]
15,448
—
—
4000
+0.0319
NAPSO [18]
15,443.77
15,443.77
15,443.77
3000
+0.0045
NAPSO1 [18]
15,449.14
15,451.04
15,452.92
3000
+0.0393
NAPSO 2 [18]
15,444.77
15,446.77
15,449.77
3000
+0.0110
NAPSO 3 [18]
15,449.08
15,450.33
15,452.78
3000
+0.0389
α-MFA [53]
15442.9
—
—
1200
−0.0011
IFA
15,443.075
15,443.12
15,443.52
300


## 5.2.2. Comparisons of Test Case2

This section compares the performance of the proposed method with other ones when dealing with a larger scale system compared to Case1. The obtained results from the proposed method and other compared methods are given in Table 8. It is obviously seen that ICSA [16] has the best optimal solution; however, the confirmation of a valid solution is not done for ICSA because there was no optimal solution reported for the case in [16]. In comparison with mean cost and maximum cost of other ones, the proposed method has higher values because we run the proposed method with NFES of 1000 while others were run by a very high value. In case that we increase NFES of the proposed method to 6000 and higher than 6000, the result is that the mean cost and the highest cost are approximately equal to the best cost. The illustration aims to demonstrate the fast convergence of the proposed method compared to others, and it is clear that the proposed method is better than this compared method for Case 2 with 15 units considering quadratic fuel cost function and POZ constraint.Table 8
Result comparisons for Case2.
Method
Best cost ($)
Mean cost ($)
Worst cost ($)
NFES
CSA [15]
32,544.9704
32,545.006
32,546.6734
8000
ICSA [16]
32,542.5593
32,543.16884
32,546.6574
6000
IGA [22]
—
32,04.736
—
250,000
IGAMU [22]
—
32,544.990
—
50,000
IFA
32544.9704
32,545.2241
32545.54
1000


## 5.2.3. Comparisons of Test Case3

In the section, a larger system with 20 units considering quadratic fuel cost function and power loss constraint is employed for comparison. The information of the best cost andNFES is shown in Table 9 for comparison with deterministic methods and metaheuristic methods. The best can see that the proposed method can obtain less cost than FA, Lamda, and HM methods in [3] and BBOA in [26] but slightly higher than cost of CSA and ICSA in [16]. The deviation implies that the proposed method is more effective than lambda, HM, and BBOA, but it is not less effective than methods in [16]. Besides, the proposed method has used much less NFES than other ones including FA. As we use a higher number of NFES equal to 10,000, the best cost obtained by the proposed method is also as good as that of methods in [16]. The numerical comparison can show a good performance of the proposed method when dealing with a system with 20 units considering quadratic fuel cost function and power loss constraint.Table 9
Result comparisons for Case3.
Method
Cost ($)
NFES
Lambda [3]
624656.639
—
HM [3]
62456.6341
—
CSA [16]
62456.633
10,000
ICSA [16]
62456.633
10,000
BBOA [26]
62456.7926
20,000
FA
62458.881
10,000
IFA
62456.634
3000


## 5.2.4. Comparisons of Test Case4

Test Case4 is more complicated than above cases due to the consideration of 40 units and nonconvex fuel cost functions. The result comparisons of the case are shown in Table 10. The values of BCI can see that the proposed method has better cost than methods with plus sign “+” of BCI such as BBO, OLCSO, KHA1, KHA2, KHA3, FA [47], and CFA [51]. In comparison with remaining methods, the proposed method has lower quality of solution with higher cost excluding FPA, which did not report optimal solution and had the surprised lowest value. The comparison of NFES also goes against the proposed method in the evaluation of performance. In fact, the proposed method uses a higher number of NFES compared to CCDE, FA, and KHA methods but it uses smaller NFES than the rest of methods including CFA, which is an improved version of FA. The analysis of best cost and NFES can point out that the application of the proposed method is not as effective as the applications with quadratic fuel cost function in Cases 1, 2, and 3 above.Table 10
Result comparisons for Case4.
Method
Best cost ($)
Mean cost ($)
Worst cost ($)
NFES
BCI (%)
CCDE [13]
121,412.58
121,413.56
121,416.53
70,000
−0.002
CSA [16]
121,412.54
121,520.41
121810.25
160,000
−0.002
BBO [26]
121,426.95
121,508.03
121,688.66
200,000
+0.010
FPA [39]
121,074.50
121,095.70
121,196.30
2000
−0.281
OLCSO [44]
121,415.82
121,504.05
121,460.78
400,000
+0.001
FA [47]
121,415.05
121,416.57
121,424.56
50,000
+0.000
KHA1 [23]
121,460.42
121,468.98
121,477.46
10,000
+0.038
KHA2 [23]
121,448.36
121,453.68
121,461.39
10,000
+0.028
KHA3 [23]
121,423.46
121,428.23
121,433.56
10,000
+0.007
KHA4 [23]
121,412.60
121,413.15
121,415.00
10,000
−0.002
CFA [51]
121,415.02
121,416.23
121,420.67
336,000
+0.00035
IFA
121,414.6
121,549.038
121,787.5
100,000
0.000


## 5.2.5. Comparisons of Test Case5

This study case considers different fuel cost functions with nonsmooth forms such as the sum of several quadratic functions for Subcase 5.1, the sum of several nonconvex functions for Subcases 5.2 and 5.3, and different complicated constraints such as POZ, ramp rate limits, active power spinning reserve, and power losses for Subcase 5.3. Subcase 5.1 is the simplest among the three subcases considering four load values, 2400, 2500, 2600, and 2700 MW. Table11 reports the comparison of the best cost and the value of NFES from methods. It can see that the proposed method can yield optimal solutions as good as others excluding FA with the worst results and the proposed method’s NFES is also the lowest value excluding comparison to ALHN and EALHN, which are not metaheuristic methods and face to the difficulties of applying for nonconvex fuel cost function. The comparison can imply that the proposed method is very promising for dealing with ten-unit system considering multifuel options where nonconvex fuel cost function or valve-point loading effects are not taken into account.Table 11
The comparison of the best cost (in $/h) for Subcase 5.1.
Method
PD (MW)
NFES
2400
2500
2600
2700
ALHN [8]
481.723
526.239
574.381
623.809
—
EALHN [9]
481.723
526.239
574.381
623.809
—
DE [11]
481.723
526.239
574.381
623.809
12,000
HDP-ICDE [12]
481.7226
526.2388
574.3808
623.809
8000
HRCGA [21]
481.7226
526.2388
574.3808
623.8092
30,000
AIS [34]
481.723
526.24
574.381
623.8092
4000
CSA [15]
—
—
—
623.8092
6000
CFA [51]
—
—
—
623.8339
156,000
FA
505.2337
580.4417
639.287
679.9525
2000
IFA
481.7226
526.2388
574.3808
623.809
1000Comparisons for Subcases 5.2 reported in Table12 can see that KHA methods in [23] are the best ones with the surprised lowest BCI compared to the proposed method; however, the validation confirmation of solutions cannot be done due to the lack of optimal solution report in [23]. Compared to popular methods such as PSO, GA, TSA (Tabu search algorithm), and FA, the proposed method is better in terms of lower cost and lower NFES. On the contrary, the proposed method is less effective than remaining methods in [13, 15, 16, 19]. As we set NFES to higher values such as 3000, 7000, and 10,000, the results obtained are also worse than these methods. Thus, we conclude that the proposed method is not much effective for the subcase considering valve-point loading effects.Table 12
The comparison of costs (in $/h) for Subcase 5.2.
Method
Best cost ($)
Mean cost ($)
Worst cost ($)
NFES
BCI (%)
ICSA [16]
623.8684
623.9495
626.3666
12,000
−0.0013
CCDE [13]
623.8288
623.8574
623.8904
7000
−0.0077
CSA [15]
623.8684
623.9495
626.3666
10,000
−0.0013
DSPSO-TSA [19]
623.8375
623.8625
623.9001
3000
−0.0063
PSO [19]
624.3045
624.5054
625.9252
3000
0.0685
GA [19]
624.505
624.7419
624.8169
3000
0.1006
TSA [19]
624.3078
635.0623
624.8285
3000
0.0690
KHA1 [23]
611.3276
613.0895
614.8293
10,000
−2.0528
KHA2 [23]
609.0768
610.3271
611.2105
10,000
−2.4299
KHA3 [23]
607.5437
608.1164.
608.5431
10,000
−2.6884
KHA4 [23]
605.7582
605.8043
605.9426
10,000
−2.9911
FA
664.5306
675.5344
679.426
2000
6.1177
IFA
623.8768
625.2704
629.2765
1000The same conclusion is also seen for Subcase 5.3 reported in Table13 since the proposed method can improve the results 0.026% compared to RCGA and 7.284% compared to FA, but the results are 0.022% worse than that of NRCGA. It is clear that the proposed method cannot be better than further improved method like NRCGA. Once again, the proposed method is not effective for system considering VPLE.Table 13
The comparison of costs (in $/h) for Subcase 5.3.
Method
Best cost ($)
Mean cost ($)
Worst cost ($)
NFES
BCI (%)
RCGA [20]
624.6605
625.9201
628.9253
3000
0.026
NRCGA [20]
624.355
624.5792
624.7541
3000
−0.022
FA
673.5544
685.2872
699.2855
6000
7.284
IFA
624.49507
625.2647
629.3951
3000


## 5.2.6. Comparisons of Test Case6

An IEEE-118 bus system is employed for the test case. The system is similar to a real system with the presence of nearly all electric components such as transformers, capacitors, transmission lines with impedance and admittance, active and reactive power of generators, voltage of load buses and generator buses, the transmission capacity of conductors, and tap setting of transformers. The system has 64 load buses, 186 transmission lines, 54 thermal generating units with quadratic fuel cost function, 14 capacitors, and 9 transformers. The results of minimum cost, mean cost, worst cost andNFES, and the best cost improvement of the proposed method compared to others are reported in Table 14. It can be seen from BCI that MSA [43] is the best method with the best cost improvement of 0.93% compared to the proposed method. On the other hand, the proposed method can improve the best cost from 1.7% (compared to COA) to 10.08% (compared to PSO); meanwhile, its NFES is equal to NFES of PSO and PG-CF-PSO, and much lower than that of all remaining methods. In summary, the proposed method is a promising method for the system with 54 units, quadratic fuel cost function, and all constraints in transmission power networks.Table 14
Result comparisons for Case6.
Method
Best cost ($)
Mean cost ($)
Worst cost ($)
NFES
BCI (%)
PSO [17]
145,520.01
158,596.173
184,686.825
10,000
10.08
PG-CF-PSO [17]
139,604.13
152,204.261
170,022.973
10,000
6.27
COA [39]
133,110.43
153,110.432
138,260.403
—
1.70
ICBOA [38]
135,121.57
—
—
450,000
3.16
CBO [38]
135,073.00
—
—
450,000
3.13
ECBO [38]
135,172.27
—
—
450,000
3.20
DE [38]
142,751.118
—
—
225,000
8.34
ABC [38]
135,145.189
—
—
450,000
3.18
BBOA [38]
135,272.196
—
—
450,000
3.27
MSA [43]
129,640.719
—
—
—
−0.93
FA
136,472.670
139,845.99
146,821.82
10,000
4.12
IFA
130,847.3503
135,391.09
138,816.24
10,000The optimal solutions obtained by the proposed method for Cases1, 2, 4, and 6 and Subcase 5.3 are given in Tables 15–17.Table 15
The optimal solution obtained by the proposed method for Case1, Case 2, and Subcase 5.3.
Case 1
Case 2
Subcase 5.3
Unit
Ps (MW)
Unit
Ps (MW)
Unit
Ps (MW)
1
447.3425
`1
450
1
221.3859
2
173.2211
2
450
2
213.1943
3
263.4287
3
130
3
285.3037
4
138.9756
4
130
4
239.1047
5
165.4668
5
335
5
260
6
87.0112
6
455
6
240.7162
7
465
7
293.1166
8
60
8
239.2326
9
25
9
438.3883
10
20
10
269.5577
11
20
12
55
13
25
14
15
15
15Table 16
The optimal solution obtained by the proposed method for Case4.
Unit
Ps (MW)
Unit
Ps (MW)
Unit
Ps (MW)
Unit
Ps (MW)
1
110.7998
11
168.7998
21
523.2794
31
190
2
110.7998
12
168.7998
22
523.2794
32
190
3
97.3999
13
214.7598
23
523.2794
33
190
4
179.7331
14
394.2794
24
523.2794
34
164.7998
5
92.7582
15
394.2794
25
523.2794
35
164.7998
6
140
16
304.5196
26
523.2794
36
164.7998
7
259.5997
17
489.2794
27
10
37
110
8
284.5997
18
489.2794
28
10
38
110
9
284.5997
19
511.2794
29
10
39
110
10
130
20
511.2794
30
87.7999
40
511.2794Table 17
Optimal solution obtained by the proposed method for Case6.
Variable
Value
Variable
Value
Variable
Value
P1 (MW)
0.8976
P103 (MW)
34.7545
V77 (pu)
1.0133
P4 (MW)
1.7628
P104 (MW)
34.628
V80 (pu)
1.0264
P6 (MW)
0.1324
P105 (MW)
19.5482
V85 (pu)
1.0033
P8 (MW)
12.2424
P107 (MW)
56.655
V87 (pu)
1.0053
P10 (MW)
395.304
P110 (MW)
6.0809
V89 (pu)
0.993
P12 (MW)
76.2481
P111 (MW)
29.5586
V90 (pu)
0.9881
P15 (MW)
10.092
P112 (MW)
15.4431
V91 (pu)
0.9933
P18 (MW)
40.2488
P113 (MW)
0.1302
V92 (pu)
0.9868
P19 (MW)
29.3122
P116 (MW)
30.0209
V99 (pu)
1.0214
P24 (MW)
0.6953
V1 (pu)
1.0489
V100 (pu)
1.0158
P25 (MW)
188.9729
V4 (pu)
1.0345
V103 (pu)
1.0278
P26 (MW)
269.9192
V6 (pu)
1.0035
V104 (pu)
1.0094
P27 (MW)
21.4182
V8 (pu)
0.9753
V105 (pu)
1.0106
P31 (MW)
7.1482
V10 (pu)
1.0454
V107 (pu)
1.0287
P32 (MW)
6.0145
V12 (pu)
1.0275
V110 (pu)
0.996
P34 (MW)
2.3381
V15 (pu)
1.06
V111 (pu)
1.0191
P36 (MW)
0.892
V18 (pu)
1.0618
V112 (pu)
0.9841
P40 (MW)
36.72
V19 (pu)
1.0581
V113 (pu)
1.0346
P42 (MW)
28.3862
V24 (pu)
1.0388
V116 (pu)
1.0382
P46 (MW)
20.5522
V25 (pu)
1.0586
Qc5 (MVAr)
−14.3
P49 (MW)
193.6641
V26 (pu)
1.0993
Qc34 (MVAr)
9.8
P54 (MW)
50.635
V27 (pu)
1.0322
Qc37 (MVAr)
−15.2
P55 (MW)
66.2727
V31 (pu)
1.035
Qc44 (MVAr)
3.9
P56 (MW)
70.6712
V32 (pu)
1.0407
Qc45 (MVAr)
5.5
P59 (MW)
131.0352
V34 (pu)
1.0137
Qc46 (MVAr)
4.7
P61 (MW)
151.3198
V36 (pu)
1.0087
Qc48 (MVAr)
8.5
P62 (MW)
0.9109
V40 (pu)
1.0247
Qc74 (MVAr)
11.3
P65 (MW)
341.1965
V42 (pu)
1.0525
Qc79 (MVAr)
11.8
P66 (MW)
295.7935
V46 (pu)
1.028
Qc82 (MVAr)
0
P70 (MW)
22.3211
V49 (pu)
1.0495
Qc83 (MVAr)
6.3
P72 (MW)
0.0596
V54 (pu)
1.0207
Qc105 (MVAr)
14.8
P73 (MW)
0.2632
V55 (pu)
1.0253
Qc107 (MVAr)
5.3
P74 (MW)
22.6539
V56 (pu)
1.0209
Qc110 (MVAr)
2.6
P76 (MW)
21.9574
V59 (pu)
1.0345
T8 (pu)
0.98
P77 (MW)
47.6506
V61 (pu)
1.0311
T32 (pu)
0.99
P80 (MW)
389.276
V62 (pu)
1.0269
T36 (pu)
1.03
P85 (MW)
0.0008
V65 (pu)
1.0395
T51 (pu)
1.1
P87 (MW)
2.0455
V66 (pu)
1.0644
T93 (pu)
0.91
P89 (MW)
442.1625
V69 (pu)
1.0447
T95 (pu)
1.06
P90 (MW)
41.6673
V70 (pu)
1.0152
T102 (pu)
0.97
P91 (MW)
17.3391
V72 (pu)
1.0264
T107 (pu)
0.9
P92 (MW)
0.0197
V73 (pu)
1.0137
T127 (pu)
1.03
P99 (MW)
0.1407
V74 (pu)
1.0168
V77 (pu)
1.0133
P100 (MW)
212.6985
V76 (pu)
1.0206


## 6. Conclusions and Future Works

In this paper, solutions of OOTGU problem have been found by using an effective proposed IFA method, which is developed based on conventional FA by performing the three following improvements.(i)
Propose a new radius between two considered solutions.
(ii)
Determine the most appropriate updated step size for each considered solution.
(iii)
Use normal distribution instead of uniform distribution in the equation of producing new solutions.Each proposed improvement’s impact on results is investigated by running on two different systems with different data and constraints. The results of two tests can indicate that among the three improvements, the first proposed improvement has the worst efficiency while the second improvement has the best performance. The combination of the three improvements has the highest impact, resulting in superiorities of IFA over FA as follows.(i)
The proposed method can find much better optimal solutions than FA.
(ii)
The proposed method’s search ability is always more stable than FA reflected via a very small standard deviation.
(iii)
The proposed method converges to an optimal solution much faster than FA because it uses a lower number of produced new solutions.On the other hand, the proposed method is also compared to other existing methods such as original methods and improved methods via testing on six cases with different types of fuel cost functions and different constraints, especially the consideration of all constraints in transmission power networks, which are similar to real power systems. The result comparisons and discussion indicate that the proposed method is an effective method when solving systems with complicated constraints such as POZ and other constraints regarding transmission power networks. However, the proposed method puts up with difficulties of fuel cost functions taking valve-point loading effects into account. Finally, it can lead to a conclusion that the proposed method outperforms popular methods such as PSO, GA, DE, RCGA, and TSA, but it works less effectively than further improved versions of these methods and some other methods.In the paper, we have applied the proposed IFA method for minimizing total electricity generation fuel cost of all available units in cases of neglecting and considering constraints of real transmission power networks. However, the considered optimization problem will be more practical and useful if other popular electricity sources such wind farms [67] and distributed generators [68] as well as maintenance policies of electricity power networks during operation process [69] are taken into account. Thus, our future work is to construct and apply the proposed IFA and other potential algorithms for solving the new problem of power system optimization operation. The objective is also to minimize the fuel cost thermal units and distributed generators while wind energy will be exploited highly and maintenance policies regarding power networks are also considered.



---

*Source: 7267593-2018-07-09.xml*