# Improved Genetic Algorithm for Phase-Balancing in Three-Phase Distribution Networks: A Master-Slave Optimization Approach

^{1}

^{2}

^{3}

^{4}

^{*}

Next Article in Journal

Next Article in Special Issue

Next Article in Special Issue

Previous Article in Journal

Previous Article in Special Issue

Previous Article in Special Issue

Facultad de Ingeniería, Universidad Distrital Francisco José de Caldas, Bogotá D.C. 11021, Colombia

Laboratorio Inteligente de Energía, Universidad Tecnológica de Bolívar, Cartagena 131001, Colombia

Facultad de Ingeniería, Universidad Tecnológica de Pereira, Pereira 660003, Colombia

Grupo GIIEN, Facultad de Ingeniería, Institución Universitaria Pascual Bravo, Medellín 050036, Colombia

Author to whom correspondence should be addressed.

Academic Editors: George Tsakalidis and Kostas Vergidis

Received: 14 May 2021 / Revised: 2 June 2021 / Accepted: 5 June 2021 / Published: 9 June 2021

(This article belongs to the Special Issue Recent Advances in Process Modeling and Optimisation)

This paper addresses the phase-balancing problem in three-phase power grids with the radial configuration from the perspective of master–slave optimization. The master stage corresponds to an improved version of the Chu and Beasley genetic algorithm, which is based on the multi-point mutation operator and the generation of solutions using a Gaussian normal distribution based on the exploration and exploitation schemes of the vortex search algorithm. The master stage is entrusted with determining the configuration of the phases by using an integer codification. In the slave stage, a power flow for imbalanced distribution grids based on the three-phase version of the successive approximation method was used to determine the costs of daily energy losses. The objective of the optimization model is to minimize the annual operative costs of the network by considering the daily active and reactive power curves. Numerical results from a modified version of the IEEE 37-node test feeder demonstrate that it is possible to reduce the annual operative costs of the network by approximately 20% by using optimal load balancing. In addition, numerical results demonstrated that the improved version of the CBGA is at least three times faster than the classical CBGA, this was obtained in the peak load case for a test feeder composed of 15 nodes; also, the improved version of the CBGA was nineteen times faster than the vortex search algorithm. Other comparisons with the sine–cosine algorithm and the black hole optimizer confirmed the efficiency of the proposed optimization method regarding running time and objective function values.

Three-phase distribution grids are responsible for providing electrical energy to residential, commercial, and industrial areas at medium- and low-voltage levels [1,2]. These grids are typically unbalanced owing to the following factors: (i) line configurations are non-symmetric because the transposition criteria are not applicable owing to the short length of the lines (few tens of kilometers) [3]; and (ii) the nature of loads can be single-, two-, and three-phase, which results in intrinsic imbalances in the currents through the lines and voltages in the nodes [4]. The main problem arising from the general unbalanced nature of three-phase networks is the large power losses in the distribution stage (distribution lines and transformers); these losses cause increments in the final billing of the electricity service for all the end-users, which is permitted by regulatory entities within an acceptable range [2]. In Colombia, the power losses of distribution systems range from 6% to 18%, depending on the geographic location of the distribution systems (i.e., rural or urban areas) as well as the level of investment in maintenance activities on the physical layer of the distribution network [5]. The regulatory entity of the electrical distribution service in Colombia, that is, the Energy and Gas Regulation Commission (CREG, based on its acronym in Spanish) allows transferring part of the energy losses to the end-users through the final billing for the electricity service; however, this value cannot exceed 8% of the total energy costs. This implies that, if the distribution company reduces the amount of energy losses to less than 8%, the following benefits can be achieved: (i) an improvement in the potential of the network to supply additional users with the same distribution grid and minimum investment costs (expansion in coverage capacity); (ii) additional economic reedits because lower energy input is needed for supplying the same group of users; and (iii) the possibility to win by billing the difference between the upper regulation bound (i.e., 8%) and the current power losses [6].

Owing to the importance of strategies for reducing the amount of energy loss in electrical distribution networks, multiple well-known approaches have been established. Some of these include (i) optimal sizing and installation of distributed generation [7]; (ii) optimal installation and sizing of reactive power compensators [5]; (iii) optimal grid reconfiguration [8]; and (iv) optimal phase-balancing strategies [2]. Each of the aforementioned methodologies can help distribution companies in significantly reducing the amount of power loss. However, the first two strategies entail large amounts of inversions due to the necessity of integrating new devices into the distribution grid, which implies that the return rate is in the order of years (5 to 15 years) [9]. The third approach, which is based on reconfiguration, requires less investment because few tie lines need to be constructed for the optimal reconfiguration of the grid [10]. The fourth approach is the most economical strategy for reducing power and energy losses in distribution grids, because it only requires a few work crews to reconfigure the load phases without investing in additional devices [11]. Here, considering the low costs associated with the phase-balancing strategy for reducing the amount of energy loss in electrical distribution networks, we propose a new master–slave optimization strategy to address this problem.

As indicated by existing literature, the problem of phase balancing has been solved using multiple optimization approaches, including the classical Chu and Beasley genetic algorithms [12,13,14,15]; tabu search algorithm [16,17]; ant colony optimization [18,19]; immune optimization algorithm [20]; branch and bound and convex methods [4,21]; bat optimization algorithm [22]; vortex search algorithm [2]; particle swarm optimization methods [8,23,24]; differential evolution algorithm [25]; and simulated annealing optimization method [26]. The main characteristics of these approaches are the hybridization of metaheuristic discrete optimization methods (with binary and integer codifications) with three-phase power flow methods, which are typically based on sweep iterative backward/forward power flow methods, and the minimization of the amount of power loss for a particular load condition, which typically corresponds to the peak load condition.

Based on the aforementioned state-of-the-art strategies, in this paper, we propose a master–slave optimization approach based on an improved version of the classical Chu and Beasley genetic algorithm (CBGA) in the master stage, in conjunction with using the three-phase version of the successive approximation power flow method in the slave stage. The master stage is entrusted by defining the phase load configuration using an integer codification between 1 and 6, which represents the six possible methods for connecting a three-phase load. The slave stage is entrusted with the evaluation of the multiperiod power flow problem, in order to determine the amount of power loss in each nodal phase configuration. Improvements in the classical CBGA are realized in the selection, mutation, and recombination procedures based on a probability criterion. If the probability is less than 50%, two individuals are selected in the tournament stage, which are then recombined using a single point. Subsequently, an adaptive multi-point mutation strategy is employed to increase the rate of convergence of the algorithm; however, if the probability exceeds 50%, an offspring is generated using a normal Gaussian distribution probability, which is employed in the evolution process of the vortex search algorithm approach [27]. The main advantage of the vortex search algorithm is the exploration and exploitation of a solution space by using adaptive hyper-ellipses based on one individual solution selected from the current population. It should be noted that this evolution criterion of the proposed CBGA is applied at each iteration; this implies that the evolution process of the algorithm is adaptive.

The main contributions of this research are as follows: (i) we propose an improved version of the CBGA based on adaptive evolution criteria by combining the classical selection, recombination, and mutation stages with the hyper-ellipses used in the vortex search algorithm, in order to explore and exploit the solution space; (ii) the solution of the phase-balancing problem in distribution grids is determined considering daily and reactive power load curves; and (iii) a three-phase version of the successive approximation power flow method, which is suitable for loads connected in Y and $\Delta $, is formulated.

It is important to highlight that the proposed improved version of the classical CBGA with the vortex search algorithm (VSA) increases the possibility of enhancing the exploration and exploitation of the solution space. This is made by introducing normal Gaussian distributions to sweep some promissory solution zones. Hence, this new version of the CBGA with the VSA method can be considered as a powered metaheuristic optimization technique that can be extended to multiple MINLP problems in future researches [27]. Furthermore, an additional important contribution presented in this study corresponds to the extension of the recently developed successive approximation power flow method [28] to three-phase unbalanced distribution networks. The resulting method could be applied to loads with $\Delta -$ and $Y-$connections including possible meshed grid configurations with the possibility of demonstrating its convergence by applying the Banach fixed-point theorem as recently proposed in [29] for the upper-triangular power flow method.

The remainder of the paper is structured as follows: Section 2 presents the mathematical formulation of the phase-balancing problem in three-phase distribution grids by highlighting the mixed-integer nonlinear structure of the optimization problem. Section 3 presents the master–slave optimization approach based on the improved version of the CBGA and the three-phase version of the successive approximation power flow method. Section 4 describes the main characteristics of the test feeders comprising 15 and 37 nodes. Section 5 presents the numerical validation of the proposed optimization approach for the IEEE 37-node system considering the daily load curves. The peak load scenario operative scenario is also compared with the classical CBGA, the black hole optimizer, the sine–cosine algorithm and the vortex search algorithm. Section 6 presents the main concluding remarks derived from this research as well as some possible future works.

The phase-balancing problem in three-phase imbalanced distribution networks can be formulated as a mixed-integer nonlinear programming (MINLP) model, in which the integer variables are related to a combination of the load phases. Moreover, the nonlinearities are produced by the power balance equations, which yield products among voltage magnitudes and trigonometric functions. Here, we present the objective function and a set of constraints that represent the phase-balancing problem, considering multiperiod analyses.

