Noninferior Solution Grey Wolf Optimizer with an Independent Local Search Mechanism for Solving Economic Load Dispatch Problems

The economic load dispatch (ELD) problem is a complex optimization problem in power systems. The main task for this optimization problem is to minimize the total fuel cost of generators while also meeting the conditional constraints of valve-point loading effects, prohibited operating zones, and nonsmooth cost functions. In this paper, a novel grey wolf optimization (GWO), abbreviated as NGWO, is proposed to solve the ELD problem by introducing an independent local search strategy and a noninferior solution neighborhood independent local search technique to the original GWO algorithm to achieve the best problem solution. A local search strategy is added to the standard GWO algorithm in the NGWO, which is called GWOI, to search the local neighborhood of the global optimal point in depth and to guarantee a better candidate. In addition, a noninferior solution neighborhood independent local search method is introduced into the GWOI algorithm to find a better solution in the noninferior solution neighborhood and ensure the high probability of jumping out of the local optimum. The feasibility of the proposed NGWO method is verified on five different power systems, and it is compared with other selected methods in terms of the solution quality, convergence rate, and robustness. The compared experimental results indicate that the proposed NGWO method can efficiently solve ELD problems with higher-quality solutions.


Introduction
Optimization problems widely exist in various fields in real-life. Some of these optimization problems are simple, while others are very complex due to nonconvex objective functions and complex model constraints. For complex optimization problems, a typical characteristic is the minimum or maximum objective function that is subject to heavy equality and/or inequality constraints. The economic load dispatch problem (ELD) is a famous complex power system operation optimization problem. ELD is a computational process, in which the total demanded generation is optimally allocated to each generation unit in operation by minimizing the selected cost criterion while also satisfying the total demand, transmission losses, and a set of physical and operational constraints imposed by the generators and system limitations [1,2]. The optimization study of the ELD problem is In the case of solving the minimum value that is shown in Figure 1, when A is the globally optimal particle, the search with only global guidance leads the particles with poor fitness to gather around A rather than A' and easily fall into the local optimum. The main reason for this phenomenon is the lack of a further in-depth search for the neighborhood of the globally optimal particle A. Therefore, this paper adds the local search strategy to this algorithm and proposed an improved version of the GWO, namely, GWOI to improve the search ability of the standard GWO algorithm. As depicted in Figure 2, after adding the local search strategy into the standard GWO algorithm, GWOI can search the local neighborhood of the global optimal point A in depth and it can easily search for a better candidate A'. The direct method for enhancing its local search ability is to replace the alpha, beta, and delta with the first three optimal individuals of the current iteration in its encircling step formulas, since the GWO algorithm search is mainly controlled by the first three best wolves. In the case of solving the minimum value that is shown in Figure 1, when A is the globally optimal particle, the search with only global guidance leads the particles with poor fitness to gather around A rather than A and easily fall into the local optimum. The main reason for this phenomenon is the lack of a further in-depth search for the neighborhood of the globally optimal particle A. Therefore, this paper adds the local search strategy to this algorithm and proposed an improved version of the GWO, namely, GWOI to improve the search ability of the standard GWO algorithm. As depicted in Figure 2, after adding the local search strategy into the standard GWO algorithm, GWOI can search the local neighborhood of the global optimal point A in depth and it can easily search for a better candidate A . The direct method for enhancing its local search ability is to replace the alpha, beta, and delta with the first three optimal individuals of the current iteration in its encircling step formulas, since the GWO algorithm search is mainly controlled by the first three best wolves. In the case of solving the minimum value that is shown in Figure 1, when A is the globally optimal particle, the search with only global guidance leads the particles with poor fitness to gather around A rather than A' and easily fall into the local optimum. The main reason for this phenomenon is the lack of a further in-depth search for the neighborhood of the globally optimal particle A. Therefore, this paper adds the local search strategy to this algorithm and proposed an improved version of the GWO, namely, GWOI to improve the search ability of the standard GWO algorithm. As depicted in Figure 2, after adding the local search strategy into the standard GWO algorithm, GWOI can search the local neighborhood of the global optimal point A in depth and it can easily search for a better candidate A'. The direct method for enhancing its local search ability is to replace the alpha, beta, and delta with the first three optimal individuals of the current iteration in its encircling step formulas, since the GWO algorithm search is mainly controlled by the first three best wolves. The improved encircling step formulas are as follows: where X a ¢  , X b ¢  , and X¢  are the first three best locally optimal individuals. Therefore, the mathematical models of the search agents update their positions, as follows: The improved encircling step formulas are as follows: where → X α , → X β , and → X are the first three best locally optimal individuals. Therefore, the mathematical models of the search agents update their positions, as follows: Energies 2019, 12, 2274 6 of 26

