A Novel Constraint Handling Approach for the Optimal Reactive Power Dispatch Problem

: This paper presents an alternative constraint handling approach within a specialized genetic algorithm (SGA) for the optimal reactive power dispatch (ORPD) problem. The ORPD is formulated as a nonlinear single-objective optimization problem aiming at minimizing power losses while keeping network constraints. The proposed constraint handling approach is based on a product of sub-functions that represents permissible limits on system variables and that includes a speciﬁc goal on power loss reduction. The main advantage of this approach is the fact that it allows a straightforward veriﬁcation of both feasibility and optimality. The SGA is examined and tested with the recommended constraint handling approach and the traditional penalization of deviations from feasible solutions. Several tests are run in the IEEE 30, 57, 118 and 300 bus test power systems. The results obtained with the proposed approach are compared to those offered by other metaheuristic techniques reported in the specialized literature. Simulation results indicate that the proposed genetic algorithm with the alternative constraint handling approach yields superior solutions when compared to other recently reported techniques.


Introduction
The optimal reactive power dispatch (ORPD) consists of scheduling available reactive power sources so that operational constraints are met while optimizing a given objective function (typical minimization of power losses or voltage deviation from a desired level). The ORPD plays an important role on the economic and secure operation of power systems. It is a complex combinatorial optimization problem involving a nonlinear objective function, nonlinear constraints and a mixture of continuous and discrete control variables [1]. The control variables of the ORPD are transformer tap settings, generator set points and reactive power compensations. Initial attempts to approach the ORPD problem resorted to linear programming [2], nonlinear programming [3], quadratic programming [4], interior point methods [5], Newton based method [6], dynamic programming [7] and mixed integer programming [8]. Although these techniques are computationally fast they do not perform well when dealing with non-convex problems and discrete variables. Also they tend to converge to local minima and have difficulties handling a large number of decision variables.
The ORPD constitutes an example of a non-convex and multi-modal optimization problem. Due to its nature, several metaheuristic optimization techniques have been tested to solve this problem in the last two decades. The main advantage of metaheuristic techniques is that they can handle both discrete and continuous variables. Furthermore, they do not require differentiability of the objective function or constraints, overcoming the disadvantages of classic optimization algorithms.
In [9] the ORPD is solved by means of Particle Swarm Optimization (PSO). This metaheuristic was first introduced by Eberhart and Kennedy in 1995 [10] and is based on the sociological behaviour associated with bird flocking. Several modifications to improve the performance of this technique have been proposed and some of them have been applied to the ORPD such as parallel vector evaluated PSO algorithm [11] and coordinated aggregation PSO algorithm [12].
Genetic and evolutionary algorithms have also been used to approach the ORPD. These techniques mimic the process of natural selection using the concepts of inheritance, mutation selection and crossover [13,14]. In [15] reactive power optimization is performed by means of a Genetic Algorithm (GA) aiming at minimizing the total support cost from generators and reactive compensators. In [16] a quantum-inspired evolutionary algorithm is developed for real and reactive power optimization. Mean-Variance Mapping Optimization (MVMO) has also been successfully applied to solve the ORPD problem [17]. The working principle of this methodology is based on a special mapping function applied for mutating the offspring on the basis of mean and variance of the set comprising the n-best solutions currently obtained in the algorithm. Hybrid approaches, combining characteristics of two or more metaheuristic techniques are reported in [18,19]. Comparisons of different solution techniques applied to the ORPD problem can be consulted in [20] and [21]. Finally, a review regarding metaheuristic techniques applied to the ORPD problem can be reviewed in [22]. Despite the current trend using novel metaheuristic techniques for solving the ORPD problem, classic metaheuristics such as GA, when properly designed, can be highly competitive.
Genetic and evolutionary algorithms are directly suited to unconstrained optimization. Therefore, the application of such type of algorithms to constraint optimization is a challenging effort. The most common method in GA to handle constraints is the use of penalty functions [23]. In this paper, two different formulations of penalty functions (also called fitness functions) are considered. The first formulation guarantees constraint enforcement by penalizing deviations from the feasible region, and it is the one commonly used in dealing with the ORPD problem. The second one consists of a product of sub-functions which gives the planner the chance to select a specific target on power losses and considers voltage and power flow limits as soft constraints. Such penalty function is devised in such a way that its maximum value is equal to one, only if all voltage magnitudes and power flows are within specified limits and the target on power losses has been achieved. This way, it permits a straightforward verification of both feasibility and optimality. Such penalty function also allows identifying the limit beyond it is not possible to reduce power losses without compromising feasibility. The contributions of this paper are twofold: (a) An alternative constraint handling approach within a specialized genetic algorithm (SGA) is presented for the ORPD problem. The proposed constraint handling approach is based on a product of sub-functions that allows a straightforward verification of both feasibility and optimality. (b) A comparison with other metaheuristic techniques is provided, showing the superiority of the proposed approach. Also, results for the IEEE 300 bus power system (not reported before for the ORPD problem) are reported with the aim of providing solutions for comparative studies in later works.
This paper is organized as follows: Section 1 presents an introduction of the ORPD problem and a literature review regarding the main techniques used to solve it. Section 2 presents the mathematical formulation of the ORPD problem. Section 3 describes the implemented SGA and the alternative constraint handling approach. Section 4 presents the results with IEEE 30, 57, 118 and 300 bus tests systems and a comparison of results with other metaheuristic techniques. Finally, Section 5 presents the conclusions.