The objective function considered for the phase-balancing problem corresponds to the minimization of costs for the daily energy loss at all the conductors in the grid. This objective function can be expressed as follows:
where ${f}_{1}$ represents the annual operative costs associated with the total energy losses in all the branches of the network; ${C}_{\mathrm{kWh}}$ is the average cost of energy; ${Y}_{kfmg}$ represents the magnitude of the admittance that relates node k at phase f with node m at phase g; ${V}_{kft}$ (${V}_{mft}$) corresponds to the voltage magnitude at node k (m) in phase f (g) at time period t; ${\delta}_{kft}$ (${\delta}_{mgt}$) is the angle of the voltage at node k (m) in phase f (g) at time period t; ${\theta}_{kfmg}$ represents the angle of the admittance that relates node k at phase f with node m at phase g; and ${\Delta}_{t}$ is the time period during which the power demands remain constant. It should be noted that $\mathcal{F}$, $\mathcal{N}$, and $\mathcal{T}$ are the sets that contain all the phases, nodes, and time periods, respectively.

$$\begin{array}{c}\hfill {f}_{1}={C}_{\mathrm{kWh}}\sum _{t\in \mathcal{T}}\sum _{k\in \mathcal{N}}\sum _{m\in \mathcal{N}}\sum _{f\in \mathcal{F}}\sum _{g\in \mathcal{F}}{Y}_{kfmg}{V}_{kft}{V}_{mgt}cos\left({\delta}_{kft}-{\delta}_{mgt}-{\theta}_{kfmg}\right){\Delta}_{t},\end{array}$$

Notably, the objective function is nonlinear and non-convex owing to the product among voltages and trigonometric functions; this implies that its minimization requires advanced optimization strategies that can efficiently deal with nonlinearities [30]. This justifies the use of the master–slave optimization strategy, as in the case of the improved version of the CBGA proposed herein, owing to its simplicity in terms of programming and its high efficiency in relation to the processing times.

The set of constraints for the phase-balancing problem is typically associated with the operative conditions of the network, that is, the power balance equations and the voltage regulation bounds [2]. These equations are as follows:
where ${P}_{kft}^{s}$ represents the active power generated by source s connected at node k in phase f at time period t; ${Q}_{kft}^{s}$ represents the reactive power generated by source s connected at node k in phase f at time period t; ${P}_{kgt}^{d}$ represents the active power demand at node k in phase g at time period t; ${Q}_{kgt}^{d}$ is the reactive power demand at node k in phase g at time period t; ${x}_{kfg}$ is a binary variable that defines the connection of the demand at node k at f in phase g; and ${V}_{min}$ and ${V}_{max}$ correspond to the minimum and maximum voltage regulation bounds allowed for all nodes of the network at each period of time, respectively.

$$\begin{array}{c}\hfill {P}_{kft}^{s}-\sum _{g\in \mathcal{F}}{x}_{kfg}{P}_{kgt}^{d}={V}_{kft}\sum _{m\in \mathcal{N}}\sum _{g\in \mathcal{F}}{Y}_{kfmg}{V}_{mgt}cos\left({\delta}_{kft}-{\delta}_{mgt}-{\theta}_{kfmg}\right),\phantom{\rule{0.277778em}{0ex}}\left\{\begin{array}{c}\forall f\in \mathcal{F}\\ \forall k\in \mathcal{N}\\ \forall t\in \mathcal{T}\end{array}\right\}\end{array}$$

$$\begin{array}{c}\hfill {Q}_{kft}^{s}-\sum _{g\in \mathcal{F}}{x}_{kfg}{Q}_{kgt}^{d}={V}_{kft}\sum _{m\in \mathcal{N}}\sum _{g\in \mathcal{F}}{Y}_{kfmg}{V}_{mgt}sin\left({\delta}_{kft}-{\delta}_{mgt}-{\theta}_{kfmg}\right),\phantom{\rule{0.277778em}{0ex}}\left\{\begin{array}{c}\forall f\in \mathcal{F}\\ \forall k\in \mathcal{N}\\ \forall t\in \mathcal{T}\end{array}\right\}\end{array}$$

$$\begin{array}{c}\hfill \sum _{g\in \mathcal{F}}{x}_{kfg}=1,\phantom{\rule{0.277778em}{0ex}}\left\{\forall f\in \mathcal{F},\phantom{\rule{0.277778em}{0ex}}\forall k\in \mathcal{N},\phantom{\rule{0.277778em}{0ex}}\forall t\in \mathcal{T}\right\}\end{array}$$

$$\begin{array}{c}\hfill \sum _{f\in \mathcal{F}}{x}_{kfg}=1,\phantom{\rule{0.277778em}{0ex}}\left\{\forall g\in \mathcal{F},\phantom{\rule{0.277778em}{0ex}}\forall k\in \mathcal{N},\phantom{\rule{0.277778em}{0ex}}\forall t\in \mathcal{T}\right\}\end{array}$$

$$\begin{array}{c}\hfill {V}_{min}\le {V}_{kft}\le {V}_{max},\phantom{\rule{0.277778em}{0ex}}\left\{\forall g\in \mathcal{F},\phantom{\rule{0.277778em}{0ex}}\forall k\in \mathcal{N},\phantom{\rule{0.277778em}{0ex}}\forall t\in \mathcal{T}\right\}\end{array}$$

Expressions (2) and (3), associated with the active and reactive power balance equations at each node, phase, and period of time, evidence the complication of the three-phase power flow problem in electrical networks, even when the connection of the loads at these nodes are well-known. This is due to the nonlinear, non-affine equations that must be solved numerically.

The following interpretations can be drawn from the optimization model described in (1) to (6). Equation (1) corresponds to the objective function of the phase-balancing problem, which is formulated as the minimization of the daily costs of energy loss in all branches of the network. Equations (2) and (3) are the active and reactive power balance constraints that must hold at every phase, node, and time period in the three-phase distribution network. These constraints originate from the three-phase power flow problem in electrical grids, which has been solved numerically using numerical methods such as the Newton–Raphson [31], backward/forward [32], and graph-based methods [33]. Equations (4) and (5) guarantee that the loads are connected to the phases in a unique form by using a matrix connection at each node (i.e., bus k), formed by the variables ${x}_{kfg}$, in which each row and each column must be equal to 1. Finally, the box-type constraint (6) defines the upper and lower voltage regulation bounds of the network, which are typically $\pm 10\%$ for medium-voltage levels [2].

The problem of optimal phase balancing in electrical distribution systems can be solved using a master–slave metaheuristic optimization approach, which represents binary variables of the optimization model (1)–(6) with an integer codification that represents the possible phase connections per node. The phase connection is defined in the master stage, which corresponds to an improved CBGA; in the slave stage, each phase configuration is evaluated to determine the total costs of daily energy losses by using the successive approximation power flow method. Details of the proposed master–slave optimization approach are presented in the following subsections.

The successive approximation power flow method is a recently developed power flow method for single-phase and direct current distribution grids [28]; however, thus far, its extension to three-phase unbalanced distribution grids with Y and $\Delta $ loads has not been reported. The main characteristics of the successive approximation power flow method are that its recursive formula does not use derivatives and the matrices that intervene in the iteration process are constant, implying that the required processing times to solve the power flow problem are in the order of milliseconds [2]; this power flow method can be formulated directly in a complex domain, which simplifies its computational implementation.

To formulate the three-phase successive approximation power flow method, we present the relation between the nodal voltages and the injected currents using the admittance matrix, as follows:
where ${\mathbb{V}}_{s3\phi ,t}$ is the vector that contains all the complex voltages in the slack sources (note that these voltages are clearly known for power flow purposes); ${\mathbb{V}}_{d3\phi ,t}$ is the vector that contains all voltages in the demand nodes (unknown variables of interest); ${\mathbb{I}}_{s3\phi ,t}$ is the vector of the net injected current in the slack nodes; ${\mathbb{I}}_{d3\phi ,t}$ represents the vector that contains all currents in the demand nodes; ${\mathbb{Y}}_{gs3\phi}$ is a component of the admittance matrix that relates the slack sources among them; ${\mathbb{Y}}_{gd3\phi}={\mathbb{Y}}_{ds3\phi}^{T}$ is the component of the admittance matrix that relates the slack and demand nodes and vice versa; and ${\mathbb{Y}}_{dd3\phi}$ is the component of the admittance matrix that relates the demand nodes among them. It should be noted that the voltages and currents in Equation (6) are defined under the three-phase condition and ordered by the nodes and phases, respectively; this equation is applied for each time period.

$$\begin{array}{c}\hfill \left[\begin{array}{c}{\mathbb{I}}_{s3\phi ,t}\\ {\mathbb{I}}_{d3\phi ,t}\end{array}\right]=\left[\begin{array}{cc}{\mathbb{Y}}_{ss3\phi}& {\mathbb{Y}}_{ds3\phi}\\ {\mathbb{Y}}_{sd3\phi}& {\mathbb{Y}}_{dd3\phi}\end{array}\right]\left[\begin{array}{c}{\mathbb{V}}_{s3\phi ,t}\\ {\mathbb{V}}_{d3\phi ,t}\end{array}\right],\phantom{\rule{0.277778em}{0ex}}\forall t\in \mathcal{T}\end{array}$$

From (6), it is evident that the second row contains the demand nodal voltages (i.e., ${\mathbb{V}}_{s3\phi ,t}$), which are the variables of interest in power flow studies; for simplicity, this equation can be rewritten as follows:

$$\begin{array}{c}\hfill {\mathbb{V}}_{d3\phi ,t}=-{\mathbb{Y}}_{dd3\phi}^{-1}\left[{\mathbb{Y}}_{gd3\phi}{\mathbb{V}}_{s3\phi ,t}-{\mathbb{I}}_{d3\phi ,t}\right].\phantom{\rule{0.277778em}{0ex}}\forall t\in \mathcal{T}\end{array}$$

Equation (8) enables the determination of all the voltage profiles at all the demand nodes; however, it is necessary to determine the form of calculating the demand current because this depends on the voltage profile and the amount of power demanded, including whether the loads have $Y-$ and $\Delta -$connections.

The case where node k has a load with a Y-connection (we assume that this is solidly grounded) is depicted in Figure 1.

In the three-phase load with a Y-connection, we assume that the voltage subjected to each load is phase-to-neutral, and each line current is calculated as follows:
which can be compacted as presented below:

$$\begin{array}{c}\hfill {\mathbb{I}}_{dka,t}=-\frac{{\mathbb{S}}_{dka,t}^{\star}}{{\mathbb{V}}_{dka,t}^{\star}},\phantom{\rule{0.277778em}{0ex}}{\mathbb{I}}_{dkb,t}=-\frac{{\mathbb{S}}_{dkb,t}^{\star}}{{\mathbb{V}}_{dkb,t}^{\star}},\phantom{\rule{0.277778em}{0ex}}{\mathbb{I}}_{dkc,t}=-\frac{{\mathbb{S}}_{dkc,t}^{\star}}{{\mathbb{V}}_{dkc,t}^{\star}},\phantom{\rule{0.277778em}{0ex}}\forall t\in \mathcal{T}\end{array}$$

$$\begin{array}{c}\hfill {\mathbb{I}}_{dk3\phi ,t}=-{\mathbf{diag}}^{-1}\left({\mathbb{V}}_{dk3\phi ,t}^{\star}\right){\mathbb{S}}_{dk3\phi ,t}^{\star}.\phantom{\rule{0.277778em}{0ex}}\forall t\in \mathcal{T}\end{array}$$

To obtain the line currents for a three-phase network with loads connected in $\Delta $, the load diagram presented in Figure 2 is considered for the k node.
which can be compacted as follows:
where the matrices $\mathbf{M}$ and $\mathbf{H}$ take the following form:

$$\begin{array}{c}\hfill {\mathbb{I}}_{dka,t}=-{\left(\frac{{\mathbb{S}}_{kab,t}}{{\mathbb{V}}_{dka,t}-{\mathbb{V}}_{dkb,t}}\right)}^{\star}+{\left(\frac{{\mathbb{S}}_{kca,t}}{{\mathbb{V}}_{dkc,t}-{\mathbb{V}}_{dka,t}}\right)}^{\star},\phantom{\rule{0.277778em}{0ex}}\forall t\in \mathcal{T}\\ \hfill {\mathbb{I}}_{dkb,t}=-{\left(\frac{{\mathbb{S}}_{kbc,t}}{{\mathbb{V}}_{dkb,t}-{\mathbb{V}}_{dkc,t}}\right)}^{\star}+{\left(\frac{{\mathbb{S}}_{kab,t}}{{\mathbb{V}}_{dka,t}-{\mathbb{V}}_{dkb,t}}\right)}^{\star},\phantom{\rule{0.277778em}{0ex}}\forall t\in \mathcal{T}\\ \hfill {\mathbb{I}}_{dkc,t}=-{\left(\frac{{\mathbb{S}}_{kca,t}}{{\mathbb{V}}_{dkc,t}-{\mathbb{V}}_{dka,t}}\right)}^{\star}+{\left(\frac{{\mathbb{S}}_{kbc,t}}{{\mathbb{V}}_{dkb,t}-{\mathbb{V}}_{dkc,t}}\right)}^{\star},\phantom{\rule{0.277778em}{0ex}}\forall t\in \mathcal{T}\end{array}$$

$$\begin{array}{c}\hfill {\mathbb{I}}_{dk3\phi ,t}=-\left({\mathbf{diag}}^{-1}\left(\mathbf{M}{\mathbb{V}}_{dk3\phi ,t}^{\star}\right)-{\mathbf{diag}}^{-1}\left({\mathbf{M}}^{T}{\mathbb{V}}_{dk3\phi ,t}^{\star}\right)\mathbf{H}\right){\mathbb{S}}_{dk3\phi ,t}^{\star},\phantom{\rule{0.277778em}{0ex}}\forall t\in \mathcal{T}\end{array}$$

$$\begin{array}{c}\hfill \mathbf{M}=\left[\begin{array}{ccc}1& -1& 0\\ 0& 1& -1\\ -1& 0& 1\end{array}\right],\phantom{\rule{0.277778em}{0ex}}\mathbf{H}=\left[\begin{array}{ccc}0& 0& 1\\ 1& 0& 0\\ 0& 1& 0\end{array}\right]\end{array}$$

Considering a case where a three-phase distribution grid has loads connected with the Y and $\Delta $ structures, for calculating the demanded currents, these connections must be considered, as presented in Equations (9) and (10), respectively. When solving the power flow problem, the main goal is to determine the daily energy losses for each configuration of the phases provided in the master stage. To this end, we employ the admittance nodal matrix, as follows:
where ${\mathbb{S}}_{\mathrm{loss}}$ represents the apparent power losses of the network.

$$\begin{array}{c}\hfill {\mathbb{S}}_{\mathrm{loss},t}={\left[{\mathbb{V}}_{s3\phi ,t},{\mathbb{V}}_{d3\phi ,t}\right]}^{T}{\mathbb{Y}}_{3\phi}^{\mathrm{bus}}\left[{\mathbb{V}}_{s3\phi ,t},{\mathbb{V}}_{d3\phi ,t}\right]\Delta t,\end{array}$$

Algorithm 1 shows the general implementation of the successive approximation power flow for unbalanced distribution networks with loads connected in Y and $\Delta $.

Note that the variable ${\mathbb{S}}_{\mathrm{loss},t}$ contains information regarding the energy losses of the system for the phase configuration provided by the master stage at each time period. To calculate the objective function defined in Equation (1), we use the value of the daily energy losses provided by the three-phase successive approximation power flow in Algorithm 1, as follows:

$$\begin{array}{c}\hfill {f}_{1}={C}_{\mathrm{kWh}}\sum _{t\in \mathcal{T}}{\mathbb{S}}_{\mathrm{loss},t}\Delta t\end{array}$$

Note that expression (12) presents the connection of the master and slave stages because the slave stage is necessary to determine the objective function value of each phase configuration for the loads provided by the proposed improved CBGA.

Algorithm 1: Solution of the three-phase power flow problem for unbalanced distributions networks with Y and $\Delta $ loads. |

The master stage is entrusted with providing the phase configuration to the successive approximation power flow approach defined in Algorithm 1. Here, we propose an improved version of the classical CBGA in the selection and mutation stages by introducing a normal Gaussian distribution, which is used in the vortex search algorithm [34]. To explain the improvements in the CBGA, we begin with the codification of the phase-balancing problem, which is based on the possible connections reported in Table 1.

The proposed codification takes the following structure:
where ${x}_{i}^{s}$ is the individual i in iteration s, which has a dimension of $1\times (n-1)$. The main advantage of this codification is that the 9 variables associated with the binary variables ${x}_{kfg}$ per node are represented solely by an integer between 1 and 6. The improved version of the CBGA is developed using probability criteria at each iteration, which defines the methodology used to generate the offspring. If the probability criterion is lower than $50\%$ in iteration s, the classical CBGA with multi-point mutation is implemented; otherwise, the vortex search algorithm is used to generate the offspring.

$$\begin{array}{c}\hfill {x}_{i}^{s}=\left[1,\phantom{\rule{0.277778em}{0ex}}5,\phantom{\rule{0.277778em}{0ex}}6,\phantom{\rule{0.277778em}{0ex}}2,\phantom{\rule{0.277778em}{0ex}}\cdots ,\phantom{\rule{0.277778em}{0ex}}3\right],\phantom{\rule{0.277778em}{0ex}}\end{array}$$

The classical approach is based on the selection of two individuals from the population, that is, ${x}_{i}^{s}$ and ${x}_{j}^{s}$ being $i\ne j$ from the current population. These two solutions are recombined using an arbitrary point randomly generated as a number between 1 and $n-2$. To illustrate this stage, suppose that a test feeder is composed of 11 nodes (1 slack and 10 demands), in which the individuals ${x}_{i}^{s}$ and ${x}_{j}^{s}$ take the following form:

$$\begin{array}{c}\hfill {x}_{i}^{s}=\left[1,\phantom{\rule{0.277778em}{0ex}}5,\phantom{\rule{0.277778em}{0ex}}6,\phantom{\rule{0.277778em}{0ex}}2,\phantom{\rule{0.277778em}{0ex}}4,\phantom{\rule{0.277778em}{0ex}}3,\phantom{\rule{0.277778em}{0ex}}4,\phantom{\rule{0.277778em}{0ex}}1,\phantom{\rule{0.277778em}{0ex}}5,\phantom{\rule{0.277778em}{0ex}}6\right],\\ \hfill {x}_{j}^{s}=\left[2,\phantom{\rule{0.277778em}{0ex}}6,\phantom{\rule{0.277778em}{0ex}}6,\phantom{\rule{0.277778em}{0ex}}4,\phantom{\rule{0.277778em}{0ex}}4,\phantom{\rule{0.277778em}{0ex}}4,\phantom{\rule{0.277778em}{0ex}}1,\phantom{\rule{0.277778em}{0ex}}5,\phantom{\rule{0.277778em}{0ex}}5,\phantom{\rule{0.277778em}{0ex}}2\right].\end{array}$$