The Compared GWOII Algorithm
Although the local search strategy is added to the GWOI algorithm to enhance its search ability, it only focuses on the local neighborhood search of the first three best particles. Once the globally optimal particles fall into the local optimum, the algorithm loses the ability to jump out of the local optimum, as shown in Figure 2. If the GWO algorithm can conduct an independent local search in the neighborhood of particles with similar fitness to the current optimal fitness, such as points B, C, D, E, and F, the probability of finding a better solution will greatly increase. This paper proposes a noninferior solution neighborhood independent local search technique based on this analysis. The main idea of this method is as follows: if the fitness value error between some particles and the current optimal particles is small, those particles are considered to be a noninferior solution and will no longer move toward the global optimal particle, but conduct a local depth search in their neighborhood to find a better solution that may exist. Figure 3 shows this search situation. In Figure 3, after introducing the noninferior solution neighborhood independent local search strategy and carrying out a certain number of iterative operations, points B, C, D, E, and F find better points B , C , D , E , and F .

The Proposed NGWO Algorithm
In this subsection, the GWOI algorithm is combined with the GWOII algorithm to form the proposed NGWO algorithm. The NGWO algorithm not only retains the strong search performance of the GWOI algorithm, but it also has the strong ability to jump out of the local optimal solution, as in the GWOII algorithm. Figure 5 shows the superiority of NGWO over GWO, GWOI, and GWOII in terms of optimization performance under certain optimization conditions. In Figure 5, the local depth search method is used to search the neighborhood of the global optimal particle A and to then find a better potential solution A'; the noninferior solution neighborhood independent local search technique is adapted to search the noninferior solution B, C, D, E, and F, and better particles B', C', D', E', and F' are found. Moreover, among these better To implement the noninferior solution neighborhood independent local search technique, the determination conditions of the noninferior solution are as follows [37]: where λ is the adjustment parameter, Fb i (t), F best (t), and Fb(t) represents the best fitness value recorded by the ith particle in the ith iteration, the best fitness value searched by the algorithm so far, and the average value of the optimal fitness value searched by each particle, respectively. If the algorithm satisfies the determination condition of Equation (9) during the iteration, then the corresponding particle is a noninferior solution; then, in the (t + 1)th iteration, the noninferior solution neighborhood independent local search method is executed, and its mathematical expression is as follows: where pbest i (t) is the best position of the ith particle obtained so far, Cauchy(t) is the Cauchy random number in the tth iteration and the reason that we chose this parameter is that it has better stability than the standard normal uniform distribution and is more conducive to the exploration of the algorithms, ub i and lb i are the ith upper boundary and ith lower boundary, respectively, of the search space, and T represents the maximum number of iterations. Therefore, the GWOII algorithm is proposed by introducing the noninferior solution independent local search strategy into the standard GWO algorithm. As described in Figure 4, this strategy enables the GWOII algorithm to find a better C in the neighborhood of noninferior solution C. If the fitness value of C is better than that of the global optimal point A, then the particles in the population will no longer move toward A, but toward C .

The Proposed NGWO Algorithm
In this subsection, the GWOI algorithm is combined with the GWOII algorithm to form the proposed NGWO algorithm. The NGWO algorithm not only retains the strong search performance of the GWOI algorithm, but it also has the strong ability to jump out of the local optimal solution, as in the GWOII algorithm. Figure 5 shows the superiority of NGWO over GWO, GWOI, and GWOII in terms of optimization performance under certain optimization conditions. In Figure 5, the local depth search method is used to search the neighborhood of the global optimal particle A and to then find a better potential solution A'; the noninferior solution neighborhood independent local search technique is adapted to search the noninferior solution B, C, D, E, and F, and better particles B', C', D', E', and F' are found. Moreover, among these better particles, the fitness of C' is better than that of the global optimal particle A'. Therefore, the particles in the population move toward particle C' instead of toward particle A', and particle C' becomes the globally optimal particle.