Problem Statement
The ORPD has traditionally been solved to reduce active power losses and improvement of voltage profile, subject to various equality and inequality constraints. Power system operators typically include the ORPD problem in operational planning studies rather than in real time applications. The mathematical formulation of the ORPD problem is as follows [22,23].

Objective Function
The objective function considered in this case is the minimization of active power losses given by Equation (1). Where P loss denotes the total active power losses of the transmission network, g k and θ ij are the line conductance and the angular difference of buses i and j, respectively; finally, N K is the total number of network branches: An alternative or complementary objective function is the minimization of absolute value of total voltage deviations (TVD) usually expressed as shown in Equation (2). In this case, N L is the number of load buses in the power system, V i is the voltage magnitude of bus i and V re f i is the voltage magnitude reference of the ith bus (usually 1.0 pu) [24][25][26]: Although the TVD is a commonly used metric to evaluate quality of solutions of the ORPD problem, it only measures the distance of the operating point to a given reference and does not consider the fact that real power systems operate within certain operative limits. In real power systems, a TVD value of zero is not achievable, since it would imply that all voltages are equal to a given reference. A more realistic way of assessing the feasibility of an operative condition is considering an operative range rather than a fixed reference. This alternative is implemented in this paper as explained in Section 3.3.2.

Equality Constraints
The equality constraints of the ORPD problem are the real and reactive power balance Equations which are given by (3) and (4), respectively. In this case, N B is the number of buses; P gi and Q gi are the active and reactive power generation in node i, respectively; P di and Q di are the active and reactive demand in node i, respectively; finally, G ij and B ij are the transfer conductance and susceptance between bus i and bus j, respectively:

Inequality Constraints
Inequality constraints of the ORPD are given by Equations (5)- (9). Superscripts min and max account for minimum and maximum limits of the respective variable.  (5) and (6). In this case, N G is the number of generators in the power system; V gi and Q gi are the voltage magnitude and reactive power of the ith generator, respectively:

Transformer Constraints
Transformers tap settings are bounded by lower and upper constraints as indicated in Equation (7), where N T is the number of transformers with tap setting in the power system:

Shunt VAR Constraints
Shunt VAR compensations are restricted as indicated in Equations (8) and (9), where N C and N L are the number of shunt capacitors and reactors, respectively; while Q ci and Q Li are the reactive power injected by the ith capacitor and reactor, correspondingly:

Security Constraints
These constraints include voltage limits in load buses and transmission line loading as indicated in Equations (10) and (11). In this case, V Li and S li are the voltage magnitude the ith bus and apparent power flow in line li, respectively:

Implemented Genetic Algorithm
Genetic Algorithms are inspired by the mechanisms of natural evolution. They offer an adaptive search based on the Darwinian principle of reproduction and survival of individuals that best adapt themselves to environmental conditions. These algorithms have been successfully applied in optimization problems of great complexity as shown in [27][28][29][30]. As well as other metaheuristic techniques, GAs are commonly used to tackle non-convex multimodal optimization problems, and do not guarantee finding global optimal solutions. The application of basic principles of genetics to mathematical optimization begins with the random or pseudo-random generation of an initial set of solutions (population). The algorithm starts by reading system data and defining the codification of solutions (chromosome). As it will be explained later, the codification was envisaged to take into account real power systems. Then, the SGA parameters are set and an initial population is generated. In this case, it is guaranteed that all candidate solutions are feasible (all control variables are within specified limits). Each individual must be read and decoded by the algorithm indicating the set points of control variables (voltage of generators, transformer taps, capacitors and reactor banks). With this information a power flow is run and power losses are computed. After that, the operators of the SGA are applied (selection, crossover and mutation) until a stopping criterion is met. Further details of

