An improved crow search algorithm applied to the phase swapping problem in asymmetric distribution systems

symmetry S S Article An Improved Crow Search Algorithm Applied to the Phase Swapping Problem in Asymmetric Distribution Systems Brandon Cortés-Caicedo 1, Laura Sofía Avellaneda-Gómez 1, Oscar Danilo Montoya 2,3,* , Lázaro Alvarado-Barrios 4and César Álvarez-Arroyo 5   Citation: Cortés-Caicedo, B.; Avellaneda-Gómez, L.S.; Montoya, O.D.; Alvarado-Barrios, L.; Álvarez-Arroyo, C. An Improved Crow Search Algorithm Applied to the Phase Swapping Problem in Asymmetric Distribution Systems. Symmetry 2021,13, 1329. https:// doi.org/10.3390/sym13081329 Academic Editors: Juan Alberto Rodríguez Velázquez and Alejandro Estrada-Moreno Received: 6 July 2021 Accepted: 21 July 2021 Published: 23 July 2021 Publisher’s Note: MDPI stays neutral with regard to jurisdictional claims in published maps and institutional affiliations. Copyright: © 2021 by the authors. Licensee MDPI, Basel, Switzerland. This article is an open access article distributed under the terms and conditions of the Creative Commons Attribution (CC BY) license (https:// creativecommons.org/licenses/by/ 4.0/). 1Ingeniería Eléctrica, Universidad Distrital Francisco José de Caldas, Bogotá 110231, Colombia; [email protected] (B.C.-C.); [email protected] (L.S.A.-G.) 2Facultad de Ingeniería, Universidad Distrital Francisco José de Caldas, Bogotá 110231, Colombia 3Laboratorio Inteligente de Energía, Universidad Tecnológica de Bolívar, Cartagena 131001, Colombia 4Department of Engineering, Universidad Loyola Andalucía, 41704 Sevilla, Spain; [email protected] 5Department of Electrical Engineering, Universidad de Sevilla, 41092 Sevilla, Spain; [email protected] *Correspondence: [email protected] Abstract: This paper discusses the power loss minimization problem in asymmetric distribution systems (ADS) based on phase swapping. This problem is presented using a mixed-integer nonlinear programming model, which is resolved by applying a master–slave methodology. The master stage consists of an improved version of the crow search algorithm. This stage is based on the generation of candidate solutions using a normal Gaussian probability distribution. The master stage is responsible for providing the connection settings for the system loads using integer coding. The slave stage uses a power flow for ADSs based on the three-phase version of the iterative sweep method, which is used to determine the network power losses for each load connection supplied by the master stage. Numerical results on the 8-, 25-, and 37-node test systems show the efficiency of the proposed approach when compared to the classical version of the crow search algorithm, the Chu and Beasley genetic algorithm, and the vortex search algorithm. All simulations were obtained using MATLAB and validated in the DigSILENT power system analysis software. Keywords: improved crow search algorithm; normal Gaussian distribution; phase swapping problem; power losses; asymmetric distribution grids; vortex search algorithm 1. Introduction Due to the economic and population growth, the dependence on electrical systems has equally grown to satisfy humanity’s basic needs, changing the habits and customs of how individuals live and work [ 1 ]. To ensure this, three-phase distribution networks are used, which are responsible for interconnecting transmission and sub-transmission networks with end-users (i.e., residential, industrial, and commercial areas) requiring medium and low voltage [ 2 , 3 ]. These systems generally operate in an asymmetric manner due to the following factors. (i) The configurations on the distribution lines are asymmetrical since the transposition criterion is not applicable due to the short length of the lines [ 4 , 5 ]. (ii) The nature of the loads may be 1 ϕ , 2 ϕ , or 3 ϕ , which generates unbalances in voltages at the nodes and in the line currents [6]. (iii) The arbitrary location of single-phase transformers on the phases of the system causes an unbalance in the currents through the lines [ 7 ]. Load unbalances in distribution systems create undesirable scenarios such as the increase of current in any phase system, which produces an increase in power losses through its constituent elements [ 8 ]. These power losses can exceed the capacity required to supply the demand, cause equipment to age, and increase investment and operating costs for network operators [7,9]. The importance of reducing power losses in distribution networks has established multiple approaches, such as (i) optimal placement and sizing of distributed generation [10] , Symmetry 2021,13, 1329. https://doi.org/10.3390/sym13081329 https://www.mdpi.com/journal/symmetry Symmetry 2021,13, 1329 2 of 20 (ii) optimal capacitor placement and sizing [ 11 ], (iii) optimal network reconfiguration [ 12 ], (iv) optimal conductor sizing in distribution networks [ 13 ], (v) optimal power system restoration [ 14 ], and (vi) optimal phase swapping [ 15 , 16 ]. These strategies can significantly help distribution companies to reduce the number of power losses. However, the first two approaches involve significant investments since they integrate new devices into the distribution network [ 17 ]. The third approach requires less investment since few distribution lines need to be constructed to realize the optimal network reconfiguration [ 18 ]. The fourth approach also requires high investment as the system conductors need to be renewed [ 19 ]. The fifth methodology is more appropriate for operation of the power system after fault isolation. The sixth strategy is the most economical as it requires few teams to reconfigure the system loads without investing in new equipment [ 20 ]. Bearing in mind the low phase swapping costs to minimize the power losses in ADSs, a new master–slave optimization strategy is proposed to solve this problem. In the specialized literature, the balance phase problem, with the minimizing power losses approach, has been solved using different optimization methods, including the Chu and Beasley genetic algorithms [ 8 , 16 , 21 – 24 ], particle swarm optimization [ 9 ], mixedinteger convex optimization [ 25 ], bat optimization algorithm [ 26 ], differential evolution algorithm [27] , simulated annealing optimizer [ 28 ], and vortex search algorithm [ 15 ], among others. The main feature of the optimization methodologies described above is that they employ the master–slave optimization scheme to solve the problem [ 15 ]. The master defines the connection of the loads to the nodes. The slave is typically a power flow tool that allows one to revise and exploit the solution space through the power losses calculation [16]. Similar to the metaheuristic optimization methods described above, a master–slave methodology is proposed in this work to solve the phase swapping problem in ADSs. The proposed optimization algorithm corresponds to an improved version of the crow optimization algorithm (CSA) to select the connection of the loads in the master stage, together with the use of the iterative sweep power flow method in its three-phase version in the slave stage. In the master stage, the connection of each load is defined using an integer encoding between 1 and 6, which represents the six possible connection forms for a three-phase charge [ 8 ]. The slave stage is responsible for evaluating the power flow to determine the total power losses for the connection set provided in the master stage [ 15 ]. Improvements in the classical CSA are carried out in the crow avoidance stage based on a probability criterion [ 29 ]. If the probability is higher than the crow knowledge probability ( Ap ), the new crow position i is provided using the classical CSA exploration proposed in [29] . Likewise, if the possibility is less than the crow knowledge probability (i.e., Ap ), the new crow position i is generated through a regular Gaussian distribution (GD) used in the process of evolution of the vortex search algorithm (VSA) [ 30 ]. The main benefit of the VSA is that the solution space can be explored and exploited through the use of hyper-spheres derived from the selection of an individual from the current population. It also clarifies that the criterion of evolution in our proposal is applied at each iteration, which implies that this process is of the adaptive type. The main contributions of our proposal are listed below. • It proposes an improved approach for the classical CSA using the VSA evolution mechanism to revise and exploit the solution space. • The interaction between the improved CSA (i.e., ICSA) and the three-phase power flow (TPPF), based on the classical iterative sweep method, allows the application of phase swapping in radial or meshed systems with connected loads, either in Yor ∆. It is relevant to mention that, upon analyzing the specialized literature, there was no evidence of the CSA application to the phase swapping problem in distribution systems, which corresponds to a research gap that this work intends to fill. In addition, the numerical results obtained in test systems of 8, 25, and 37 nodes prove the quality of the algorithm when compared with classical metaheuristic optimization methodologies. Symmetry 2021,13, 1329 3 of 20 The structure of the remainder of the document takes the following form: Section 2 presents the optimal phase swapping problem representation in ADSs, Section 3presents the ICSA incorporated with the TPPF method, Section 4portrays the electrical networks used in this research, and Section 5represents the obtained results for the connections set and the grid power losses. Finally, Section 6states the conclusions drawn from the development of this article. 2. Mathematical Formulation The optimal phase swapping problem in ADSs is represented by a mixed-integer nonlinear programming (MINLP) model [ 16 ]. The binary variables are the decision variables, which correspond to the connections set for each load presented in the system [ 8 ]. Additionally, the power flow formulation provides the continuous part of the decision variables. The nonlinear nature of the products appears between the different voltage magnitudes at the nodes and the trigonometric functions [ 25 ]. Next, the objective function and the set of constraints representing the phase swapping problem are presented. 2.1. Formulation of the Objective Function The phase swapping problem has an objective function associated with the minimization of total active power losses of the ADS, as presented in Equation (1): min z=∑ n∈N ∑ m∈N ∑ f∈F ∑ g∈F Yn f mgVn f Vmgcos(δn f −δmg −θnm f g), (1) where z defines the value of the objective function. Further, Yn f mg is the admittance magnitude associated with node n in the electrical phase f with node m in the electrical phase g , Vn f (Vmg) corresponds to the voltage magnitude at node n(m) in the electrical phase f(g) , δn f (δmg) represents the angle of the voltage at node n(m) in the electrical phase f(g) , and θnm f g represents the admittance angle associated with node n in the electrical phase f with node m in the electrical phase g . It is relevant to mention that F and N are the sets containing all phases and nodes, respectively. Remark 1. The product between the magnitudes of the voltages and the trigonometric functions makes the objective function nonlinear and nonconvex [ 16 ]. The structure of the objective function makes advanced numerical optimization techniques necessary to minimize it efficiently [ 15 ]. A master–slave methodology is proposed, as is the case of the improved version of the developed CSA, due to its simplicity in programming terms. 2.2. Set of Constraints The phase swapping problem has a set of constraints that corresponds to the different operating limitations in an ADS [15]. These are shown from Equation (2) to Equation (6): Ps n f −∑ g∈F xn f gPd ng =Vn f ∑ m∈N ∑ g∈F Yn f mgVmgcos(δn f −δmg −θn f mg),∀f∈ F ∀n∈ N , (2) Qs n f −∑ g∈F xn f gQd ng =Vn f ∑ m∈N ∑ g∈F Yn f mgVmgsin(δn f −δmg −θn f mg),∀f∈ F ∀n∈ N , (3) ∑ g∈F xn f g =1, ∀f∈ F,∀n∈ N , (4) ∑ f∈F xn f g =1, ∀g∈ F,∀n∈ N , (5) Vmin ≤Vn f ≤Vmax,∀g∈ F,∀n∈ N , (6) where Ps n f is the variable associated with the active power produced at generator s connected to node n in the electrical phase f , and Qs n f is the generated reactive power at Symmetry 2021,13, 1329 4 of 20 generator s sited at node n in the electrical phase f . Pd ng indicates the active demanded power at node n in the electrical phase g , while Qd ng describes the reactive power needed at node n in the electrical phase g . xn f g is a binary variable defining the configuration of demand node n at f in the electrical phase g . Finally, Vmin and Vmax correspond to the allowable limits of voltage regulation for all the system nodes. Remark 2. Equations (2) and (3) are nonlinear and nonconvex. These features highlight the complexity of the TPPF problem for electrical systems, making it necessary to use numerical methods, such as the iterative sweep method, to solve it. Note that Equation (1) defines the form of the objective function for the phase swapping problem formulated as the minimization of power losses under a given demand condition in all network sections of the system. Equations (2) and (3) define the apparent power balance constraints maintained at each phase and node of the ADS. Equations (4) and (5) ensure that loads take a unique connection form by using a matrix of connections (i.e., xn f g ) at each node [ 5 ]. Finally, the constraint presented in (6) defines the allowable limits of voltage regulation for all nodes of the system [15]. 3. Methodology Proposed To solve the optimal phase swapping problem in ADS with the objective of minimize power losses for a specific demand condition, this paper proposes to use an ICSA [ 29 ] as the master stage in conjunction with the iterative swept TPPF as the slave stage [ 15 ]. The master stage defines the phase configurations at each system demand node to achieve the most balanced system possible, while the slave stage evaluates the power flow constraints defined in (2) and (3). The section below will describe each component of the proposed master–slave methodology. 3.1. Slave Stage: TPPF Method The iterative sweep power flow method is a numerical method typically used for single-phase distribution networks [ 31 ]. Nevertheless, this method has been adapted for three-phase ADSs with wye (i.e., Y ) and delta (i.e., ∆ ) loads [ 15 ]. This method is derived from graph theory, where the topology of the network is represented by an incidence matrix, which relates the nodes and the links of the system [ 31 ]. First, Kirchhoff’s first law is used to calculate currents in the system nodes, starting from the terminal nodes and until the source node, which corresponds to the implementation of the backward sweep stage. Then, Kirchhoff’s second law, which corresponds to the implementation of the forward sweep stage, is used to calculate the voltage drops in the network sections from the slack node to the terminal nodes [31]. One of the most important aspects of the iterative sweep power flow method is that it is derivative-free. Likewise, the matrices involved in the calculations are constant, which implies that the computing times required to obtain a solution are in milliseconds [15]. In order to expose the iterative swept TPPF method developed by [ 15 ], in any n− node system, we used the relationship between the nodal voltage and the injected current that is presented using the equivalent between the admittance matrix and the incidence matrix [32], as shown in Equation (7). Ig3ϕ Id3ϕ="Ag3ϕZ−1 r3ϕAT g3ϕAg3ϕZ−1 r3ϕAT d3ϕ Ad3ϕZ−1 r3ϕAT g3ϕAd3ϕZ−1 r3ϕAT d3ϕ#Vg3ϕ Vd3ϕ, (7) where Vg3ϕ is the vector containing all the voltages at the slack node that are known for power flow purposes [ 21 ]. Vd3ϕ is the vector containing all the unknown variables of interest, i.e., the demand voltages. Further, Ig3ϕ represents the vector with the net current injections at the slack node, Id3ϕ represents the vector that involves all the currents at the nodes of consumption, and Zr3ϕ is the matrix that contains all the impedance matrices of Symmetry 2021,13, 1329 5 of 20 the distribution lines present in the system. Ag3ϕ is the component of the incidence matrix that associates the source nodes to each other, while Ad3ϕ is the component of the incidence matrix relating the nodes of consumption to each other. Remark 3. The voltage variables and current ones in Equation (7) are organized by nodes and phases according to the three-phase condition. From Equation (7), it can be observed that the second row has the nodal voltages at the demand nodes (i.e., Vd3ϕ ), which are the unknown variables in power flow studies [ 16 ]. Equation (7) can then be rewritten as follows: (8). Vd3ϕ=−Y−1 dd3ϕYdg3ϕVg3ϕ−Id3ϕ, (8) where Ydg3ϕ=Ad3ϕZ−1 r3ϕAT g3ϕand Ydd3ϕ=Ad3ϕZ−1 r3ϕAT d3ϕ. Equation (8) allows the determination of all the nodal demand voltages per phase. However, it is necessary to consider the type of load, either Y or ∆ , to establish the demand current (i.e., Id3ϕ). In the case where node m has a constant power load with a Y structure (assuming it is solidly earthed [ 33 ]), the demand current can be shown as in Equation (9), as reported in [21]. Idm3ϕ=−diag−1(V∗ dm3ϕ)S∗ dm3ϕ(9) If node m has a load with a connection ∆ , the demand current can be as shown in Equation (10), as reported in [21]. Idm3ϕ=−(diag−1(MV∗ dm3ϕ)−diag−1(MTV∗ dm3ϕ)H)S∗ dm3ϕ, (10) where the Hand Mmatrices are defined as follows: H=  001 100 010 ,M=  1−1 0 0 1 −1 −1 0 1   Through a t iteration counter, the solution of Equation (8) is obtained if maxn||Vt+1 d3ϕ|−|Vt d3ϕ||o≤e , where e is the maximum tolerance suggested as 1×10−10 [34]. When solving the TPPF, the main objective is to evaluate the power losses for the phase connection set established in the master stage. For this purpose, Equation (11), as described in [32], is used. Ploss =realnJ3ϕ TZr3ϕJ3ϕo, (11) where Ploss describes the total system effective power losses and J3ϕ represents the current per phase flowing through the system branches expressed as shown in Equation (12) through Ohm’s Law, as reported in [15]. J3ϕ=Z−1 r3ϕE3ϕ, (12) where E3ϕ represents the voltage drop per phase in the system branches, which can be written in terms of the generation and demand using the three-phase incidence matrix as shown in Equation (13), as reported in [15]. E3ϕ=AT g3ϕVg3ϕ+AT d3ϕVd3ϕ(13) Algorithm 1shows the general implementation of the TPPF method by an iterative sweep for ADSs with connected loads in Yand ∆. Symmetry 2021,13, 1329 6 of 20 Algorithm 1: Solution of the TPPF problem for ADSs with Yand ∆loads Define the characteristics of the unbalanced three-phase system under study; Obtain the per-unit equivalent of the system; Generate the 3ϕincidence matrix A3ϕ; Extract the components Ag3ϕand Ad3ϕ; Define Zr3ϕ; Calculate Ydg3ϕand Ydd3ϕ; Select tmax; Define e; Chose the voltages per phase of the Slack node: Vg3ϕ=1∠0, 1∠−2π 3, 1∠2π 3T; Do t=0; for m≥n−1do Do Vt dm3ϕ=Vg3ϕ; end for t≤tmax do Define m=1; for m≥n−1do if node load m is connected at Y then Calculate It dm3ϕusing Equation (9); else Calculate It dm3ϕusing Equation (10); end end Calculate the new voltages at the demand nodes Vt d3ϕusing Equation (8); if maxn||Vt+1 d3ϕ|−|Vt d3ϕ||o≤ethen Report the nodal voltages as V3ϕ=Vg3ϕ,Vd3ϕT; Calculate the voltage drop across the branches of the system using Equation (13); Calculate the current flowing in the system branches using Equation (12); Calculate the power losses using Equation (11) ; break; else Do Vt d3ϕ=Vt+1 d3ϕ; end end Remark 4. The convergence of the matrix iterative sweep method can be demonstrated with the Banach fixed-point theorem, as reported in [ 31 ]. So, if the system is far enough from the stress collapse point, it can be guaranteed that the solution of Equations (2) and (3) obtain any combination of nodal loads provided by the master stage. 3.2. Master Stage: ICSA The master stage is responsible for providing the nodal connections set for evaluation in the iterative sweep TPPF presented in Algorithm 2. This paper proposes an improved version of the classical CSA modifying the solution space exploration by introducing a normal GD employed by the VSA optimization method [ 35 ]. Before explaining the improvements made to the CSA, the encoding of the phase swapping problem in ADSs according to the possible configurations is shown in Table 1. Symmetry 2021,13, 1329 7 of 20 Table 1. Possible connection types for the system loads [8]. Connection Type Phases Sequence 1 UVW 2 VWU No vary 3 WUV 4 UWV 5 WVU Vary 6 VUW The coding used to represent an individual iin iteration kis presented as follows: Xi,k=[6, 2, 4, 5, ..., c, ..., 1] Here, the dimensions are 1 ×(n− 1 ) and c is an integer that defines the type of connection (see Table 1) [ 15 ]. The ICSA is developed using a probability criterion at each iteration, which will define the methodology to be used to define the new site of crow i . If the probability criterion is greater than the crow knowledge probability (i.e., AP ) at iteration k, the traditional classical CSA search is employed [29]; otherwise, it will use the GD of the VSA to generate the new crow position i. 3.2.1. Classical Approach: CSA CSA is a metaheuristic optimization technique inspired by the intelligent performance of crows [ 29 ]. Crows are considered to be the smartest birds in nature; they possess a brain much larger in relation to the size of their bodies [ 36 ]. In groups, crows show notable traits of intelligence. They can learn and remember faces, use tools, communicate using sophisticated manners, and manage their food throughout the seasons because they hide their excess food in certain places (caches) in their environment and retrieve it when necessary [29]. Crows are ambitious birds as they pursue each other for better food reserves, observe where other birds hide their food, and steal it once they have left [ 37 ]. If a crow has stolen something, it takes additional cautions such as relocating its hiding places to prevent becoming a future victim [ 36 ]. They use their experience in theft to predict another robber’s behavior and to determine the quickest course to protect their caches from being stolen [ 38 ]. The main bases of CSA are listed below [ 29 ]: (i) crows live in flocks, (ii) crows remember the site of their hiding places, (iii) crows pursue others to commit robberies, and (iv) crows shield their caches from being robbed using probability. The following is an example of how the CSA mechanism works. Crows explore and exploit their environment, which is the solution space. Each environment cache corresponds to a feasible solution, the quality of the food source is the objective function. The best food source in the environment is the solution to the problem, Thus, the CSA seeks to simulate the intelligent behavior of crows to find the optimal solution [29]. The CSA is also based on a population. The population size (i.e., flock) consists of N individuals (number of crows), which belong to a d dimensional solution space, where d= (n− 1 ) . The position Xi,k of crow i at iteration k is described in Equation (14) and represents a feasible problem solution. Xi,k= [x1 i,k,x2 i,k,...., xn i,k], (14) where i= 1,2,..., N and k= 1,2,..., tmax . Further, tmax is the maximum number of iterations of the exploration and exploitation process of the solution space. Each crow (individual) can memorize the position of its hiding place. At iteration k , the position of crow i hiding place is represented as Mi,k , being the best position that crow i has obtained so far. Of course, the position of its best experience is memorized because crows move in their environment and search for the best food source (i.e., hiding places). Symmetry 2021,13, 1329 8 of 20 At the start of the optimization process, it is assumed that at iteration k , crow j wants to visit its hideout, Mj,k . In this iteration, crow i decided to follow crow j to approach crow j’s hideout. In this case, two situations can occur. Situation 1: Search Crow j does not know that crow i is following it. As a result, crow i approaches crow j ’s hiding place. In this case, through Equation (15), the new position of crow i is obtained. Xi,k+1=Xi,k+rif li,kMj,k−Xi,k, (15) where ri a random number between 0 and 1, and f li,k is the flight length of the crow i at each iteration k . According to [ 29 ], small values of f l allow a local exploration of the solution space (close to Xi,k ), while large values of f l allow a global solution space exploration (far from Xi,k). Situation 2: Evasion Crow j knows that crow i is following it. As a result, crow j tries to deceive crow i by heading towards another position in the solution space to protect its caches from being stolen. Either way, situations 1 and 2 can be represented as shown in Equation (16). Xi,k+1=Xi,k+rif li,k(Mj,k−Xi,k)if rj≥AP random otherwise, (16) where rj is a random number between 0 and 1, and AP denotes the probability of crow j ’s knowledge of crow i. Once the crows’ positions are modified, the memory of each crow is updated based on the objective function values of the new spots. So, if the objective function of the new location is better than the objective function of the memorized position, the crows update their memory to the new area, as shown in Equation (17). Mi,k+1=F(Xi,k+1)if F(Xi,k+1)<F(Mi,k) Mi,kotherwise , (17) where F(·)represents the minimized objective function. 3.2.2. Proposed ICSA The CSA has confirmed its capability to reach the optimal solution for particular solution space configurations [ 29 , 39 , 40 ]. Nevertheless, its convergence is not ensured due to the inefficient exploration of its search strategy [ 37 ]. Therefore, it presents difficulties when facing high-dimensional problems [ 37 , 41 ]. In the original CSA method [ 29 ], there are two elements responsible for the search process: knowledge probability (i.e., AP ) and random motion (i.e., Situation 2: Avoidance) [41]. The value of AP is entrusted with providing an adequate equilibrium between diversification and intensification [ 29 ]. With small AP values, a local solution space search is obtained, increasing intensification [ 29 ]. On the other hand, with large AP values, a global solution space search is obtained, which increases diversification [ 29 ]. Since metaheuristic algorithms require a balance between diversification and intensification to find a globally optimal solution when solving problems with large dimensions [ 42 , 43 ], AP is taken as an intermediate value, i.e., 50%. Moreover, the random motion specifically impacts the CSA search mechanism since it resets the candidate solutions, deviating them from the current best solution and delaying the convergence of the problem [ 37 ]. In the proposed ICSA, the random motion is reformulated, as shown below. Symmetry 2021,13, 1329 9 of 20 3.2.3. Improved Approach: VSA for Random Movement The evasion behavior is simulated by a random motion implementation computed through a uniformly distributed random value [ 29 ]. In the proposed ICSA, to have a better solution space diversification, the possibility of generated candidate solutions based on the VSA evolution criteria is added to the classical CSA [ 15 ], considering that the main feature of the VSA is the use of a regular GD to generate neighbors around the current best solution named as the center of the hyper-sphere µ [ 44 ]. In this paper, the hyper-sphere center is selected as the best solution in the current population during iteration kas µk=Xbest,k. The set of candidate solutions Cv,k(y) is generated using a random GD in the space d around the best solution Xbest,k, as shown in Equation (18). Cv,k(y) = p(y|µk,Σ) = 1 p(2π)d|Σ|)exp−1 2(y−µk)TΣ−1(y−µk), (18) where d is the solution space dimension, y∈Rd×1 corresponds to a vector of random variables with values between zero and one, and Σ∈Rd×d is the matrix of covariance. If, in Σ , the on-diagonal elements (variance) are equally defined and if the off-diagonal components (covariance) are chosen as zero, then the GD will generate hyper-spheres in the d− dimensional space [ 35 ]. Equation (19) displays a simple way to calculate Σ , considering equal variance and zero covariances. Σ=σ2Id×d, (19) where Id×d is an identity matrix and σ represents the variance of the GD. Note that the standard deviation of the GD can be defined as shown in Equation (20). σ0=max{ymax}−minymin 2, (20) where ymax and ymin are vectors of dimension d× 1 that define the upper and lower bounds of the decision variables of the optimization problem, respectively. Here, σ0 can also be considered as the initial radius r0 of the hyper-sphere [ 35 ]. To achieve a proper exploration of the solution space, initially, σ0is the largest possible hyper-sphere. The candidate solutions obtained and contained in the set Cv,k(y) must guarantee that the results lie within the bounds of the solution space; so, Equation (21) is employed. Cv,k(y) = Cv,k(y)ymin ≤y≤ymax ymin + (ymax −ymin)rand otherwise , (21) where rand generates random numbers between 0 and 1. Once verified the limits, the best solution obtained in the set Cv,k(y)is selected to be the new position of the crow i. It is necessary to mention that the radius of the hyper-sphere decreases as the iteration process progresses using an inverse incomplete gamma function, as reported in [ 45 ] and as shown in Equation (22). rk=σ0γ−1(y,ak)(22) The inverse incomplete gamma function for the variable radius calculation can be calculated in MATLAB®, as shown in Equation (23) [35]. rk=σ01 ygammaincinv(y,ak), (23) where ak is a parameter defined as ak= 1 −k tmax . Moreover, the parameter y is chosen as 0.1, as recommended in [35]. Algorithm 2summarizes the implementation of the CSA, considering the VSA evolution criterion for Situation 2, to solve the phase swapping problem in ADSs. Symmetry 2021,13, 1329 16 of 20 point for future approaches that can be proposed to solve the phase swapping problem in ADSs. CBGA VSA CSA ICSA 19.000 19.050 19.100 19.150 19.200 19.250 19.120 19.249 19.018 19.252 Reduction (%) Figure 8. Reduction of the power losses in the 37-bus test feeder. 5.4. Complementary Results The following can be concluded from the results obtained in Section 5. X The optimal solution achieved with the enhanced version of the CSA for each test system is equal to the one reported in [ 15 ] for the 8- and 25-node systems. However, it obtained a better solution for the 37-node system. Note that the obtained solutions (power losses) in this research are 10.5869 kW, 72.2888 kW, and 61.4781 kW for the 8-, 25-, and 37-node systems, respectively, while the best solutions reported in [ 15 ] are 10.5869 kW, 72.2888 kW, and 61.4801 kW, respectively. It is relevant to highlight that the ICSA requires a longer processing time. Nevertheless, these times do not exceed 6 minutes, which is not significant considering the terms for optimization problems with solution spaces with hundreds of thousands of combinations. It ensures excellent quality solutions and even better ones than the values reported in the literature review in some cases (see results for the 37-node system). However, simulations in tests feeders with large number of nodes will be required to ensure that, in all of the cases, the processing times spent by the proposed ICSA will be compensated with optimal solutions better than the solutions provided by the VSA. X The standard deviations reported in Tables 3–5for the 8-, 25-, and 37-node systems, respectively, for ICSA are lower than those reported in [ 15 ]. In addition, a standard deviation of 0.1344 kW for the 37-node test system demonstrates the repeatability properties of the ICSA for solving the phasing problem, considering that the dimensions of the solution space are higher than 1 × 10 28 . Regarding metaheuristic optimization methods, the main preoccupation in the literature is associated with the ability of these methods to find the same optimal solution at each simulation. However, this is not possible due to the random nature of the exploration and exploitation criteria inside of each metaheuristic optimizer. Nevertheless, when an optimization methodology exhibits low standard deviations, all the solutions are contained inside of a small hyper-sphere close to the global optimum. This improves the optimization method of the proposed ICSA when compared with a family of metaheuristic optimizers. X When comparing the base cases of each test system with the proposed methodologies, as shown in Table 3and Figures 5and 7, the reductions in energy losses resulting from ICSA in the 8, 25, and 37-node systems are 24.34%, 4.152%, and 19.252%, respectively. The results show that better results are obtained when applying the proposed methodology for the 37-node system, in contrast to the CSA and the VSA [ 15 ]. Likewise, it is observed that the proposed enhancement for the CSA effectively explores and exploits the solution space for large systems higher than 25 nodes in the case of distribution power systems. Symmetry 2021,13, 1329 17 of 20 X The comparison of the effect in the energy losses redistribution at each phase using the classic and improved CSA methods is presented in Table 6. Table 6. Comparison between the CSA and the proposed ICSA regarding power losses among phases. Method Phase a(kW) Phase b(kW) Phase b(kW) Total (kW) 8-bus system Benchmark case 1.7158 2.3305 9.9462 13.9925 CSA 3.7617 2.7295 4.0957 10.5869 ICSA 2.7412 3.9930 3.8464 10.5869 25-bus system Benchmark case 36.8801 14.7837 23.7570 75.4207 CSA 24.9590 26.2079 21.1627 72.3296 ICSA 25.6645 26.1613 20.4630 72.2888 IEEE 37-bus system Benchmark case 27.1532 11.9143 37.0683 76.1357 CSA 19.1177 20.5926 21.1627 61.6565 ICSA 20.3761 21.8327 19.2694 61.4781 With these results, we can make several observations. (i) In the 8-bus system, the phase c presents power losses higher than 4 kW, while the ICSA does not support this value; however, both solutions are indeed optimal since the total power losses is the same for both methods. (ii) In the 25-bus system both methods, i.e., the CSA and its improved version, the phase c has been identified to present higher power losses surpassing 26 kW; however, regarding the final power losses, the ICSA presents better load redistribution, since the amount of power losses is about 0.0408 kW. Finally, (iii) in the 37-bus system, the CSA method presents a difference between phases a and b power losses of 2.0450 kW, while the proposed ICSA has a minor difference between both phases with a value of 1.1067 kW. This directly impacts the final grid results with a general improvement of about 0.1784 kW in favor of the proposed ICSA, which is indeed the best global optimum reported in the current literature for the IEEE 37-bus system. 6. Conclusions and Future Works This paper presented a master–slave methodology to solve the phase equilibrium problem in ADSs. For this purpose, an ICSA was proposed using the VSA evolution mechanism. The master stage determined the set of phase configurations for the system three-phase charges through the ICSA, using an integer encoding between 1 and 6 representing the load connections in the three-phase system. On the other hand, the slave stage evaluated each of the load connections in the TPPF, which correspond to the extended version of the iterative sweep power flow method for unbalanced three-phase systems. The numerical results on the 8-, 25-, and 37-node test systems showed that the proposed ICSA, compared to the VSA, achieves identical solutions for the 8-node and 25-node systems. However, for the 37-node system, the ISCA obtains a better optimal solution when compared with the current report employing the VSA method. The solutions obtained by the ICSA are 10.5869 kW, 72.2888 kW, and 61.4781 kW for the 8-, 25-, and 37-node systems, representing a reduction of 24.34%, 4.152%, and 19.252%, respectively. In the same way, the solutions for the VSA are 10.5869 kW, 72.2888 kW, and 61.4801 kW, respectively. In addition, the proposed methodology has minor standard deviations for solving the phase swapping problem for the 8-, 25-, and IEEE 37-node test systems, which were 1.099 ×10−5kW , 0.0116 kW, and 0.1344 kW, respectively. This demonstrates better repeatability properties of the improved algorithm, since all the solutions are contained inside small hyper-spheres around the global optimum. Further, the solution space for the phase equilibrium problem potentially increases as a function of the demand nodes. Thus, for the Symmetry 2021,13, 1329 18 of 20 IEEE 37-node system, there is a solution space bigger than 1 × 10 28 , which is by far higher than 100 billion combinations. This result confirmed the effectiveness of the proposed methodology for solving complex MINLP models such as the phase equilibrium problem in ADSs as well as in being better than the optimal solution reported in the current literature using the VSA. In future work, it is possible to accomplish the following: (i) combine the proposed phase swapping with the optimal placement of distributed generators to reduce power losses in ADSs; (ii) employ typical active and reactive power demand curves to solve the phase swapping problem using the ICSA to reduce power losses; and (iii) use the ICSA to solve the optimal reconfiguration problem in three-phase radial distribution networks. Author Contributions: Conceptualization, methodology, software, and writing—review and editing, B.C.-C., L.S.A.-G., O.D.M., L.A.-B. and C.Á.-A. All authors have read and agreed to the published version of the manuscript. Funding: This work was partially supported in part by the Laboratorio de Simulación Hardware-in- the-loop para Sistemas Ciberfísicos de la Universidad Loyola Andalucía. Institutional Review Board Statement: Not applicable. Informed Consent Statement: Not applicable. Data Availability Statement: No new data were created or analyzed in this study. Data sharing is not applicable to this article. Acknowledgments: This work was supported in part by the Centro de Investigación y Desarrollo Científico de la Universidad Distrital Francisco José de Caldas under grant 1643-12-2020 associated with the project: “Desarrollo de una metodología de optimización para la gestión óptima de recursos energéticos distribuidos en redes de distribución de energía eléctrica.” and in part by the Dirección de Investigaciones de la Universidad Tecnológica de Bolívar under grant PS2020002 associated with the project: “Ubicación óptima de bancos de capacitores de paso fijo en redes eléctricas de distribución para reducción de costos y pérdidas de energía: Aplicación de métodos exactos y metaheurísticos.” Conflicts of Interest: The authors declare no conflict of interest. References 1. Montoya, O.D.; Serra, F.M.; De Angelo, C.H. On the efficiency in electrical networks with ac and dc operation technologies: A comparative study at the distribution stage. Electronics 2020,9, 1352. [CrossRef] 2. Montoya, O.D.; Gil-González, W.; Hernández, J.C. Efficient Operative Cost Reduction in Distribution Grids Considering the Optimal Placement and Sizing of D-STATCOMs Using a Discrete-Continuous VSA. Appl. Sci. 2021,11, 2175. [CrossRef] 3. Nassar, M.E.; Salama, M. A novel branch-based power flow algorithm for islanded AC microgrids. Electr. Power Syst. Res. 2017 , 146, 51–62. [CrossRef] 4. Aboshady, F.; Thomas, D.W.; Sumner, M. A wideband single end fault location scheme for active untransposed distribution systems. IEEE Trans. Smart Grid 2019,11, 2115–2124. [CrossRef] 5. Montoya, O.D.; Arias-Londoño, A.; Grisales-Noreña, L.F.; Barrios, J.Á.; Chamorro, H.R. Optimal Demand Reconfiguration in Three-Phase Distribution Grids Using an MI-Convex Model. Symmetry 2021,13, 1124. [CrossRef] 6. Arias, J.; Calle, M.; Turizo, D.; Guerrero, J.; Candelo-Becerra, J.E. Historical load balance in distribution systems using the branch and bound algorithm. Energies 2019,12, 1219. [CrossRef] 7. Alvarado-Barrios, L.; Álvarez-Arroyo, C.; Escaño, J.M.; Gonzalez-Longatt, F.M.; Martinez-Ramos, J.L. Two-Level Optimisation and Control Strategy for Unbalanced Active Distribution Systems Management. IEEE Access 2020 ,8, 197992–198009. [CrossRef] 8. Granada Echeverri, M.; Gallego Rendón, R.A.; López Lezama, J.M. Optimal phase balancing planning for loss reduction in distribution systems using a specialized genetic algorithm. Ing. Cienc. 2012,8, 121–140. [CrossRef] 9. Hooshmand, R.; Soltani, S. Simultaneous optimization of phase balancing and reconfiguration in distribution networks using BF–NM algorithm. Int. J. Electr. Power Energy Syst. 2012,41, 76–86. [CrossRef] 10. Ogunsina, A.A.; Petinrin, M.O.; Petinrin, O.O.; Offornedo, E.N.; Petinrin, J.O.; Asaolu, G.O. Optimal distributed generation location and sizing for loss minimization and voltage profile optimization using ant colony algorithm. SN Appl. Sci. 2021 ,3, 1–10. [CrossRef] 11. Asadi, M.; Shokouhandeh, H.; Rahmani, F.; Hamzehnia, S.M.; Harikandeh, M.N.; Lamouki, H.G.; Asghari, F. Optimal placement and sizing of capacitor banks in harmonic polluted distribution network. In Proceedings of the 2021 IEEE Texas Power and Energy Conference (TPEC), College Station, TX, USA, 2–5 February 2021; pp. 1–6. Symmetry 2021,13, 1329 19 of 20 12. Jakus, D.; ˇ Ca ¯ denovi´c, R.; Vasilj, J.; Sarajˇcev, P. Optimal reconfiguration of distribution networks using hybrid heuristic-genetic algorithm. Energies 2020,13, 1544. [CrossRef] 13. Montoya, O.D.; Garces, A.; Castro, C.A. Optimal conductor size selection in radial distribution networks using a mixed-integer non-linear programming formulation. IEEE Lat. Am. Trans. 2018,16, 2213–2220. [CrossRef] 14. Łukaszewski, A.; Nogal, Ł.; Robak, S. Weight Calculation Alternative Methods in Prime’s Algorithm Dedicated for Power System Restoration Strategies. Energies 2020,13, 6063. [CrossRef] 15. Cortés-Caicedo, B.; Avellaneda-Gómez, L.S.; Montoya, O.D.; Alvarado-Barrios, L.; Chamorro, H.R. Application of the Vortex Search Algorithm to the Phase-Balancing Problem in Distribution Systems. Energies 2021,14, 1282. [CrossRef] 16. Montoya, O.D.; Molina-Cabrera, A.; Grisales-Noreña, L.F.; Hincapié, R.A.; Granada, M. Improved Genetic Algorithm for Phase-Balancing in Three-Phase Distribution Networks: A Master-Slave Optimization Approach. Computation 2021 ,9, 67. [CrossRef] 17. Saad Al-Sumaiti, A.; Kavousi-Fard, A.; Salama, M.; Pourbehzadi, M.; Reddy, S.; Rasheed, M.B. Economic Assessment of Distributed Generation Technologies: A Feasibility Study and Comparison with the Literature. Energies 2020 ,13, 2764. [CrossRef] 18. Rajaram, R.; Kumar, K.S.; Rajasekar, N. Power system reconfiguration in a radial distribution network for reducing losses and to improve voltage profile using modified plant growth simulation algorithm with Distributed Generation (DG). Energy Rep. 2015 , 1, 116–122. [CrossRef] 19. Ahshan, R. Analysis of Loss Reduction Techniques for Low Voltage Distribution Network. J. Eng. Res. [TJER] 2020 ,17, 100–111. 20. Grigora s , , G.; Neagu, B.C.; Gavrila s , , M.; Tri s , tiu, I.; Bulac, C. Optimal phase load balancing in low voltage distribution networks using a smart meter data-based algorithm. Mathematics 2020,8, 549. [CrossRef] 21. Montoya, O.D.; Giraldo, J.S.; Grisales-Noreña, L.F.; Chamorro, H.R.; Alvarado-Barrios, L. Accurate and Efficient Derivative-Free Three-Phase Power Flow Method for Unbalanced Distribution Networks. Computation 2021,9, 61. [CrossRef] 22. Taghipour Boroujeni, S.; Mardaneh, M.; Hashemi, Z. A dynamic and heuristic phase balancing method for LV feeders. Appl. Comput. Intell. Soft Comput. 2016,2016, doi:10.1155/2016/6928080. [CrossRef] 23. Gandomkar, M. Phase balancing using genetic algorithm. In Proceedings of the 39th International Universities Power Engineering Conference, 2004, UPEC 2004, Bristol, UK, 6–8 September 2004; Volume 1, pp. 377–379. 24. Rios, M.A.; Castaño, J.C.; Garcés, A.; Molina-Cabrera, A. Phase Balancing in Power Distribution Systems: A heuristic approach based on group-theory. In Proceedings of the 2019 IEEE Milan PowerTech, Milan, Italy, 23–27 June 2019; pp. 1–6. 25. Garces, A.; Gil-González, W.; Montoya, O.D.; Chamorro, H.R.; Alvarado-Barrios, L. A Mixed-Integer Quadratic Formulation of the Phase-Balancing Problem in Residential Microgrids. Appl. Sci. 2021,11, 1972. [CrossRef] 26. Amon, D.A.; Adeyemi, A. A modified bat algorithm for power loss reduction in electrical distribution system. Indones. J. Electr. Eng. Comput. Sci. (IJEECS) 2015,14, 55–61. 27. Sathiskumar, M.; kumar, A.N.; Lakshminarasimman, L.; Thiruvenkadam, S. A self adaptive hybrid differential evolution algorithm for phase balancing of unbalanced distribution system. Int. J. Electr. Power Energy Syst. 2012,42, 91–97. [CrossRef] 28. Zhu, J.; Bilbro, G.; Chow, M.Y. Phase balancing using simulated annealing. IEEE Trans. Power Syst. 1999 ,14, 1508–1513. [CrossRef] 29. Askarzadeh, A. A novel metaheuristic method for solving constrained engineering optimization problems: Crow search algorithm. Comput. Struct. 2016,169, 1–12. [CrossRef] 30. Dogan, B. A Modified Vortex Search Algorithm for Numerical Function Optimization. Int. J. Artif. Intell. Appl. 2016 ,7, 37–54. [CrossRef] 31. Shen, T.; Li, Y.; Xiang, J. A graph-based power flow method for balanced distribution systems. Energies 2018 ,11, 511. [CrossRef] 32. Herrera-Briñez, M.C.; Montoya, O.D.; Alvarado-Barrios, L.; Chamorro, H.R. The Equivalence between Successive Approximations and Matricial Load Flow Formulations. Appl. Sci. 2021,11, 2905. [CrossRef] 33. Łukaszewski, A.; Nogal, Ł. Influence of lightning current surge shape and peak value on grounding parameters. Bull. Pol. Acad. Sci. Tech. Sci. 2021,69, e136730. [CrossRef] 34. Gil-González, W.; Montoya, O.D.; Rajagopalan, A.; Grisales-Noreña, L.F.; Hernández, J.C. Optimal Selection and Location of Fixed-Step Capacitor Banks in Distribution Networks Using a Discrete Version of the Vortex Search Algorithm. Energies 2020 , 13, 4914. [CrossRef] 35. Do˘gan, B.; Ölmez, T. A new metaheuristic for numerical function optimization: Vortex Search algorithm. Inf. Sci. 2015 , 293, 125–145. [CrossRef] 36. Hussien, A.G.; Amin, M.; Wang, M.; Liang, G.; Alsanad, A.; Gumaei, A.; Chen, H. Crow Search Algorithm: Theory, Recent Advances, and Applications. IEEE Access 2020,8, 173548–173565. [CrossRef] 37. Jain, M.; Rani, A.; Singh, V. An improved Crow Search Algorithm for high-dimensional problems. J. Intell. Fuzzy Syst. 2017 , 33, 3597–3614. [CrossRef] 38. Clayton, N.; Emery, N. Corvid cognition. Curr. Biol. 2005,15, R80–R81. [CrossRef] 39. Askarzadeh, A. Capacitor placement in distribution systems for power loss reduction and voltage improvement: A new methodology. IET Gener. Transm. Distrib. 2016,10, 3631–3638. [CrossRef] 40. Rajput, S.; Parashar, M.; Dubey, H.M.; Pandit, M. Optimization of benchmark functions and practical problems using Crow Search Algorithm. In Proceedings of the 2016 Fifth International Conference on Eco-friendly Computing and Communication Systems (ICECCS), Bhopal, India, 8–9 December 2016; pp. 73–78. Symmetry 2021,13, 1329 20 of 20 41. Díaz, P.; Pérez-Cisneros, M.; Cuevas, E.; Avalos, O.; Gálvez, J.; Hinojosa, S.; Zaldivar, D. An improved crow search algorithm applied to energy problems. Energies 2018,11, 571. [CrossRef] 42. Li, Z.; Wang, W.; Yan, Y.; Li, Z. PS–ABC: A hybrid algorithm based on particle swarm and artificial bee colony for high-dimensional optimization problems. Expert Syst. Appl. 2015,42, 8881–8895. [CrossRef] 43. Sayed, G.I.; Hassanien, A.E.; Azar, A.T. Feature selection via a novel chaotic crow search algorithm. Neural Comput. Appl. 2019 , 31, 171–188. [CrossRef] 44. Li, P.; Zhao, Y. A quantum-inspired vortex search algorithm with application to function optimization. Nat. Comput. 2019 , 18, 647–674. [CrossRef] 45. Montoya, O.D.; Grisales-Noreña, L.F.; Amin, W.T.; Rojas, L.A.; Campillo, J. Vortex Search Algorithm for Optimal Sizing of Distributed Generators in AC Distribution Networks with Radial Topology. In Workshop on Engineering Applications; Springer: Berlin/Heidelberg, Germany, 2019; pp. 235–249. [CrossRef]