The Proposed NGWO Algorithm
In this subsection, the GWOI algorithm is combined with the GWOII algorithm to form the proposed NGWO algorithm. The NGWO algorithm not only retains the strong search performance of the GWOI algorithm, but it also has the strong ability to jump out of the local optimal solution, as in the GWOII algorithm. Figure 5 shows the superiority of NGWO over GWO, GWOI, and GWOII in terms of optimization performance under certain optimization conditions.
In Figure 5, the local depth search method is used to search the neighborhood of the global optimal particle A and to then find a better potential solution A ; the noninferior solution neighborhood independent local search technique is adapted to search the noninferior solution B, C, D, E, and F, and better particles B , C , D , E , and F are found. Moreover, among these better particles, the fitness of C is better than that of the global optimal particle A . Therefore, the particles in the population move toward particle C instead of toward particle A , and particle C becomes the globally optimal particle.

The Proposed NGWO Algorithm
In this subsection, the GWOI algorithm is combined with the GWOII algorithm to form the proposed NGWO algorithm. The NGWO algorithm not only retains the strong search performance of the GWOI algorithm, but it also has the strong ability to jump out of the local optimal solution, as in the GWOII algorithm. Figure 5 shows the superiority of NGWO over GWO, GWOI, and GWOII in terms of optimization performance under certain optimization conditions.
In Figure 5, the local depth search method is used to search the neighborhood of the global optimal particle A and to then find a better potential solution A'; the noninferior solution neighborhood independent local search technique is adapted to search the noninferior solution B, C, D, E, and F, and better particles B', C', D', E', and F' are found. Moreover, among these better particles, the fitness of C' is better than that of the global optimal particle A'. Therefore, the particles in the population move toward particle C' instead of toward particle A', and particle C' becomes the globally optimal particle.

Economic Load Dispatch Formulations
The ELD problem can be described as an optimization problem to minimize the total fuel cost of the individual dispatchable generating power while being subject to different constraints. We adopt the problem descriptions and formulations from refs. [38,39].

Objective Function
The ELD problem sums all the costs of the committed generators. Mathematically, this problem can be modeled in Equation (12), as: where F is the total cost function of n committed generating units, ( ) i i F P is the generating cost function of the ith generator with the generation output Pi, and ai, bi, and ci are the smooth cost fuel coefficients of the ith generator, which are constants.

Begin
Initialize the prey wolf population X i (i = 1, 2, . . . , n) and set the maximum number of iterations T. Initialize a, A, and C Calculate fitness function value of each search agent f (X i ) X α = the global best search agent; X α = the local best search agent X β = the global second search agent; X β = the local second search agent X δ = the global third search agent; X δ = the local third search agent while t < T do for each search agent Update the position of the current search agent by Equation (8) end for Update a, A, and C Calculate the fitness function value of all search agents for each search agent if Equation (9) is ture Update the position of the current search agent by Equation (11) end if end for Update X α , X β , X δ , X α , X β and X δ t = t + 1 end while return X α and f (X α )

Economic Load Dispatch Formulations
The ELD problem can be described as an optimization problem to minimize the total fuel cost of the individual dispatchable generating power while being subject to different constraints. We adopt the problem descriptions and formulations from refs. [38,39].

Objective Function
The ELD problem sums all the costs of the committed generators. Mathematically, this problem can be modeled in Equation (12), as: where F is the total cost function of n committed generating units, F i (P i ) is the generating cost function of the ith generator with the generation output P i , and a i , b i , and c i are the smooth cost fuel coefficients of the ith generator, which are constants.
In real-life, the valve-point loading effects are modeled by adding a higher-order nonlinearity rectified sinusoid contribution to the power generating systems and they are represented using Equation (13), as follows: where e i and f i represent the nonsmooth cost fuel coefficients of the ith generator; and, P min i is the minimum generating capacity of the ith generator. According to the description above, the objective function of the ELD problem with the valve-point effect can be formulated as:

Power Balance Constraints and Variables
The whole power demand must equal to the total power generated by available units minus the total transmission loss, which can be modeled as: where P loss and P demand are the value of power demand and whole transmission loss, respectively, in the system. Generally, P loss is calculated by Kron's loss formula, as shown in Equation (17).
where B ij , B i0 , and B 00 are the loss coefficients, which are assumed to be constants under normal circumstances.

Generating Capacity Limits and Variables
The actual output P i that is generated by the ith available unit should range between its minimum generation capacity and maximum generation capacity: where P min i and P max i are the minimum and maximum generating capacity of the ith generator, respectively.