Codification
The proposed codification was devised to be suitable for real power systems. Transformers taps as well as capacitor and reactor banks are discretized based on system data when this one is available or using default parameters when it is not. Figure 2 illustrates the representation of a potential solution to the ORPD problem. It consists of a vector with the discretization of all control variables. Such variables are the setpoints of generators, transformers taps and reactive power injections (from both capacitors and reactors). Control variables are discretized as follows: Voltage setpoints of generators are in a typical range of [0.95, 1.1] pu, coded between discrete values in the range [−100, 100]. However, any other range limit can be considered (depending on specific system data).
Each capacitor is coded using the limits and step size reported for each power system test case. The number of steps for a given capacitor bank is computed using its capacity and step size (if provided). In this way, each capacitor might be coded differently. For example, in the IEEE 30 bus test system all capacitor banks have a maximum capacity of 5 MW; however, the step size is not provided in the original data; in this case, the step was set by default at 0.05 MVAR. Transformers taps vary within the range [ , ] that may be different for every transformer. If the limits of the tap setting and step size are provided in the system data, this information is used in the codification. The number of steps is calculated as the integer number that

Codification
The proposed codification was devised to be suitable for real power systems. Transformers taps as well as capacitor and reactor banks are discretized based on system data when this one is available or using default parameters when it is not. Figure 2 illustrates the representation of a potential solution to the ORPD problem. It consists of a vector with the discretization of all control variables. Such variables are the setpoints of generators, transformers taps and reactive power injections (from both capacitors and reactors). Control variables are discretized as follows: Voltage setpoints of generators are in a typical range of [0.95, 1.1] pu, coded between discrete values in the range [−100, 100]. However, any other range limit can be considered (depending on specific system data).
Each capacitor is coded using the limits and step size reported for each power system test case. The number of steps for a given capacitor bank is computed using its capacity and step size (if provided). In this way, each capacitor might be coded differently. For example, in the IEEE 30 bus test system all capacitor banks have a maximum capacity of 5 MW; however, the step size is not provided in the original data; in this case, the step was set by default at 0.05 MVAR.

Codification
The proposed codification was devised to be suitable for real power systems. Transformers taps as well as capacitor and reactor banks are discretized based on system data when this one is available or using default parameters when it is not. Figure 2 illustrates the representation of a potential solution to the ORPD problem. It consists of a vector with the discretization of all control variables. Such variables are the setpoints of generators, transformers taps and reactive power injections (from both capacitors and reactors). Control variables are discretized as follows: Voltage setpoints of generators are in a typical range of [0.95, 1.1] pu, coded between discrete values in the range [−100, 100]. However, any other range limit can be considered (depending on specific system data).
Each capacitor is coded using the limits and step size reported for each power system test case. The number of steps for a given capacitor bank is computed using its capacity and step size (if provided). In this way, each capacitor might be coded differently. For example, in the IEEE 30 bus test system all capacitor banks have a maximum capacity of 5 MW; however, the step size is not provided in the original data; in this case, the step was set by default at 0.05 MVAR. Transformers taps vary within the range [ , ] that may be different for every transformer. If the limits of the tap setting and step size are provided in the system data, this information is used in the codification. The number of steps is calculated as the integer number that  If the limits of the tap setting and step size are provided in the system data, this information is used in the codification. The number of steps is calculated as the integer number that results from dividing the tap settings range (T max i − T min i ) by the step size; otherwise, a default range of [−10, 10] with steps of 1% is considered.

GA Operators
Initial population is randomly generated within the specified limits of each control variable. This is done to guarantee feasible candidate solutions. After that, the corresponding fitness function of each candidate solution is evaluated. In order to compute power losses, it is necessary to decode and run a power flow for each candidate solution. This is done with the software Matpower 5.1 [31]. Once the fitness function of each candidate solution is calculated, the selection operator is carried out. In this case, selection is performed by tournament method. The recombination or crossover stage combines the information of selected individuals in every subset of control variables (multipoint crossover); Figure 3 illustrates this stage. Mutation rate is dynamic (starts with a high rate and decreases steadily in every generation) and can be applied differently to every subset of control variables. For example, the mutation rate is lower for the subset of capacitors than it is for the subset of transformers tap; consequently, at the end of the evolution process, there is a greater probability of change in transformers taps than in capacitor banks. In this case, the mutated element takes a random value within its limits. This is done to conserve the feasibility of candidate solutions.
In every cycle or generation, the offspring replace the parents only if they represent solutions with better fitness functions. The process of selection, crossover and mutation is repeated until the SGA reaches a specific stopping criterion. Such stopping criterion is determined by a maximum number of generations or when a target on fitness function has been achieved without any violation of system constraints. results from dividing the tap settings range ( − ) by the step size; otherwise, a default range of [−10, 10] with steps of 1% is considered.