Considering that the recombination point for the iteration s is 4, the offspring takes the following form:

$$\begin{array}{c}\hfill {y}_{i}^{s}=\left[1,\phantom{\rule{0.277778em}{0ex}}5,\phantom{\rule{0.277778em}{0ex}}6,\phantom{\rule{0.277778em}{0ex}}2,\phantom{\rule{0.277778em}{0ex}}4,\phantom{\rule{0.277778em}{0ex}}4,\phantom{\rule{0.277778em}{0ex}}1,\phantom{\rule{0.277778em}{0ex}}5,\phantom{\rule{0.277778em}{0ex}}5,\phantom{\rule{0.277778em}{0ex}}2\right],\\ \hfill {y}_{j}^{s}=\left[2,\phantom{\rule{0.277778em}{0ex}}6,\phantom{\rule{0.277778em}{0ex}}6,\phantom{\rule{0.277778em}{0ex}}4,\phantom{\rule{0.277778em}{0ex}}4,\phantom{\rule{0.277778em}{0ex}}3,\phantom{\rule{0.277778em}{0ex}}4,\phantom{\rule{0.277778em}{0ex}}1,\phantom{\rule{0.277778em}{0ex}}5,\phantom{\rule{0.277778em}{0ex}}6\right].\end{array}$$

For each descending individual, the multi-point mutation criteria are applied as follows: (i) the number of mutation points is randomly selected from 1 to $1+20\%(n-1)rand\left(1\right)$ of the population size (i.e., 1 and 3). Let us suppose that, for individual ${y}_{i}^{s}$, the number of mutation points is 2, whereas for individual ${y}_{j}^{s}$, the number of mutation points is 3. Thus, the selection of mutation points is randomly generated for each descending individual; here, we consider that these points are 1 and 8 for ${y}_{i}^{s}$ and 2, 5, and 7 for ${y}_{j}^{s}$. A random number between 1 and 6 is generated at each of these points to define the phase connection at this node. This procedure yields the following structure in the current offspring:

$$\begin{array}{c}\hfill {y}_{i}^{s}=\left[3,\phantom{\rule{0.277778em}{0ex}}5,\phantom{\rule{0.277778em}{0ex}}6,\phantom{\rule{0.277778em}{0ex}}2,\phantom{\rule{0.277778em}{0ex}}4,\phantom{\rule{0.277778em}{0ex}}4,\phantom{\rule{0.277778em}{0ex}}1,\phantom{\rule{0.277778em}{0ex}}6,\phantom{\rule{0.277778em}{0ex}}5,\phantom{\rule{0.277778em}{0ex}}2\right],\\ \hfill {y}_{j}^{s}=\left[2,\phantom{\rule{0.277778em}{0ex}}1,\phantom{\rule{0.277778em}{0ex}}6,\phantom{\rule{0.277778em}{0ex}}4,\phantom{\rule{0.277778em}{0ex}}5,\phantom{\rule{0.277778em}{0ex}}3,\phantom{\rule{0.277778em}{0ex}}3,\phantom{\rule{0.277778em}{0ex}}1,\phantom{\rule{0.277778em}{0ex}}5,\phantom{\rule{0.277778em}{0ex}}6\right].\end{array}$$

The remaining steps involve the evaluation of the objective function (12) to select which among the individuals ${y}_{i}^{s}$ and ${y}_{j}^{s}$ has the opportunity to replace the current population, considering the case where the objective function is better than the worst individual in the population and is different from all the other ones (diversity criteria).

In this stage, the main goal is to add to the classical CBGA the possibility of generating an offspring population based on the evolution criteria of the VSA [2]. The main characteristic of the VSA is the usage of a normal Gaussian distribution to generate a neighborhood around the current solution, called the center of the hyper-ellipse $\mu $[35]. Here, we select the center of the hyper-ellipse as a random individual in the current population during iteration s as ${\mu}_{S}={x}_{i}^{s}$. This center is generated using the following distribution probability:
where $d=(n-1)$ represents the dimension, y is a random vector featuring values between zero and one with dimension $d\times 1$, and $\Sigma $ is the covariance matrix. If the values of the diagonal elements of $\Sigma $ are equal and those outside the diagonal are zero, then the resulting shape of the distribution will be a hyper-ellipse; therefore, equal variances with zero covariance can be calculated as the value of $\Sigma $, as presented in Equation (15).
where $\sigma $ represents the variance of the distribution, and I represents the identity matrix of dimension $d\times d$. Equation (16) can be used to calculate the current standard departure ${\sigma}_{s}$ of the distribution for iteration s.
where ${\sigma}_{s}$ is also considered as the radius ${r}_{s}$ of the current hyper-ellipse. As a first step, by multiplying the radius, the region of possible solutions for the offspring is generated.

$${C}_{i}^{s}\left(y\right)=p\left(y|\mu ,\Sigma \right)=\frac{1}{\sqrt{{\left(2\pi \right)}^{d}\left|\Sigma \right|}}exp\left\{-\frac{1}{2}{\left(y-{\mu}_{s}\right)}^{T}{\Sigma}^{-1}\left(y-{\mu}_{s}\right)\right\}$$

$$\Sigma ={\sigma}^{2}\xb7{\left[I\right]}_{d\times d}$$

$${\sigma}_{s}=\frac{max\left\{{x}^{max}\right\}-min\left\{{x}^{min}\right\}}{2}$$

The set of neighboring solutions obtained and contained in ${C}_{i}^{s}\left(y\right)$, which are verified using the algorithm and by considering the lower and upper limits set in ${\mu}_{s}$, is defined by Equation (17).
where $rand$ is a function that generates a random number between 1 and 0 with a normal distribution. After verifying the limits, at this point, the exploiter stage is initiated. It selects the optimal response within the solution space of ${C}_{i}^{s}$ obtained in (17), which is selected as the potential individual to be part of the current population of the CBGA, provided it is better than the worst solution in the current population and fulfills the diversity criteria.