Ramp Rate Limits and Variables
In real circumstances, the operating range of each unit is restricted by its ramp-rate limit constraint: where P and P 0 i are the current and previous power output, respectively, and UR i and DR i are the ramp-up and ramp-down limits of generator i, respectively.

Prohibited Operating Zones Constraints and Variables
In the actual situation, the valve-point loading effects affect the power system, and each generator contains some discontinuous POZs where the generator cannot work. Therefore, the feasible operating zones of each unit should be avoided in these prohibited zones and they can be demonstrated, as follows: are the lower and upper bounds, respectively, of the jth POZs of the ith generating unit, where j ∈ [1, n i ], and n i is the total number of POZs of unit number i.

Implementation of NGWO Method in Solving the ELD Problem
In this subsection, the connection between the NGWO algorithm and the ELD problem was developed to obtain an efficient and high-quality solution. The NGWO algorithm was primarily employed to determine the optimal power generation for each unit that was operational during a particular period to minimize the total power generation cost. Two following definitions should be described in detail before using the proposed NGWO method to solve the ELD problem.

Constraints Handling in ELD Problems with NGWO Approache
The key point in applying the NGWO method to optimizing the ELD problem is how the NGWO algorithm handles the constraints that exist in the problem. In general, most of the researchers are more likely to employ the penalty function methods to treat the constrained optimization problems [44]. The introduction of a penalty function can transform a constrained problem into an unconstrained problem and build a single objective function, so using an unconstrained optimization method can minimize it. When using the NGWO algorithm to solve a constrained ELD problem, it is common to handle constraints using principles of penalty functions, as follows [44]: where, penalty(P i ) is 0 if no constraint is violated; otherwise it is positive value, Fl indicates the feasible region.

ELD Problem without the Valve-Point Loading Effects
In our work, when using the NGWO algorithm to handle the ELD problem without considering the valve-point loading effects, the map methods is built, as in Equation (22).
where, x i is a value between 0 and 1 obtained by the NGWO method, and the meanings of P min i and P max i are shown in Section 3.2. After establishing the map method, Equation (13) is rewritten as: where, q is a positive constants (penalty factors).

ELD Problem with Considering the Valve-Point Loading Effects
In this article, the map method used for handling the valve-point loading effects in NGWO approaches is according to Equations (22)-(24), as follows: where, t is the current iteration, the meanings of P min i , P max i , P 0 i , DR i , and UR i are shown in Section 3.2, x i is a value between 0 and 1 that is obtained by the NGWO method.
Next, if the equality Equation (16) and inequality Equations (19) and (20) are not solved, then Equation (15) is rewritten as: where, q is a positive constant (penalty factors), the meanings of P loss and P demand are shown in Section 3.2.

Implementation Steps of NGWO to ELD Problem
This work presents a quick solution to the ELD problem while utilizing the NGWO algorithm to obtain global optimal or near global optimal generation quantity of each generator unit. The development steps of the proposed technique to solve the ELD problem were detailed, as below.
Step 1: Initialize the population size N and the maximum number of iteration T, and randomly generate the grey individuals of the population between 0 and 1.
Step 2: Calculate the fitness value.
Step 2.1 If it is the ELD problem without the valve-point loading effects, map these initialized grey wolf individuals to the feasible domain of the practical operation constraints according to Equation (22). Calculate the total cost function of n committed generating units by using Equation (23) as the fitness value.
Step 2.2 If the valve-point loading effects are considered in the ELD problem, then map these initialized grey wolf individuals to the feasible domain of the practical operation constraints according to Equation (26), and employ the loss coefficients B, B 0 , and B 00 to calculate the transmission loss P loss while using Equation (17). Calculate the total cost function of n committed generation units using Equation (27) as the fitness value.
Step 3: Compare each individual's fitness value to find out X a , X β , X δ , X α , X β , and X δ .
Step 4: Update the position of the current search agent by Equation (8).
Step 5: Update parameters a, A, and C.
Step 6: Calculate Equation (9). If Equation (9) is satisfied, then update the position of the current search agent by Equation (11).
Step 7: Use Step 2 to calculate fitness function value for each individual search agents.
Step 9: If the number of iterations t reaches the maximum T, then go to Step 10. Otherwise, go to Step 3.
Step 10: The latest generated individual X α is the optimal and then maps X α according to Step 2 to obtain P(X α ). P(X α ) is the optimal generation power of each unit, and its fitness value F(P(X α )) is the minimum total generation cost.
Based on above analysis, the pseudocode of the NGWO algorithm employed to solve the ELD problem is shown in Algorithm 2.