GA Operators
Initial population is randomly generated within the specified limits of each control variable. This is done to guarantee feasible candidate solutions. After that, the corresponding fitness function of each candidate solution is evaluated. In order to compute power losses, it is necessary to decode and run a power flow for each candidate solution. This is done with the software Matpower 5.1 [31]. Once the fitness function of each candidate solution is calculated, the selection operator is carried out. In this case, selection is performed by tournament method. The recombination or crossover stage combines the information of selected individuals in every subset of control variables (multipoint crossover); Figure 3 illustrates this stage. Mutation rate is dynamic (starts with a high rate and decreases steadily in every generation) and can be applied differently to every subset of control variables. For example, the mutation rate is lower for the subset of capacitors than it is for the subset of transformers tap; consequently, at the end of the evolution process, there is a greater probability of change in transformers taps than in capacitor banks. In this case, the mutated element takes a random value within its limits. This is done to conserve the feasibility of candidate solutions.
In every cycle or generation, the offspring replace the parents only if they represent solutions with better fitness functions. The process of selection, crossover and mutation is repeated until the SGA reaches a specific stopping criterion. Such stopping criterion is determined by a maximum number of generations or when a target on fitness function has been achieved without any violation of system constraints.

Constraint Handling Approaches
Evolutionary algorithms usually perform unconstrained searches, and thus require additional mechanisms to handle constraints. In the ORPD problem, equality constraints (3) and (4) are met by the load flow solution while constraints on control variables can be handled directly in the problem codification. The remaining constraints to be enforced are voltage magnitudes in load buses and power flow limits in lines (security constraints given by Equations (10) and (11)). These constraints are commonly enforced by some sort of penalty function. Two penalty functions are explored in this paper as detailed below.

Traditional Penalty Function Approach
A penalty function guarantees constraint enforcement by penalizing deviations of candidate solutions from the feasible region of the problem. There are different ways of forming a penalty

Constraint Handling Approaches
Evolutionary algorithms usually perform unconstrained searches, and thus require additional mechanisms to handle constraints. In the ORPD problem, equality constraints (3) and (4) are met by the load flow solution while constraints on control variables can be handled directly in the problem codification. The remaining constraints to be enforced are voltage magnitudes in load buses and power flow limits in lines (security constraints given by Equations (10) and (11)). These constraints are commonly enforced by some sort of penalty function. Two penalty functions are explored in this paper as detailed below.

Traditional Penalty Function Approach
A penalty function guarantees constraint enforcement by penalizing deviations of candidate solutions from the feasible region of the problem. There are different ways of forming a penalty function and several versions of them have been applied in the ORPD problem as reported in [22] and [32]. For comparative purposes, the penalty function approach shown in Equation (12) was selected, which is named as F f 1 (fitness function 1): (12) In this case, x is the general representation of the optimization variables. In (12) the second and third terms correspond to the traditional penalty function approach, where V(x) and P f (x) represent constraint violations on voltage magnitudes in buses and power flows in lines, respectively. µ V and µ P f are penalty constants. Both V(x) and P f (x) are sub-functions that represent the distance to the feasible region of the problem; each of these are expressed in general form as D(x) in Equation (13): where x j , x minj and x maxj represent the optimization variables and their operational limits, respectively. The fitness function used here primarily aims to control the voltage profile and power flow limits. However, it can also be used to handle constraints on other variables such as voltage levels of particular nodes (which cannot operate within conventional ranges), lines with special load capability, etc. Note that V(x) is an alternative way of representing TVD. In this case, this expression considers the fact that voltages operate within a given range and are not compared to a fixed reference.