$${C}_{i}^{s}\left(y\right)=\left\{\begin{array}{c}{C}_{i}^{s}\left(y\right),\phantom{\rule{4pt}{0ex}}\phantom{\rule{4pt}{0ex}}\phantom{\rule{4pt}{0ex}}\phantom{\rule{4pt}{0ex}}\phantom{\rule{4pt}{0ex}}{x}^{min}\le x\le {x}^{max}\\ {x}^{min}+\left({x}^{max}-{x}^{min}\right)rand,\phantom{\rule{4pt}{0ex}}\phantom{\rule{4pt}{0ex}}\phantom{\rule{4pt}{0ex}}\phantom{\rule{4pt}{0ex}}\phantom{\rule{4pt}{0ex}}\mathrm{otherwise}\end{array}\right.$$

It should be noted that the radius of the hyper-ellipse is gradually reduced as the optimization process progresses. Here, we selected a linear rule to decrease the radius as follows:
where ${s}_{max}$ is the maximum number of iterations of the CBGA.

$$\begin{array}{c}\hfill {r}_{s}=1-\frac{s}{{s}_{max}},\end{array}$$

The size of the offspring in ${C}_{i}^{s}\left(y\right)$ is selected as 20% of the population size in order to reduce the number of evaluations required in the slave stage for determining the costs of daily energy losses.

The implementation of the proposed improved CBGA considering multi-point mutation and offspring generation with the VSA is presented in Algorithm 2.

Algorithm 2: Improved CBGA for solving the phase-balancing problem in three-phase unbalanced distribution networks. |

Two test feeders are considered to validate the proposed optimization approach based on the improved version of the CBGA with the multi-point mutation strategy and the VSA for generating offspring. The main characteristics of the test feeders are discussed below.

The 15-bus test feeder is a three-phase unbalanced distribution network composed of 15 nodes and 14 distribution lines, which are operated with a phase voltage of 13.2 kV at the substation (typical voltage in Colombian distribution grids). The electrical configuration of the test feeder is depicted in Figure 3.

This electric distribution grid is a portion of an underground power system in California, USA, and is composed of 37 buses and 36 lines, as depicted in Figure 4. In this test system, the slack source is assigned at node 1 with a line-to-line voltage of 4.8 kV. The total power consumptions for phases a, b, and c are 727 kW and 357 kvar, 639 kW and 314 kvar, and 1091 kW and 530 kvar, respectively. Note that the IEEE 37-bus test system employed in this research is the adaptation proposed in [2], where the voltage regulator between nodes 1 and 2 is replaced by a single three-phase line with a length of 1850 feet, and the transformers located at nodes 10 and 24 are removed, including bus 24, as recommended by [13].

This section presents the numerical validation of the proposed approach in the 15- and IEEE 37-bus test feeders, considering the following scenarios: (i) application of the improved CBGA and comparative methodologies to the 15-bus test system under the peak load condition, and (ii) usage of active and reactive demand load curves in the IEEE 37-bus test feeder to minimize the annual operative costs of the network.

In this simulation case, we present the positive effects of the power loss reduction after applying a phase-balancing approach in the 15-bus test feeder, considering the peak load scenario. Table 6 presents the optimal solutions obtained via the classical CBGA [27], the sine-cosine algorithm (SCA) [36], the black hole optimizer (BHO) [37], the VSA [2], and the improved CBGA. Notably, the proposed and comparative algorithms were set with 10 individuals in the population, 1000 iterations, and 100 consecutive evaluations to determine the average processing time.

The results in Table 6 show that (i) all the methodologies enable a reduction of more than 18% in the amount of power loss in the peak load case. The worst result was yielded by the BHO, with a reduction of $18.06\%$, whereas the best result was achieved by the improved CBGA, with a reduction of $18.66\%$. All the compared methods, including the proposed approach, require less than 40 s to identify an adequate solution by using the phase-balancing configuration. Nevertheless, the proposed approach is the fastest, with an average processing time of $2.0760$ s. (iii) The classical CBGA and the VSA afford the second and third best results, respectively. Thus, the combination of the best characteristics of both these methods in the improved CBGA enables excellent numerical results for the phase-balancing problem.

Note the following in Table 6: our proposal-the proposed improved version of the CBGA-is at least three times faster than the classical CBGA, which is the second faster method; in addition, the VSA method proposed in [2] is more than nineteen times slower than our proposed optimization approach. Those comparisons confirms the effectiveness of the approach here presented to solve complex MINLP problems with solution spaces with large dimensions.

To verify the effectiveness of the redistribution of all the loads in the distribution network, Figure 5 presents the amount of power loss per phase for the benchmark case, as well as for the solution obtained using the improved CBGA.

The variations in the power losses per phase, as shown in Figure 5, indicate that, after applying phase balancing with the improved version of the CBGA, the losses of phases b and c increase by $9.1521\phantom{\rule{3.33333pt}{0ex}}$kW and $21.8832$ kW, respectively. By contrast, the power loss of phase c decreases by approximately $56.0845$ kW, which clearly compensates for the increments experienced in phases a and b. In addition, the power losses per phase are balanced at approximately 36 kW, with small differences lower than 2 kW; this shows that the phase-balancing approach using the improved CBGA effectively redistributes the loads in all the phases as uniformly as possible.

In the case where the 15-bus system has all the loads connected in $\Delta $, the power loss in the benchmark case is approximately $110.7841$ kW, which can be reduced to $107.0944$ kW by using the proposed improved version of the CBGA. The codification that presents this solution is $\left\{5,6,6,2,3,2,3,5,5,2,4,3,5,6\right\}$. The difference between the loads connected in Y and $\Delta $ is notable in terms of the amount of power loss in the benchmark case, because the difference between both connections is approximately $23.4631$ kW. In addition, the percentage of power loss reduction in the case of $\Delta -$connections is just $3.33\%$, which is constrained by the $Y-$connection case, where this reduction exceeds $18\%$. The average processing time in the case of $\Delta -$connections is $1.6199\phantom{\rule{3.33333pt}{0ex}}$s, which is $0.4563$ s less than that for the $Y-$connections.

To demonstrate the effectiveness of the improved version of the CBGA for performing the phase-balancing task in three-phase distribution networks, the IEEE-37 node test feeder was evaluated under a daily power flow environment, considering active and reactive power demand curves. The demand curves are presented in Figure 6 and also reported in Table 7, for a comparison of our results with future research (note that the scaling factor of the data reported in Table 7 for active and reactive power demands is 2).

To evaluate the objective function in terms of the annual operative costs, the cost of the energy is assumed in US$/kWh $0.1390$ (this value is the average cost of the energy in Bogotá, Colombia in May 2019 [38]); the number of days is considered as 365, and the length of the power flow period, ${\Delta}_{t}$, is 0.5 h.

Table 8 presents the numerical performance of the improved version of the CBGA when 10 individuals, 1000 iterations, and 100 consecutive evaluations are adopted for the IEEE 37-bus system. Based on the results reported in Table 8, it is evident that (i) the best numerical solution (see solution 1) enables a reduction in the annual operative costs by approximately US$/year $8121.7220$, which corresponds to a reduction of $18.79\%$ with respect to the benchmark case; and (ii) the difference between solution 1 and 10 is just US$/year $75.1586$, which implies that the 10 solutions reported in the final population of the improved CBGA are adequate to solve the phase-balancing problem of three-phase networks. The implementation of any one of these depends exclusively on the practices of the utility.

To demonstrate the effectiveness of the proposed optimization approach, Figure 7 presents the percentages of annual reduction in the energy costs for each of the ten solutions reported in Table 8. Note that the differences among all solutions are lower than $0.18\%$, thereby confirming that the proposed approach can efficiently solve the phase-balancing problem of unbalanced electrical three-phase grids. In addition, it is important to mention that the average processing time for solving the phase-balancing problem in the IEEE 37-bus system is $136.3901\phantom{\rule{3.33333pt}{0ex}}$s. This is a short time, considering that, for each combination of the phases provided by the CBGA, 48 power flows are evaluated, which implies that, for a population of 10 individuals, 1000 iterations with two offspring individuals per iteration are evaluated as approximately 96,480 power flows.

To elucidate the effect of phase balancing on the amount of daily energy loss, shown in Figure 8, we present the daily energy losses per phase before and after applying the phase-balancing strategy (solution 1 in Table 8).

The results presented in Figure 8 indicate the following: (i) the benchmark case shows that the difference in the daily energy losses between phases b and c is approximately $562.2398$ kWh/day, which confirms the higher imbalance in the IEEE 37-bus system; (ii) after applying the phase-balancing strategy, the difference between phases b and c is just $127.1357$ kWh/day, which represents a reduction of $77.3876\%$ with respect to the benchmark case; and (iii) the phase-balancing approach leads to better load distribution per phase, which balances the power losses among them with differences lower than 130 kWh/day in the worst case.

In this research, the problem of phase balancing in three-phase distribution grids was addressed from the master–slave optimization perspective. To this end, an improved version of the CBGA using the VSA was proposed. The master stage provides the configuration of the load phases through the CBGA by using an integer codification between 1 and 6, which represents all the possible load connections in a three-phase grid. Alternatively, the slave stage evaluates each of the load configurations in the three-phase power flow tool, corresponding to an extended version of the successive approximation power flow method for three-phase unbalanced networks.

Numerical results for the 15-bus test system demonstrate that the improved version of the CBGA proposed herein attains the highest quality solution with a reduction of $18.66\%$, as compared with the classical CBGA ($18.64\%$), BHO ($18.06\%$), SCA ($18.51\%$), and VSA ($18.57\%$). In addition, the proposed approach exhibited the fastest processing time ($2.0762\phantom{\rule{3.33333pt}{0ex}}$s) for solving the phase-balancing problem in the 15-bus system, followed by the classical CBGA ($6.9435\phantom{\rule{3.33333pt}{0ex}}$s). These processing times are important because the dimension space of the phase-balancing problem potentially increases, that is, ${6}^{\left(n-1\right)}$, implying that the 15-bus system corresponds to ${6}^{14}$ = 78,364,164,096, that is, more than 78,000 million combinations.

In the IEEE 37-bus system, the improved version of the CBGA reduces the annual energy loss cost by approximately $18.79\%$, with an average processing time of $136.3901$ s. In addition, the difference between the ten best solutions reported by the proposed CBGA is less than $0.18\%$, which is less than 76 dollars per year of operation. This confirms the effectiveness of the proposed approach in solving complex MINLP models, such as the phase-balancing problem in three-phase networks.

In future works, it will be possible to (i) combine the proposed three-phase-balancing approach with the optimal location of capacitor banks to reduce the annual operating costs of three-phase grids; (ii) study the optimal integration of static distribution compensators in three-phase networks by using a hybrid master-slave approach employing the successive approximation power flow method in the slave stage; and (iii) use the improved version of the CBGA to solve the problem of the optimal selection of conductors in three-phase networks.

Conceptualization, O.D.M., A.M.-C. and L.F.G.-N.; methodology, O.D.M., A.M.-C. and L.F.G.-N.; investigation, O.D.M., A.M.-C. and L.F.G.-N.; writing—review and editing, O.D.M., A.M.-C., L.F.G.-N., R.A.H. and M.G. All authors have read and agreed to the published version of the manuscript.

This research received no external funding.

Not applicable.

Not applicable.

No new data were created or analyzed in this study. Data sharing is not applicable to this article.

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.”

The authors declare no conflict of interest.

Here is presented a list with the main abbreviations used along within this document.

CBGA | Chu & Beasley genetic algorithm |

VSA | Vortex search algorithm |

MINLP | Mixed-integer nonlinear programming |

CREG | Energy and Gas Regulation Commission |

SCA | Sine-cosine algorithm |

BHO | Black-hole optimizer |

${\Delta}_{t}$ | Time period during which the power demands remain constant (h). |

${\delta}_{kft}$ | Angle of the voltage at node k in phase f at time period t (rad). |

${\delta}_{mgt}$ | Angle of the voltage at node m in phase g at time period t (rad). |

${\mathbb{I}}_{d3\phi ,t}$ | Complex three-phase current at demand node d in time period t (A). |

${\mathbb{I}}_{dk3\phi ,t}$ | Complex demanded three-phase current at k in time period t (A). |

${\mathbb{I}}_{dka,t}$ | Complex demanded current at node k in phase a at time period t (A). |

${\mathbb{I}}_{dkb,t}$ | Complex demanded current at node k in phase b at time period t (A). |

${\mathbb{I}}_{dkc,t}$ | Complex demanded current at node k in phase c at time period t (A). |

${\mathbb{I}}_{s3\phi ,t}$ | Complex three-phase current at source node s in time period t (A). |

${\mathbb{S}}_{\mathrm{loss}}$ | Apparent power losses of the network (VA) |

${\mathbb{S}}_{dk3\phi ,t}$ | Complex demanded three-phase power at k in time period t (VA). |

${\mathbb{S}}_{dka,t}$ | Complex demanded power at node k in phase a at time period t (VA). |

${\mathbb{S}}_{dkab,t}$ | Complex demanded power at node k between phases a and b at time period |

t (VA). | |

${\mathbb{S}}_{dkb,t}$ | Complex demanded power at node k in phase b at time period t (VA). |

${\mathbb{S}}_{dkbc,t}$ | Complex demanded power at node k between phases b and c at time period |

t (VA). | |

${\mathbb{S}}_{dkc,t}$ | Complex demanded power at node k in phase c at time period t (VA). |

${\mathbb{S}}_{dkca,t}$ | Complex demanded power at node k between phases c and a at time period |

t (VA). | |

${\mathbb{V}}_{d3\phi ,t}$ | Complex three-phase voltage at demand node d in time period t (V). |

${\mathbb{V}}_{dk3\phi ,t}$ | Complex demanded three-phase voltage at k in time period t (V). |

${\mathbb{V}}_{dka,t}$ | Complex demanded voltage at node k in phase a at time period t (V). |

${\mathbb{V}}_{dkb,t}$ | Complex demanded voltage at node k in phase b at time period t (V). |

${\mathbb{V}}_{dkc,t}$ | Complex demanded voltage at node k in phase c at time period t (V). |

${\mathbb{V}}_{s3\phi ,t}$ | Complex three-phase voltage at source node s in time period t (V). |

${\mathbb{Y}}_{dd3\phi}$ | Three-phase submatrix of the admittance nodal matrix that relates demand |

nodes among them (S). | |

${\mathbb{Y}}_{ds3\phi}$ | Three-phase submatrix of the admittance nodal matrix that relates demand |

and source nodes among them (S). | |

${\mathbb{Y}}_{sd3\phi}$ | Three-phase submatrix of the admittance nodal matrix that relates source |

and demand nodes among them (S). | |

${\mathbb{Y}}_{ss3\phi}$ | Three-phase submatrix of the admittance nodal matrix that relates source |

nodes among them (S). | |

$\mathcal{F}$ | Set containing all the phases of the system. |

$\mathcal{N}$ | Set containing all the nodes of the network. |

$\mathcal{T}$ | Set containing all the time periods of the operation horizon. |

$\mathbf{M}$ | Matrix to calculate the line-to-line voltages from line-to-ground voltages. |

$\mathbf{H}$ | Matrix of demand rotation. |

${\theta}_{kfmg}$ | Angle of the admittance that relates node k at phase f with node m at phase |

g (rad). | |

${C}_{\mathrm{kWh}}$ | Average cost of energy (US$/kWh). |

${f}_{1}$ | Cost of daily energy losses (US$/day). |

${P}_{kft}^{s}$ | Active power generated by source s connected at node k in phase f at time |

period t (W). | |

${P}_{kgt}^{d}$ | Active power demanded at node k in phase g at time |

period t (W). | |

${Q}_{kft}^{s}$ | Reactive power generated by source s connected at node k in phase f at time |

period t (var). | |

${Q}_{kgt}^{d}$ | Reactive power demanded at node k in phase g at time |

period t (var). | |

${V}_{max}$ | Maximum voltage regulation bound (V). |

${V}_{min}$ | Minimum voltage regulation bound (V). |

${V}_{kft}$ | Voltage magnitude at node k in phase f at time period t (V). |

${V}_{mgt}$ | Voltage magnitude at node m in phase g at time period t (V). |

${x}_{kfg}$ | Binary variable that defines the connection of the demand in node k at f in |

phase g. | |

${Y}_{kfmg}$ | Magnitude of admittance that relates node k at phase f with node m at |

phase g (S). |

- Temiz, A.; Almalki, A.M.; Kahraman, Ö.; Alshahrani, S.S.; Sönmez, E.B.; Almutairi, S.S.; Nadar, A.; Smiai, M.S.; Alabduljabbar, A.A. Investigation of MV Distribution Networks with High-Penetration Distributed PVs: Study for an Urban Area. Energy Procedia
**2017**, 141, 517–524. [Google Scholar] [CrossRef] - 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. [Google Scholar] [CrossRef] - Aboshady, F.M.; Thomas, D.W.P.; Sumner, M. A Wideband Single End Fault Location Scheme for Active Untransposed Distribution Systems. IEEE Trans. Smart Grid
**2020**, 11, 2115–2124. [Google Scholar] [CrossRef] - Arias, J.; Calle, M.; Turizo, D.; Guerrero, J.; Candelo-Becerra, J. Historical Load Balance in Distribution Systems Using the Branch and Bound Algorithm. Energies
**2019**, 12, 1219. [Google Scholar] [CrossRef] - 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. [Google Scholar] [CrossRef] - Cabrera, J.B.; Veiga, M.F.; Morales, D.X.; Medina, R. Reducing Power Losses in Smart Grids with Cooperative Game Theory. In Advanced Communication and Control Methods for Future Smartgrids; IntechOpen: London, UK, 2019. [Google Scholar] [CrossRef]
- 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. [Google Scholar] [CrossRef] - 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. [Google Scholar] [CrossRef] - Al-Sumaiti, A.S.; 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. [Google Scholar] [CrossRef] - 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. [Google Scholar] [CrossRef] - Grigoraș, G.; Neagu, B.C.; Gavrilaș, M.; Triștiu, I.; Bulac, C. Optimal Phase Load Balancing in Low Voltage Distribution Networks Using a Smart Meter Data-Based Algorithm. Mathematics
**2020**, 8, 549. [Google Scholar] [CrossRef] - Boroujeni, S.T.; Mardaneh, M.; Hashemi, Z. A Dynamic and Heuristic Phase Balancing Method for LV Feeders. Appl. Comput. Intell. Soft Comput.
**2016**, 2016, 1–8. [Google Scholar] [CrossRef] - 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. [Google Scholar] [CrossRef] - Montoya, O.D.; Grajales, A.; Hincapié, R.A.; Granada, M. A new approach to solve the distribution system planning problem considering automatic reclosers. Ingeniare. Rev. Chil. Ing.
**2017**, 25, 415–429. [Google Scholar] [CrossRef] - Garcés, A.; Castaño, J.C.; Rios, M.A. Phase Balancing in Power Distribution Grids: A Genetic Algorithm with a Group-Based Codification. In Energy Systems; Springer International Publishing: Cham, Switzerland, 2020; pp. 325–342. [Google Scholar] [CrossRef]
- Darmawan, I.; Kuspriyanto, Y.; Priyan, M.I.J. Integration of Genetic and Tabu Search algorithm based load balancing for heterogenous grid computing. In Proceedings of the 2013 International Conference on Computer, Control, Informatics and Its Applications (IC3INA), Jakarta, Indonesia, 19–21 November 2013. [Google Scholar] [CrossRef]
- Li, D.; Qiang, R.; Yoon, S.W. Balanced Adaptive Tabu Search Algorithm to Optimize Dual-gantry Pick-and-place Assembly. Procedia Manuf.
**2017**, 11, 1892–1899. [Google Scholar] [CrossRef] - Sim, K.M.; Sun, W.H. Ant colony optimization for routing and load-balancing: Survey and new directions. IEEE Trans. Syst. Man Cybern. Part A Syst. Hum.
**2003**, 33, 560–572. [Google Scholar] [CrossRef] - Keskinturk, T.; Yildirim, M.B.; Barut, M. An ant colony optimization algorithm for load balancing in parallel machines with sequence-dependent setup times. Comput. Oper. Res.
**2012**, 39, 1225–1235. [Google Scholar] [CrossRef] - Huang, M.Y.; Chen, C.S.; Lin, C.H.; Kang, M.S.; Chuang, H.J.; Huang, C.W. Three-phase balancing of distribution feeders using immune algorithm. IET Gener. Transm. Distrib.
**2008**, 2, 383. [Google Scholar] [CrossRef] - 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. [Google Scholar] [CrossRef] - Amon, D.A. A Modified Bat Algorithm for Power Loss Reduction in Electrical Distribution System. Telkomnika Indones. J. Electr. Eng.
**2015**, 14. [Google Scholar] [CrossRef] - Toma, N.; Ivanov, O.; Neagu, B.; Gavrila, M. A PSO Algorithm for Phase Load Balancing in Low Voltage Distribution Networks. In Proceedings of the 2018 International Conference and Exposition on Electrical Furthermore, Power Engineering (EPE), Iasi, Romania, 18–19 October 2018. [Google Scholar] [CrossRef]
- Schweickardt, G.; Alvarez, J.M.G.; Casanova, C. Metaheuristics approaches to solve combinatorial optimization problems in distribution power systems. An application to Phase Balancing in low voltage three-phase networks. Int. J. Electr. Power Energy Syst.
**2016**, 76, 1–10. [Google Scholar] [CrossRef] - 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. [Google Scholar] [CrossRef] - Zhu, J.; Bilbro, G.; Chow, M.Y. Phase balancing using simulated annealing. IEEE Trans. Power Syst.
**1999**, 14, 1508–1513. [Google Scholar] [CrossRef] - Montoya, O.D.; Gil-González, W.; Orozco-Henao, C. Vortex search and Chu-Beasley genetic algorithms for optimal location and sizing of distributed generators in distribution networks: A novel hybrid approach. Eng. Sci. Technol. Int. J.
**2020**, 23, 1351–1363. [Google Scholar] [CrossRef] - Montoya, O.D.; Gil-González, W. On the numerical analysis based on successive approximations for power flow problems in AC distribution systems. Electr. Power Syst. Res.
**2020**, 187, 106454. [Google Scholar] [CrossRef] - 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. [Google Scholar] [CrossRef] - Da Cunha, A.S.; Peixoto, F.C.; Prata, D.M. Robust data reconciliation in chemical reactors. Comput. Chem. Eng.
**2021**, 145, 107170. [Google Scholar] [CrossRef] - Sereeter, B.; Vuik, K.; Witteveen, C. Newton Power Flow Methods for Unbalanced Three-Phase Distribution Networks. Energies
**2017**, 10, 1658. [Google Scholar] [CrossRef] - Shen, T.; Li, Y.; Xiang, J. A Graph-Based Power Flow Method for Balanced Distribution Systems. Energies
**2018**, 11, 511. [Google Scholar] [CrossRef] - Oliveira-De-Jesus, P.M.D.; Rojas, A.A.; Gonzalez-Longatt, F.M. Unbalanced Power Flow Analysis in Distribution Systems Using TRX Matrix: Implementation Using DIgSILENT Programming Language. In PowerFactory Applications for Power System Analysis; Springer International Publishing: Cham, Switzerland, 2014; pp. 85–110. [Google Scholar] [CrossRef]
- Doğan, B.; Ölmez, T. Vortex search algorithm for the analog active filter component selection problem. AEU Int. J. Electron. Commun.
**2015**, 69, 1243–1253. [Google Scholar] [CrossRef] - Li, P.; Zhao, Y. A quantum-inspired vortex search algorithm with application to function optimization. Nat. Comput.
**2018**, 18, 647–674. [Google Scholar] [CrossRef] - Montoya, O.D.; Molina-Cabrera, A.; Chamorro, H.R.; Alvarado-Barrios, L.; Rivas-Trujillo, E. A Hybrid Approach Based on SOCP and the Discrete Version of the SCA for Optimal Placement and Sizing DGs in AC Distribution Networks. Electronics
**2020**, 10, 26. [Google Scholar] [CrossRef] - Deeb, H.; Sarangi, A.; Mishra, D.; Sarangi, S.K. Improved Black Hole optimization algorithm for data clustering. J. King Saud Univ. Comput. Inf. Sci.
**2020**. [Google Scholar] [CrossRef] - Montoya, O.D.; Gil-González, W.; Grisales-Noreña, L.; Orozco-Henao, C.; Serra, F. Economic Dispatch of BESS and Renewable Generators in DC Microgrids Using Voltage-Dependent Load Models. Energies
**2019**, 12, 4494. [Google Scholar] [CrossRef]

Connection Type | Phases | Sequence |
---|---|---|

1 | ABC | |

2 | CAB | No change |

3 | BCA | |

4 | ACB | |

5 | BAC | Change |

6 | CBA |

Line | Node i | Node j | Cond. | Length (ft) | ${\mathit{P}}_{\mathbf{ja}}$ | ${\mathit{Q}}_{\mathbf{ja}}$ | ${\mathit{P}}_{\mathbf{jb}}$ | ${\mathit{Q}}_{\mathbf{jb}}$ | ${\mathit{P}}_{\mathbf{jc}}$ | ${\mathit{Q}}_{\mathbf{jc}}$ |
---|---|---|---|---|---|---|---|---|---|---|

1 | 1 | 2 | 1 | 603 | 0 | 0 | 725 | 300 | 1100 | 600 |

2 | 2 | 3 | 2 | 776 | 480 | 220 | 720 | 600 | 1040 | 558 |

3 | 3 | 4 | 3 | 825 | 2250 | 1610 | 0 | 0 | 0 | 0 |

4 | 4 | 5 | 3 | 1182 | 700 | 225 | 0 | 0 | 996 | 765 |

5 | 5 | 6 | 4 | 350 | 0 | 0 | 820 | 700 | 1220 | 1050 |

6 | 2 | 7 | 5 | 691 | 2500 | 1200 | 0 | 0 | 0 | 0 |

7 | 7 | 8 | 6 | 539 | 0 | 0 | 960 | 540 | 0 | 0 |

8 | 8 | 9 | 6 | 225 | 0 | 0 | 0 | 0 | 2035 | 1104 |

9 | 9 | 10 | 6 | 1050 | 1519 | 1250 | 1259 | 1200 | 0 | 0 |

10 | 3 | 11 | 3 | 837 | 0 | 0 | 259 | 126 | 1486 | 1235 |

11 | 11 | 12 | 4 | 414 | 0 | 0 | 0 | 0 | 1924 | 1857 |

12 | 12 | 13 | 5 | 925 | 1670 | 486 | 0 | 0 | 726 | 509 |

13 | 6 | 14 | 4 | 386 | 0 | 0 | 850 | 752 | 1450 | 1100 |

14 | 14 | 15 | 2 | 401 | 486 | 235 | 887 | 722 | 0 | 0 |

Conductor | Impedance Matrix ($\mathsf{\Omega}/$mi) | ||
---|---|---|---|

$0.3686+j0.6852$ | $0.0169+j0.1515$ | $0.0155+j0.1098$ | |

1 | $0.0169+j0.1515$ | $0.3757+j0.6715$ | $0.0188+j0.2072$ |

$0.0155+j0.1098$ | $0.0188+j0.2072$ | $0.3723+j0.6782$ | |

$0.9775+j0.8717$ | $0.0167+j0.1697$ | $0.0152+j0.1264$ | |

2 | $0.0167+j0.1697$ | $0.9844+j0.8654$ | $0.0186+j0.2275$ |

$0.0152+j0.1264$ | $0.0186+j0.2275$ | $0.9810+j0.8648$ | |

$1.9280+j1.4194$ | $0.0161+j0.1183$ | $0.0161+j0.1183$ | |

3 | $0.0161+j0.1183$ | $1.9308+j1.4215$ | $0.0161+j0.1183$ |

$0.0161+j0.1183$ | $0.0161+j0.1183$ | $1.9337+j1.4236$ |

Line | Node i | Node j | Cond. | Length (ft) | ${\mathit{P}}_{\mathbf{ja}}$ | ${\mathit{Q}}_{\mathbf{ja}}$ | ${\mathit{P}}_{\mathbf{jb}}$ | ${\mathit{Q}}_{\mathbf{jb}}$ | ${\mathit{P}}_{\mathbf{jc}}$ | ${\mathit{Q}}_{\mathbf{jc}}$ |
---|---|---|---|---|---|---|---|---|---|---|

1 | 1 | 2 | 1 | 1850 | 140 | 70 | 140 | 70 | 350 | 175 |

2 | 2 | 3 | 2 | 960 | 0 | 0 | 0 | 0 | 0 | 0 |

3 | 3 | 24 | 4 | 400 | 0 | 0 | 0 | 0 | 0 | 0 |

4 | 3 | 27 | 3 | 360 | 0 | 0 | 0 | 0 | 85 | 40 |

5 | 3 | 4 | 2 | 1320 | 0 | 0 | 0 | 0 | 0 | 0 |

6 | 4 | 5 | 4 | 240 | 0 | 0 | 0 | 0 | 42 | 21 |

7 | 4 | 9 | 3 | 600 | 0 | 0 | 0 | 0 | 85 | 40 |

8 | 5 | 6 | 3 | 280 | 42 | 21 | 0 | 0 | 0 | 0 |

9 | 6 | 7 | 4 | 200 | 42 | 21 | 42 | 21 | 42 | 21 |

10 | 6 | 8 | 4 | 280 | 42 | 21 | 0 | 0 | 0 | 0 |

11 | 9 | 10 | 3 | 200 | 0 | 0 | 0 | 0 | 0 | 0 |

12 | 10 | 23 | 3 | 600 | 0 | 0 | 85 | 40 | 0 | 0 |

13 | 10 | 11 | 3 | 320 | 0 | 0 | 0 | 0 | 0 | 0 |

14 | 11 | 13 | 3 | 320 | 85 | 40 | 0 | 0 | 0 | 0 |

15 | 11 | 12 | 4 | 320 | 0 | 0 | 0 | 0 | 42 | 21 |

16 | 13 | 14 | 3 | 560 | 0 | 0 | 0 | 0 | 42 | 21 |

17 | 14 | 18 | 3 | 640 | 140 | 70 | 0 | 0 | 0 | 0 |

18 | 14 | 15 | 4 | 520 | 0 | 0 | 0 | 0 | 0 | 0 |

19 | 15 | 16 | 4 | 200 | 0 | 0 | 0 | 0 | 85 | 40 |

20 | 15 | 17 | 4 | 1280 | 0 | 0 | 42 | 21 | 0 | 0 |

21 | 18 | 19 | 3 | 400 | 126 | 62 | 0 | 0 | 0 | 0 |

22 | 19 | 20 | 3 | 400 | 0 | 0 | 0 | 0 | 0 | 0 |

23 | 20 | 22 | 3 | 400 | 0 | 0 | 0 | 0 | 42 | 21 |

24 | 20 | 21 | 4 | 200 | 0 | 0 | 0 | 0 | 85 | 40 |

25 | 24 | 26 | 4 | 320 | 8 | 4 | 85 | 40 | 0 | 0 |

26 | 24 | 25 | 4 | 240 | 0 | 0 | 0 | 0 | 85 | 40 |

27 | 27 | 28 | 3 | 520 | 0 | 0 | 0 | 0 | 0 | 0 |

28 | 28 | 29 | 4 | 80 | 17 | 8 | 21 | 10 | 0 | 0 |

29 | 28 | 31 | 3 | 800 | 0 | 0 | 0 | 0 | 85 | 40 |

30 | 29 | 30 | 4 | 520 | 85 | 40 | 0 | 0 | 0 | 0 |

31 | 31 | 34 | 4 | 920 | 0 | 0 | 0 | 0 | 0 | 0 |

32 | 31 | 32 | 3 | 600 | 0 | 0 | 0 | 0 | 0 | 0 |

33 | 32 | 33 | 4 | 280 | 0 | 0 | 42 | 21 | 0 | 0 |

34 | 34 | 36 | 4 | 760 | 0 | 0 | 42 | 21 | 0 | 0 |

35 | 34 | 35 | 4 | 120 | 0 | 0 | 140 | 70 | 21 | 10 |

Conductor | Impedance Matrix ($\mathsf{\Omega}/$mi) | ||
---|---|---|---|

$0.2926+j0.1973$ | $0.0673-j0.0368$ | $0.0337-j0.0417$ | |

1 | $0.0673-j0.0368$ | $0.2646+j0.1900$ | $0.0673-j0.0368$ |

$0.0337-j0.0417$ | $0.0673-j0.0368$ | $0.2926+j0.1973$ | |

$0.4751+j0.2973$ | $0.1629-j0.0326$ | $0.1234-j0.0607$ | |

2 | $0.1629-j0.0326$ | $0.4488+j0.2678$ | $0.1629-j0.0326$ |

$0.1234-j0.0607$ | $0.1629-j0.0326$ | $0.4751+j0.2973$ | |

$1.2936+j0.6713$ | $0.4871+j0.2111$ | $0.4585+j0.1521$ | |

3 | $0.4871+j0.2111$ | $1.3022+j0.6326$ | $0.4871+j0.2111$ |

$0.4585+j0.1521$ | $0.4871+j0.2111$ | $1.2936+j0.6713$ | |

$2.0952+j0.7758$ | $0.5204+j0.2738$ | $0.4926+j0.2123$ | |

4 | $0.5204+j0.2738$ | $2.1068+j0.7398$ | $0.5204+j0.2738$ |

$0.4926+j0.2123$ | $0.5204+j0.2738$ | $2.0952+j0.7758$ |

Method | Solution | Losses (kW) | Reduction (%) | Proc. Time (s) |
---|---|---|---|---|

Benchmark case | $\left\{1,1,1,1,1,1,1,1,1,1,1,1,1,1\right\}$ | 134.2472 | 0.00 | - |

Classical CBGA | $\left\{4,1,2,6,1,4,2,4,3,5,6,5,2,1\right\}$ | 109.2236 | 18.64 | 6.9435 |

BHO | $\left\{6,6,6,4,5,6,6,6,6,6,4,6,6,6\right\}$ | 110.0025 | 18.06 | 8.4850 |

SCA | $\left\{5,1,2,1,3,2,5,1,1,5,2,5,2,2\right\}$ | 109.3973 | 18.51 | 34.2865 |

VSA | $\left\{5,6,3,2,3,2,3,5,5,6,1,2,4,3\right\}$ | 109.3217 | 18.57 | 39.6831 |

Improved CBGA | $\left\{2,2,6,2,4,5,4,2,2,6,4,6,2,1\right\}$ | 109.1980 | 18.66 | 2.0762 |

Period | Act. (pu) | React. (pu) | Period | Act. (pu) | React. (pu) |
---|---|---|---|---|---|

1 | 0.1700 | 0.1477 | 25 | 0.4700 | 0.3382 |

2 | 0.1400 | 0.1119 | 26 | 0.4700 | 0.3614 |

3 | 0.1100 | 0.0982 | 27 | 0.4500 | 0.3877 |

4 | 0.1100 | 0.0833 | 28 | 0.4200 | 0.3434 |

5 | 0.1100 | 0.0739 | 29 | 0.4300 | 0.3771 |

6 | 0.1000 | 0.0827 | 30 | 0.4500 | 0.4269 |

7 | 0.0900 | 0.0831 | 31 | 0.4500 | 0.4224 |

8 | 0.0900 | 0.0637 | 32 | 0.4500 | 0.3647 |

9 | 0.0900 | 0.0702 | 33 | 0.4500 | 0.4226 |

10 | 0.1000 | 0.0875 | 34 | 0.4500 | 0.3081 |

11 | 0.1100 | 0.0728 | 35 | 0.4500 | 0.2994 |

12 | 0.1300 | 0.1214 | 36 | 0.4500 | 0.3336 |

13 | 0.1400 | 0.1231 | 37 | 0.4300 | 0.3543 |

14 | 0.1700 | 0.1390 | 38 | 0.4200 | 0.3399 |

15 | 0.2000 | 0.1410 | 39 | 0.4600 | 0.4234 |

16 | 0.2500 | 0.1998 | 40 | 0.5000 | 0.4061 |

17 | 0.3100 | 0.2497 | 41 | 0.4900 | 0.3820 |

18 | 0.3400 | 0.3224 | 42 | 0.4700 | 0.3820 |

19 | 0.3600 | 0.3263 | 43 | 0.4500 | 0.3887 |

20 | 0.3900 | 0.3661 | 44 | 0.4200 | 0.2751 |

21 | 0.4200 | 0.3585 | 45 | 0.3800 | 0.3383 |

22 | 0.4300 | 0.3316 | 46 | 0.3400 | 0.2355 |

23 | 0.4500 | 0.4187 | 47 | 0.2900 | 0.2301 |

24 | 0.4600 | 0.3652 | 48 | 0.2500 | 0.1818 |

Sol. Number | Solution | Losses (US$) |
---|---|---|

Benchmark case | $\left\{1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1,1\right\}$ | 43,226.9376 |

Solution 1 | $\left\{4,4,5,2,5,2,6,3,2,3,6,3,5,3,2,1,2,3,6,2,4,3,1,1,5,3,4,5,6,4,6,4,2,3,4\right\}$ | 35,105.2156 |

Solution 2 | $\left\{4,3,1,3,1,1,6,1,6,4,5,1,6,3,1,6,5,4,1,1,3,6,4,2,2,6,4,3,5,2,2,1,3,2,1\right\}$ | 35,127.0109 |

Solution 3 | $\left\{3,3,5,3,5,6,3,6,6,5,2,2,5,6,6,4,5,6,1,6,4,2,4,5,5,2,4,3,5,3,5,5,2,4,5\right\}$ | 35,133.1683 |

Solution 4 | $\left\{3,1,1,4,1,6,5,1,2,3,6,2,1,4,3,5,6,1,1,4,1,5,2,1,5,4,3,2,4,3,6,3,2,4,5\right\}$ | 35,137.2049 |

Solution 5 | $\left\{4,1,1,4,3,4,1,4,4,1,1,3,5,2,4,4,4,5,6,5,6,5,6,1,5,6,1,3,2,3,1,3,2,4,5\right\}$ | 35,140.4392 |

Solution 6 | $\left\{4,2,4,3,3,3,1,2,1,2,5,4,3,2,3,4,3,4,1,4,1,6,5,6,2,5,1,3,2,5,2,1,6,5,6\right\}$ | 35,154.0686 |

Solution 7 | $\left\{3,1,4,2,2,6,3,4,6,1,4,4,1,5,4,2,4,6,1,4,2,4,6,3,3,1,2,3,2,1,2,6,3,3,4\right\}$ | 35,157.9344 |

Solution 8 | $\left\{4,1,1,3,1,5,6,2,1,2,5,3,4,2,6,2,2,4,1,1,4,1,6,6,2,5,3,4,5,4,3,2,4,3,2\right\}$ | 35,169.6288 |

Solution 9 | $\left\{4,1,3,3,4,5,6,3,4,2,3,1,5,2,1,3,5,4,4,1,3,2,6,6,6,5,5,3,1,5,5,5,3,6,4\right\}$ | 35,179.4097 |

Solution 10 | $\left\{4,6,2,3,6,2,1,3,2,6,1,6,4,2,4,2,4,3,1,3,2,5,2,4,4,2,3,1,3,3,3,4,5,3,2\right\}$ | 35,180.3742 |

Publisher’s Note: MDPI stays neutral with regard to jurisdictional claims in published maps and institutional affiliations. |

© 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/).