Begin
Initialize the prey wolf population X i (I = 1, 2, . . . , n) and set the maximum number of iterations T. Input the relevant constraint parameters of the generator unit. If the ELD problem has no valve-point loading effects, then map the initialized grey wolf individuals to the feasible domain according to Equation (22) to obtain P(X i ). Otherwise, map the initialized grey wolf individuals to the feasible domain according to Equation (26) to obtain P(X i ), and then calculate the transmission loss P loss by using Equation (17). Initialize a, A, and C Calculate fitness function value of each search agent F(P(X i )) according to Equation (27) for considering the valve-point loading effects. Otherwise, calculate F(P(X i )) by using Equation (23). X α = the global best search agent; X α = the local best search agent X β = the global second search agent; X β = the local second search agent X δ = the global third search agent; X δ = the local third search agent while t < T do for each search agent Update the position of the current search agent by Equation (8) (23). for each search agent if Equation (9) is satisfied Update the position of the current search agent by Equation (11) end if end for Update X α , X β , X δ , X α , X β and X δ t = t + 1 end while return P(X α ) and F(P(X α ))

Numerical Simulation Results and Analysis
To verify the applicability of NGWO for solving the ELD problem, the performance of the basic GWO, the compared GWOI, the compared GWOII, and the proposed NGWO algorithms are assessed on the following ELD cases: Case I. A 3-generator system for load demand of 850 MW, and valve-point loading effects are considered. Case II. A 13-generator system for a load demand of 2520 MW, and valve-point loading effects are considered. Case III. A 40-generator system for a load demand of 10500 MW, and valve-point loading effects are considered. Case IV. A 6-generator system with a quadratic cost function, POZs and transmission loss, and a load demand of 1263 MW. Case V. A 15-generator system with a quadratic cost function, POZs and transmission loss, and a load demand of 2630 MW.
In this paper, the parameters set for each case study mentioned above are listed below. Each optimization technique was coded in MATLAB 2015a and executed on a Windows 10, 4-GHz, 2-GB RAM processor. In addition, the numbers of 50 independent runs were recorded for the compared GWOI, GWOII, and the proposed NGWO algorithms to validate the robustness of the proposed optimization technique. Furthermore, the population size was set to 30 for each case. Finally, for each ELD case, the maximum number of iterations was set to 500.

Case I: 3-Generator System
This case study consists of three generating units with quadratic cost functions and the effects of valve-point loadings are considered [45]. Table 1 provides the data of the generating units. In this case study, NGWO is compared with GWO, GWOI, GWOII, GA, and PSO [46], CJAYA and MP-CJAYA [37], and EP [47], in terms of the mean (F mean ) and the best (F best ) total generation cost. Table 2 records the comparison results. Figure 7 shows the convergence curve of the total generation cost for the mean solution, and Figure 8 explicitly shows the robustness of GWO, GWOI, GWOII, and NGWO in 30 trials. As shown in Table 2, GWO, GWOI, GWOII, and NGWO continuously decrease the values of F best and F mean , and NGWO achieves the very competitive minimum value of 8223.104 $/h relative to that of GA of 8222.07 $/h, as well as the very close minimum F mean value 8233.567 $/h relative to that of MP-CJAYA of 8232.06 $/h. Therefore, the NGWO could obtain the second best results when compared to the above eight mentioned algorithms. The convergence curve in Figure 7 shows that the convergence rates of GWO, GWOI, GWOII, and NGWO continuously improved, and NGWO had the fastest convergence rate. Figure 8 confirms that NGWO achieved the best robustness. In this paper, the parameters set for each case study mentioned above are listed below. Each optimization technique was coded in MATLAB 2015a and executed on a Windows 10, 4-GHz, 2-GB RAM processor. In addition, the numbers of 50 independent runs were recorded for the compared GWOI, GWOII, and the proposed NGWO algorithms to validate the robustness of the proposed optimization technique. Furthermore, the population size was set to 30 for each case. Finally, for each ELD case, the maximum number of iterations was set to 500.