Alternative Constraint Handling Approach
An alternative constraint handling approach is based on the fitness function shown in Equation (14), named as F f 2 (fitness function 2). This function is an adaptation of the one proposed in [33] which was first introduced in the context of expansion planning for congestion management. In this case f V N (i), f CR (j) and f loss represent sub-functions for voltage magnitude in load bus i, power flow in line j and power losses assessment, respectively. As regards voltage magnitudes F f 2 is used to guarantee that they remain within acceptable limits, rather than comparing them to a fixed reference as done by TVD given by (2). Figure 4 depicts the sub-functions under consideration. Note that f loss allows the planner to set a goal on power loss minimization. In this case, it is assumed that the system operator has a reasonable estimation of the network power losses. Also note that if all quantities are given in per unit, the maximum value of F f 2 is equal to one, independently of the number of constraints. This represents an advantage over traditional penalty functions since it allows the algorithm to stop when the optimal solution (previously selected by the planner) is achieved. It also permits to quickly assess the quality of a given solution which is given by how close F f 2 is to its maximum value. This way, the verification of both feasibility and optimality of candidate solutions is straightforward. Since F f 2 is a multiplication of sub-functions, with only one variable or limit out of bounds the fitness function is penalized indicating that the solution associated to such operational point is neither optimal nor feasible: Energies  The mathematical expressions for the sub-functions depicted in Figure 4 are given by Equations (15)- (17).
and are the maximum and minimum voltage limit on node i, respectively; and are the maximum power flow limit on line j and its actual value, respectively; and represents the goal on system real power losses which is compared to actual power losses . The lambdas in every sub-function determine the hardness of the constraint. Smaller values of lambda indicate softer constraints (see Figure 5).
Note that within sub-function ( ) and ( ) it is possible to set specific voltage limits per node and specific power flow limits per line. This characteristic allows defining nodes and lines with special limits for the ORPD problem.
(a) (b)   The mathematical expressions for the sub-functions depicted in Figure 4 are given by Equations (15)- (17).

Tests and Results
where Vmax i and Vmin i are the maximum and minimum voltage limit on node i, respectively; Load Rjmax and Load Rj are the maximum power flow limit on line j and its actual value, respectively; and loss re f represents the goal on system real power losses which is compared to actual power losses P loss . The lambdas in every sub-function determine the hardness of the constraint. Smaller values of lambda indicate softer constraints (see Figure 5). Note that within sub-function f V N (i) and f CR (j) it is possible to set specific voltage limits per node and specific power flow limits per line. This characteristic allows defining nodes and lines with special limits for the ORPD problem.  The mathematical expressions for the sub-functions depicted in Figure 4 are given by Equations (15)- (17).
and are the maximum and minimum voltage limit on node i, respectively; and are the maximum power flow limit on line j and its actual value, respectively; and represents the goal on system real power losses which is compared to actual power losses . The lambdas in every sub-function determine the hardness of the constraint. Smaller values of lambda indicate softer constraints (see Figure 5).
Note that within sub-function ( ) and ( ) it is possible to set specific voltage limits per node and specific power flow limits per line. This characteristic allows defining nodes and lines with special limits for the ORPD problem.

Tests and Results
To show the applicability of the proposed approach several tests were performed on the IEEE 30, 57, 118 and 300 bus test systems. Specialized literature regarding the ORPD problem usually reports solutions on IEEE 30, 57 and 118 bus test systems. Also, several tests were performed using the IEEE 300 bus test system with the aim of providing solutions for comparative purposes in later works. All tests were carried out on a personal computer equipped with an Intel Core i7 (Quadcore) 3.6 GHz processor and 8 GB of RAM memory. Test system data can be consulted in [34][35][36]. Active and reactive power generation limits as well as active generator settings (except for the swing generator) are taken from [36]. A summary of the test systems data is presented in Table 1.

Input Parameters
Parameters of the SGA used for all simulations are described in Table 2. For simplicity purposes, these set of parameters were tuned to be used with all tests systems. As regards fitness function 2 (given by (14)), it is necessary to set a goal on power losses for every test system. Such goal must be set by the system planner taking into account the particularities of the network. An ambitious goal on power loss reduction might result in unfeasible solutions while a conservative one might result in sub-optimal solutions. Different goals on power loss reduction were tested and those that resulted in feasible solutions are reported in Table 3.  Note that for fitness function 1 (given by (12)) there is no need of setting a specific goal on power losses; since in this case, the algorithm always aims at minimizing losses even at the expense of not fully enforcing security constraints. For both fitness functions voltage limits on load buses were set as V max Li = 1.1 and V min Li = 0.9, µ V and µ P f are 10,000 and 1000, respectively. This values were found through several simulation trials. It was found that lower penalization values would not enforce voltage and power flow limits. Lambdas for fitness function 2 are: Maximum and minimum limits of control variables for the IEEE test cases, with a base of 100 MVA, are given in Table 4 [34][35][36]. Note that the main differences among these cases are maximum limits of capacitor banks and reactors.

Results with the IEEE 30 Bus Power System
The IEEE 30 bus power system comprises nineteen control variables: six generator voltage magnitudes (at buses 1, 2, 5, 8, 11 and 13), four tap changing transformers (at branches 6-9, 6-10, 4-12 and 28-27) and nine shunt capacitor devices (at buses 10, 12, 15,17, 20, 21, 23, 24 and 29). The total system demand is 283.4 MW [34][35][36]. As it is well known, power losses are greatly affected by maximum voltage limits of generators. Allowing higher voltage limits results in lower power losses and vice versa. In regards to the IEEE 30 bus power system, some studies consider upper voltage limits of 1.1 pu while some others consider 1.05 pu. Therefore, several tests were performed considering both limits for comparative purposes. Table 5a,b present the comparison of results when the upper voltage limit of generators is set to 1.1 pu. In Table 5a ABC, FA and HFA stand for artificial bee colony, firefly algorithm and hybrid firefly algorithm, correspondingly; while CLPSO, DE and BBO stand for comprehensive learning particle swarm optimization, differential evolution and biogeography-based optimization, respectively. The solutions obtained with the proposed methodology, using F f 1 and F f 2 , are presented in the last two columns of Table 5b. In this case GSA, MFO and IGSA-CSS stand for gravitational search algorithm, moth-flame optimization and improved gravitational search algorithm with conditional selection strategies. Finally, FAHLCPSO stands for fuzzy adaptive heterogeneous comprehensive-learning particle swarm optimization. Power losses and total voltage deviation (TVD) given by Equations (1) and (2), correspondingly, are also computed for other metaheuristics using the reported values of control variables with the software Matpower 5.1 [31]. Also V(x) and P f (x) are computed as given by Equation (13). Note that both expressions represent the distance to the feasible region for voltage and power low limits, respectively. Power losses obtained with the proposed SGA were 4.5399 MW and 4.5692 MW with total voltage deviations of 2.0105 pu and 1.8333 pu for F f 1 and F f 2 , correspondingly. However, both V(x) and P f (x) are zero, which indicates that the solution found by the SGA guarantees the operation of the system within feasible ranges. When the SGA is implemented with F f 1 , it obtains lower power losses but higher voltage deviations. Note that the SGA outperforms other metaheuristic techniques reported in Table 5a,b when using F f 1 ; however, DE, MFO and BBO obtain slightly better results (with less than 1% of difference) than the proposed methodology when applying F f 2 . Table 5c presents the comparison of results when the upper voltage limit of generators is set to 1.05 pu. In this case, OGSA, ALC-PSO, KHA, CKHA, and NGBWCA stand for opposition-based gravitational search algorithm, particle swarm optimization with an aging leader and challengers, krill heard algorithm, chaotic krill heard algorithm and Gaussian bare-bones water cycle algorithm, respectively. Power losses obtained with the proposed SGA were 5.072 MW using both objective functions. As expected, power losses in this case are higher than those obtained considering higher voltage limits (see Table 5a,b). Nevertheless, the proposed SGA was able to obtain better solutions than those obtained with other metaheuristics.   Figure 6 depicts the convergence of the algorithm for both objective functions for four independent runs (considering 1.1 pu as voltage limit of generators). Note that when using F f 2 the algorithm requires fewer generations to reach convergence, which has a positive impact in computational time.  Figure 6 depicts the convergence of the algorithm for both objective functions for four independent runs (considering 1.1 pu as voltage limit of generators). Note that when using the algorithm requires fewer generations to reach convergence, which has a positive impact in computational time.

Results with the IEEE 57 Bus Power System
The IEEE 57 bus power system consists of eighty branches (lines and transformers), seven generators, fifteen transformers (available for tap changing), and three shunt capacitor devices (at buses 18, 25 and 53). The total system demand is 1250.8 MW [36]. A comparison of the best solutions found with different metaheuristics for the ORPD problem applied to this power system is reported in Table 6a,b with a base of 100 MVA. In this case SOA stands for seeker optimization algorithm. The maximum voltage limit of generators was set to 1.06 pu for comparative purposes. Power losses, voltage deviations, as well as ( ) and ( ) were computed for other metaheuristics using the reported values of control variables. Note that the proposed SGA was able to obtain better results than the other metaheuristics, especially when using ; however, at a expense of higher TVD. Furthermore, the values obtained with the SGA for ( ) and ( ) are approximately zero, meaning that the solution found meets the operational constraints defined for this system, which is not always the case for the other reported metaheuristics. Figure 7 depicts the convergence of the algorithm for both objective functions considering four independent runs. Note that fewer generations are required to reach optimality when is implemented.