Case I: 3-Generator System
This case study consists of three generating units with quadratic cost functions and the effects of valve-point loadings are considered [45]. Table 1 provides the data of the generating units. In this case study, NGWO is compared with GWO, GWOI, GWOII, GA, and PSO [46], CJAYA and MP-CJAYA [37], and EP [47], in terms of the mean (Fmean) and the best (Fbest) total generation cost. Table 2 records the comparison results. Figure 7 shows the convergence curve of the total generation cost for the mean solution, and Figure 8 explicitly shows the robustness of GWO, GWOI, GWOII, and NGWO in 30 trials. As shown in Table 2, GWO, GWOI, GWOII, and NGWO continuously decrease the values of Fbest and Fmean, and NGWO achieves the very competitive minimum value of 8223.104 $/h relative to that of GA of 8222.07 $/h, as well as the very close minimum Fmean value 8233.567 $/h relative to that of MP-CJAYA of 8232.06 $/h. Therefore, the NGWO could obtain the second best results when compared to the above eight mentioned algorithms. The convergence curve in Figure 7 shows that the convergence rates of GWO, GWOI, GWOII, and NGWO continuously improved, and NGWO had the fastest convergence rate. Figure 8 confirms that NGWO achieved the best robustness.

Case II: 13-Generator System
This case study is the second system with both valve-point loading effects and multiple fuel options and it comprises 13 generating units with quadratic cost functions. Table 3 shows all detailed data, which are taken from refs. [38,43]. The total power demand is 2520 MW. With an increasing number of generators, this system becomes more nonlinear and complex when compared to Case I. The results obtained by GWO, GWOI, GWOII, and NGWO were compared with GA [46], CPSO [48], JAYA and CJAYA [37], and SA [48], as listed in Table 4. The comparison results confirm that NGWO achieved the optimum total generation cost of both Fmean and Fbest among the other algorithms, which were 24,366.12 $/h and 24,185.45 $/h, respectively. Figure 9 shows the convergence curves of GWO, GWOI, GWOII, and NGWO for the mean value of the total generation cost in 30 trials. From Figure 9, GWO and GWOI have better convergence rates in the early iteration than GWOII and NGWO, but they easily fall into the local optimal solution, and GWOII and NGWO more easily obtain the global optimal solution in the later iteration. Figure 10 describes the distribution of the total generation cost of GWO, GWOI, GWOII, and NGWO in 30 trials. NGWO was more robust than GWO, GWOI, and GWOII.

Case II: 13-Generator System
This case study is the second system with both valve-point loading effects and multiple fuel options and it comprises 13 generating units with quadratic cost functions. Table 3 shows all detailed data, which are taken from refs. [38,43]. The total power demand is 2520 MW. With an increasing number of generators, this system becomes more nonlinear and complex when compared to Case I. The results obtained by GWO, GWOI, GWOII, and NGWO were compared with GA [46], CPSO [48], JAYA and CJAYA [37], and SA [48], as listed in Table 4. The comparison results confirm that NGWO achieved the optimum total generation cost of both F mean and F best among the other algorithms, which were 24,366.12 $/h and 24,185.45 $/h, respectively. Figure 9 shows the convergence curves of GWO, GWOI, GWOII, and NGWO for the mean value of the total generation cost in 30 trials. From Figure 9, GWO and GWOI have better convergence rates in the early iteration than GWOII and NGWO, but they easily fall into the local optimal solution, and GWOII and NGWO more easily obtain the global optimal solution in the later iteration. Figure 10 describes the distribution of the total generation cost of GWO, GWOI, GWOII, and NGWO in 30 trials. NGWO was more robust than GWO, GWOI, and GWOII.  "--" indicates the cost value is missing.

Case III: 40-Generator System
In this subsection, the largest ELD problem system consisting of 40 generators, in which the valve-point effect is considered, with a total load demand of 10,500 MW, is selected to investigate the effectiveness of the NGWO algorithm. Table 5 provides the test data for this case study, in which the valve-point loading has also been included in the fuel cost functions [47]. Due to the large number of generator units in this test system, it has a much more complex solution when compared to the previous solution; therefore, the test system is very suitable for testing the difference in the optimization performance of different improvement strategies for the same algorithm. Table 6 reports the comparison results of the proposed NGWO method with NPSO [49], PSO-LRS [37], MPSO [13], CJAYA [37], IGA [50], GWO, GWOI, and GWOII. The table shows that NGWO solved