Results with the IEEE 57 Bus Power System
The IEEE 57 bus power system consists of eighty branches (lines and transformers), seven generators, fifteen transformers (available for tap changing), and three shunt capacitor devices (at buses 18, 25 and 53). The total system demand is 1250.8 MW [36]. A comparison of the best solutions found with different metaheuristics for the ORPD problem applied to this power system is reported in Table 6a,b with a base of 100 MVA. In this case SOA stands for seeker optimization algorithm. The maximum voltage limit of generators was set to 1.06 pu for comparative purposes. Power losses, voltage deviations, as well as V(x) and P f (x) were computed for other metaheuristics using the reported values of control variables. Note that the proposed SGA was able to obtain better results than the other metaheuristics, especially when using F f 1 ; however, at a expense of higher TVD. Furthermore, the values obtained with the SGA for P f (x) and V(x) are approximately zero, meaning that the solution found meets the operational constraints defined for this system, which is not always the case for the other reported metaheuristics. Figure 7 depicts the convergence of the algorithm for both objective functions considering four independent runs. Note that fewer generations are required to reach optimality when F f 2 is implemented.

Results with the IEEE 118 Bus Power System
The IEEE 118 bus test system has seventy-seven control variables; these consist of fifty-four generator buses, nine tap changing transformers, twelve capacitor devices and two reactor devices. The total system demand is 4242 MW [36]. The optimal settings of control variables are presented in Table 7a

Results with the IEEE 118 Bus Power System
The IEEE 118 bus test system has seventy-seven control variables; these consist of fifty-four generator buses, nine tap changing transformers, twelve capacitor devices and two reactor devices. The total system demand is 4242 MW [36]. The optimal settings of control variables are presented in Table 7a,b; power losses, voltage deviations, V(x) and P f (x) were computed for other metaheuristics using the reported values of control variables. In this case, power losses are given with a base of 100 MVA. Note that the solutions obtained with the proposed approach are better than those reported with other metaheuristics. Furthermore, the values obtained with the SGA for P f (x) and V(x) are zero, which means that the solution found meets all operational constraints. Figure 8 depicts the convergence of the algorithm for both objective functions considering four independent runs. Note that in general fewer generations are needed to reach optimality when using F f 2 .

Results with the IEEE 300 Bus Power System
The IEEE 300 bus test system has one hundred ninety control variables. These consist of sixtynine generator buses, one hundred seven tap changing transformers, eight capacitor devices and six reactor devices. The total system demand is 23525.85 MW [36]. So far, no results for the ORPD problem applied to this system have been reported in the specialized literature. The best control variable settings obtained with the SGA are presented in Table 8.
Power losses are given with a base of 100 MVA and voltage limits on generators are set to 1.1 pu. In this case, a reduction of 9.9% in power losses is obtained when using . Also, note that the values of ( ) and ( ) are approximately zero, which means that the solution found meets the operational constraints. The solution reported in Table 8 can be used for further comparisons in future research. Figure 9 depicts the convergence of the algorithm for both objective functions considering four independent runs. Note that as with the previous systems, the SGA reaches optimality with fewer generations when is implemented.

Results with the IEEE 300 Bus Power System
The IEEE 300 bus test system has one hundred ninety control variables. These consist of sixty-nine generator buses, one hundred seven tap changing transformers, eight capacitor devices and six reactor devices. The total system demand is 23525.85 MW [36]. So far, no results for the ORPD problem applied to this system have been reported in the specialized literature. The best control variable settings obtained with the SGA are presented in Table 8.
Power losses are given with a base of 100 MVA and voltage limits on generators are set to 1.1 pu. In this case, a reduction of 9.9% in power losses is obtained when using F f 2 . Also, note that the values of P f (x) and V(x) are approximately zero, which means that the solution found meets the operational constraints. The solution reported in Table 8 can be used for further comparisons in future research. Figure 9 depicts the convergence of the algorithm for both objective functions considering four independent runs. Note that as with the previous systems, the SGA reaches optimality with fewer generations when F f 2 is implemented.