Case III: 40-Generator System
In this subsection, the largest ELD problem system consisting of 40 generators, in which the valve-point effect is considered, with a total load demand of 10,500 MW, is selected to investigate the effectiveness of the NGWO algorithm. Table 5 provides the test data for this case study, in which the valve-point loading has also been included in the fuel cost functions [47]. Due to the large number of generator units in this test system, it has a much more complex solution when compared to the previous solution; therefore, the test system is very suitable for testing the difference in the optimization performance of different improvement strategies for the same algorithm. Table 6 reports the comparison results of the proposed NGWO method with NPSO [49], PSO-LRS [37], MPSO [13], CJAYA [37], IGA [50], GWO, GWOI, and GWOII. The table shows that NGWO solved the large-scale ELD problem with a high-quality optimum, and its optimization ability was slightly worse than that of CJAYA, but was better than that of the other methods. In these eight comparison algorithms, NGWO achieved the second optimization results. In Figure 11, the convergence curves of four different GWO versions are compared, and it can be observed that NGWO and GWOII dramatically accelerated the convergence rate. However, GWO and GWOI were easily trapped in the local optimum. Figure 12 is the distribution of the optimal total generation cost values that were obtained by GWO, GWOI, and GWOII in 30 runs. This figure demonstrates that these three algorithms show poor stability when optimizing this large-scale power system. However, the stability of the NGWO algorithm is relatively better. All of the above comparisons provide strong evidence demonstrating the effectiveness of NGWO in solving large-scale ELD problems.

Case III: 40-Generator System
In this subsection, the largest ELD problem system consisting of 40 generators, in which the valve-point effect is considered, with a total load demand of 10,500 MW, is selected to investigate the effectiveness of the NGWO algorithm. Table 5 provides the test data for this case study, in which the valve-point loading has also been included in the fuel cost functions [47]. Due to the large number of generator units in this test system, it has a much more complex solution when compared to the previous solution; therefore, the test system is very suitable for testing the difference in the optimization performance of different improvement strategies for the same algorithm. Table 6 reports the comparison results of the proposed NGWO method with NPSO [49], PSO-LRS [37], MPSO [13], CJAYA [37], IGA [50], GWO, GWOI, and GWOII. The table shows that NGWO solved the large-scale ELD problem with a high-quality optimum, and its optimization ability was slightly worse than that of CJAYA, but was better than that of the other methods. In these eight comparison algorithms, NGWO achieved the second optimization results. In Figure 11, the convergence curves of four different GWO versions are compared, and it can be observed that NGWO and GWOII dramatically accelerated the convergence rate. However, GWO and GWOI were easily trapped in the local optimum. Figure 12 is the distribution of the optimal total generation cost values that were obtained by GWO, GWOI, and GWOII in 30 runs. This figure demonstrates that these three algorithms show poor stability when optimizing this large-scale power system. However, the stability of the NGWO algorithm is relatively better. All of the above comparisons provide strong evidence demonstrating the effectiveness of NGWO in solving large-scale ELD problems.

Case IV: 6-Generator System
This case study comprises six generating units with constraints of transmission loss, together with two POZs with ramp-up and ramp-down. The generator data and two POZs are recorded in Table 7, and Table 8 shows the loss B-coefficients [11]. The algorithm used to find the global optima for this problem always encounters challenging complexity, owing to the decision spaces being nonconvex and the cost functions being convex and represented by quadratic functions.
To validate the effectiveness of our method on this test system, NGWO is compared with SA [46], GA [51], MTS [52], NPSO [49], PSO [46], JAYA [37], GWO, GWOI, and GWOII in terms of the total generation cost. Table 9 provides the comparison confirming that GWOI obtained the lowest best total generation cost (F best ) of 15,443.25 $/h among all of the techniques, while JAYA and NGWO achieved the second and third lowest F best of 15,447.09 $/h and 15,449.17 $/h, respectively. In addition, Table 10 summarizes the best total generation cost (F best ), worst total generation cost (F worst ), and mean total generation cost (F mean ) of the four versions of GWO. From  Figure 13 plots the distribution of the total generation cost for the mean solution, which shows that NGWO is the fastest among the four versions of GWO in terms of the convergence rate and it approaches the global optimum. In addition, NGWO has the most robust characteristics, as described in Figure 14.