Results with the IEEE 300 Bus Power System
The IEEE 300 bus test system has one hundred ninety control variables. These consist of sixtynine generator buses, one hundred seven tap changing transformers, eight capacitor devices and six reactor devices. The total system demand is 23525.85 MW [36]. So far, no results for the ORPD problem applied to this system have been reported in the specialized literature. The best control variable settings obtained with the SGA are presented in Table 8.
Power losses are given with a base of 100 MVA and voltage limits on generators are set to 1.1 pu. In this case, a reduction of 9.9% in power losses is obtained when using . Also, note that the values of ( ) and ( ) are approximately zero, which means that the solution found meets the operational constraints. The solution reported in Table 8 can be used for further comparisons in future research. Figure 9 depicts the convergence of the algorithm for both objective functions considering four independent runs. Note that as with the previous systems, the SGA reaches optimality with fewer generations when is implemented.
(a) (b)    Table 9 presents a statistical description of the results obtained with the SGA for all test cases over one hundred runs. Note that using F f 1 yields better results than using F f 2 . The maximum difference of both fitness functions regarding the reduction on power losses is about 2.66%; but significantly different on computation time. Faster results are obtained when using F f 2 . For the IEEE 300 bus test system the reduction on computation time is about 21.25%; however, for the IEEE 57 bus test system this reduction is about 90.37%. This advantage is due to the fact that using F f 2 allows a straightforward verification of both feasibility and optimality. Consequently, the SGA can stop the process when the optimal solution is found even if the maximum number of iterations has not been reached.

Comparison of Fitness Functions Performance
The standard deviations results using F f 2 are smaller than those obtained using F f 1 ; this means that the reproducibility of results is higher when the SGA uses F f 2 . On the other hand, success rate is an indicator of the percentage of runs in which a feasible operational point is obtained before the SGA reaches the maximum number of generations (number of times the algorithm obtains feasible optimal solutions). For example, the last column in Table 9 indicates that in 97 of the 100 runs F f 2 reaches its optimal value before completing the maximum of generations. Fitness function  Figure 10 presents a comparison of system power losses (base case and optimized case) for the test systems under study. Note that for both fitness functions a similar reduction of power losses is achieved, being slightly higher when the SGA is run with F f 1 . These power losses are computed as a percentage of the current active power generation in each test power system.  Table 9 presents a statistical description of the results obtained with the SGA for all test cases over one hundred runs. Note that using yields better results than using . The maximum difference of both fitness functions regarding the reduction on power losses is about 2.66%; but significantly different on computation time. Faster results are obtained when using . For the IEEE 300 bus test system the reduction on computation time is about 21.25%; however, for the IEEE 57 bus test system this reduction is about 90.37%. This advantage is due to the fact that using allows a straightforward verification of both feasibility and optimality. Consequently, the SGA can stop the process when the optimal solution is found even if the maximum number of iterations has not been reached.

Comparison of Fitness Functions Performance
The standard deviations results using are smaller than those obtained using ; this means that the reproducibility of results is higher when the SGA uses . On the other hand, success rate is an indicator of the percentage of runs in which a feasible operational point is obtained before the SGA reaches the maximum number of generations (number of times the algorithm obtains feasible optimal solutions). For example, the last column in Table 9 indicates that in 97 of the 100 runs reaches its optimal value before completing the maximum of generations.  Figure 10 presents a comparison of system power losses (base case and optimized case) for the test systems under study. Note that for both fitness functions a similar reduction of power losses is achieved, being slightly higher when the SGA is run with . These power losses are computed as a percentage of the current active power generation in each test power system.

Conclusions
This paper presented an assessment of two different fitness functions applied to the ORPD within a SGA framework. Such fitness functions represent the classic approach of penalization by adding terms to the fitness function, and a novel approach that consists of the multiplication of different sub-functions representing operative system limits and a goal on power system losses. Although the first approach results in slightly better solutions, it was found that the latter approach not only guarantees the enforcement of network limits but also contributes to a significant reduction of computing time. The main advantage of the proposed fitness function relays on the fact that the optimal solution is known in advance, which was used as stopping criteria for the GA. This fitness function can also be adapted to account for other type of system constraints; such as stability criteria or specific voltage or power flow limits in a given bus or branch. Several tests were performed on the IEEE 30, 57, 118 and 300 bus power systems showing the effectiveness and robustness of the proposed approach. Also, comparisons with different metaheuristic techniques were performed showing the superiority of the proposed approach in terms of quality of solutions.
Author Contributions: All authors contributed to the paper. J.A.V.-V. was the project leader. W.M.V.-A. was responsible for the programing of the GA and running of tests; he also wrote the initial layout of the manuscript; J.M.L.-L. was the advisor in the optimization section and completed the writing of the manuscript. All the authors were responsible for organizing and revising the whole paper.