Case V: 15-Generator System
This case study comprises a larger 15-unit system with quadratic cost functions and it has the same constraints as those in case IV. Table 11 shows the generator data and POZs, in which three POZs exist in generators 2, 5, and 6, and generator 12 has two POZs. The transmission loss coefficient data are taken directly from Ref. [53]. The best output results that were achieved by GWO, GWOI, GWOII, and NGWO are compared with those of SA [46], GA [47], MTS [52], TSA [54], PSO [46], and AIS [55], as recorded in Table 12. The table shows that NGWO obtained the lowest best generation cost among all of the abovementioned methods. Table 13 compares the Fbest, Fworst, and Fmean of the four versions of GWO with those of the techniques that are listed above. As shown in Table 13, NGWO achieved the best Fbest and Fmean and the third best Fworst. GWOII provided the best Fworst, and MTS obtained the second best Fbest, Fworst, and Fmean. Figure 15 shows the convergence curves of average evaluation values of the 15-generator systems while using the four versions of the GWO method. The NGWO was the fastest algorithm to converge to the global optimal solution as can be seen from the simulation. Figure 16 displays the distribution outline of the best solution in 30 runs and, in most of the trials, the NGWO method obtained a better-quality solution and strong, robust characteristics.

Case V: 15-Generator System
This case study comprises a larger 15-unit system with quadratic cost functions and it has the same constraints as those in case IV. Table 11 shows the generator data and POZs, in which three POZs exist in generators 2, 5, and 6, and generator 12 has two POZs. The transmission loss coefficient data are taken directly from Ref. [53]. The best output results that were achieved by GWO, GWOI, GWOII, and NGWO are compared with those of SA [46], GA [47], MTS [52], TSA [54], PSO [46], and AIS [55], as recorded in Table 12. The table shows that NGWO obtained the lowest best generation cost among all of the abovementioned methods. Table 13 compares the Fbest, Fworst, and Fmean of the four versions of GWO with those of the techniques that are listed above. As shown in Table 13, NGWO achieved the best Fbest and Fmean and the third best Fworst. GWOII provided the best Fworst, and MTS obtained the second best Fbest, Fworst, and Fmean. Figure 15 shows the convergence curves of average evaluation values of the 15-generator systems while using the four versions of the GWO method. The NGWO was the fastest algorithm to converge to the global optimal solution as can be seen from the simulation. Figure 16 displays the distribution outline of the best solution in 30 runs and, in most of the trials, the NGWO method obtained a better-quality solution and strong, robust characteristics.

Case V: 15-Generator System
This case study comprises a larger 15-unit system with quadratic cost functions and it has the same constraints as those in case IV. Table 11 shows the generator data and POZs, in which three POZs exist in generators 2, 5, and 6, and generator 12 has two POZs. The transmission loss coefficient data are taken directly from Ref. [53]. The best output results that were achieved by GWO, GWOI, GWOII, and NGWO are compared with those of SA [46], GA [47], MTS [52], TSA [54], PSO [46], and AIS [55], as recorded in Table 12. The table shows that NGWO obtained the lowest best generation cost among all of the abovementioned methods. Table 13 compares the F best , F worst , and F mean of the four versions of GWO with those of the techniques that are listed above. As shown in Table 13, NGWO achieved the best F best and F mean and the third best F worst . GWOII provided the best F worst , and MTS obtained the second best F best , F worst , and F mean . Figure 15 shows the convergence curves of average evaluation values of the 15-generator systems while using the four versions of the GWO method. The NGWO was the fastest algorithm to converge to the global optimal solution as can be seen from the simulation. Figure 16 displays the distribution outline of the best solution in 30 runs and, in most of the trials, the NGWO method obtained a better-quality solution and strong, robust characteristics.

Conclusions and Future Work
In this paper, we successfully applied the proposed NGWO technique to solve the ELD problem by considering the ramp rate limits, POZ constraints, and nonsmooth cost functions. The NGWO algorithm was validated to improve both the global exploration capability and the convergence rates and it had the best robustness when compared to the other three versions of the GWO algorithms. Furthermore, the results, when compared to all the other compared algorithms in five cases, demonstrated the outstanding superiority of the NGWO method in solving the ELD problem.
The superiority of the NGWO algorithm in solving ELD problems was proven. Although our research has not been yet applied to any utility companies or energy providers, we believe that the application of this algorithm in the future will surely improve the operational level of these companies. Next, an interesting application would be to apply the algorithm to training neural networks, optimizing restrictive engineering structures, and solving multiobjective optimization problems. However, for the large-scale power systems of Cases II, III, and V, the robustness of the NGWO algorithm in solving these problems is not as perfect as in solving Cases I and IV, so the optimization performance of the NGWO algorithm can still be improved.