Optimal Load Redistribution in Distribution Systems Using a Mixed-Integer Convex Model Based on Electrical Momentum

: This paper addresses the problem concerning the efﬁcient minimization of power losses in asymmetric distribution grids from the perspective of convex optimization. This research’s main objective is to propose an approximation optimization model to reduce the total power losses in a three-phase network using the concept of electrical momentum. To obtain a mixed-integer convex formulation, the voltage variables at each node are relaxed by assuming them to be equal to those at the substation bus. With this assumption, the power balance constraints are reduced to ﬂow restrictions, allowing us to formulate a set of linear rules. The objective function is formulated as a strictly convex objective function by applying the concept of average electrical momentum, by representing the current ﬂows in distribution lines as the active and reactive power variables. To solve the relaxed MIQC model, the GAMS software (Version 28.1.2) and its CPLEX, SBB, and XPRESS solvers are used. In order to validate the effectiveness of load redistribution in power loss minimization, the initial and ﬁnal grid conﬁgurations are tested with the triangular-based power ﬂow method for asymmetric distribution networks. Numerical results show that the proposed mixed-integer model allows for reductions of 24.34%, 18.64%, and 4.14% for the 8, 15-, and 25-node test feeders, respectively, in comparison with the benchmark case. The sine–cosine algorithm and the black hole optimization method are also used for comparison, demonstrating the efﬁciency of the MIQC approach in minimizing the expected grid power losses for three-phase unbalanced networks.


Introduction
Distribution systems in Colombia have voltage ranges from 1 kV to 57.5 kV [1], connecting transmission and sub-transmission networks with end users through substations in urban and rural areas, forming medium-and low-voltage networks [2].These are radially designed to reduce the costs of structural investment and protection systems, as well as to simplify coordination protection schemes, with the main characteristic that power flows are typically unidirectional (from the substation bus to the load nodes) [3,4].However, there is a trend towards an active distribution network connecting consumers and small generator units [5].
In distribution systems (DS), there are two types of electrical energy losses.Technical losses correspond to the energy that is dissipated and cannot be recovered in any way [6].Still, they can be reduced to acceptable values according to plans established by the distribution system's operators (DSOs) for this purpose and depend on the resistance of the conductors [7].Non-technical losses are due to electricity theft, represented by the users' completely or partially fraudulent consumption [6].This problem leads to high costs and decreased electrical energy efficiency and quality in the distribution network.
One challenge for DSOs is reducing electrical losses in distribution networks with the lowest possible investment costs.For example, in Colombia, DS losses range from 6% to 18% of the total losses [8].As an incentive to reduce losses in Colombia, recognition is given to DSOs that apply a plan to reduce losses according to an efficient loss limit regarding voltage levels [1].
In distribution networks, three-phase imbalances might happen for the reasons listed below: i.
The transposition condition does not apply given the short length of the lines and the distribution line designs are asymmetrical.ii.Line current and voltage imbalances are caused by the nature of the loads, which can be 1φ, 2φ, or 3φ, as well as by their uneven distribution and unpredictable behavior.iii.An imbalance in the current flowing through the lines is caused by the single-phase transformers' haphazard placement in the system's phases.iv.The use of new low-carbon emission technologies [7], which result in inefficient use of low-voltage network resources, increased energy losses, and the spread of imbalances to other connected nodes in the network [9,10].
The main cause of power losses is an increase in current through the neutral wire which results in anomalous load behavior or outages, lower task-related line security, and equipment overload in the power system [11][12][13].
The use of capacitors [14] has limitations when the three-phase imbalance is very severe.Consequently, reactive power compensation will not entirely solve the three-phase imbalance problem [15,16].Conductor classification is the use of mathematical models [17] or metaheuristic methods [18] to analyze the economics of improving conductors in specific areas of the network.In addition to reducing losses, this enhances the network's operational features.Feeder reconfiguration, as discussed in [19], is pursued for the hourly cost optimization of distribution systems, ranging from power loss costs to switching costs with an hourly analysis within a 24-h horizon, allowing for a network reconfiguration that adapts to the needs and operational conditions offered by various system actors, as well as to the economic and operating interaction of renewable generators, energy storage systems, and network mesh structures.Its key feature is the way in which relationships are presented, i.e., using a mixed integer non-linear programming (MINLP) problem, which is recast as a linear problem but contains solutions that are based on stochastic techniques.
On the other hand, a meta-heuristic methodology for network reconfiguration with demand variation is presented in [20,21], aiming to minimize network losses while observing the operating constraints of system voltages, currents, and the final topology.This methodology includes, among other things, particle swarm optimization and an improved coyote optimization algorithm.Network reconfiguration can be employed as a reaction to anomalous events (faults) in order to restore service to users not affected directly by the fault.The goal of implementing genetic algorithms is to produce results in shorter processing times while observing operating constraints.Important aspects in the assessment of various solution methodologies for minimizing losses in distribution networks through network reconfiguration, such as convergence times and compliance with operating constraints, are discussed in [22].Additionally, case examples for unbalanced three-phase distribution systems are described in [23].For optimization issues in these systems, exact approaches, linear and nonlinear programming models [24], and distributed generation [25,26] have also been used.However, implementing these methods in the distribution system is an expensive task.
In order to reduce distribution network losses, load redistribution optimization techniques must first identify the load of each phase [27] and calculate the network losses [8,28].These losses result in harmonics [29], overloading of the network in one or two specific phases, conductor overheating (which shortens the insulation life cycle), increased voltage drops that interfere with proper equipment operation due to increased current, and increased current through the neutral wire (which degrades protection performance).
In the scientific literature, the issues of ideal load redistribution, phase-balancing, or phase-swapping in three-phase asymmetric distribution networks have received considerable attention.Some of the relevant publications in this field are outlined below.
The authors of [30] present a solution methodology for the phase-balancing problem in low-voltage distribution networks by applying a hybrid optimization methodology that combines fuzzy logic with the power flow solution via the Newton-Raphson algorithm.This study's key contribution is that it suggests installing automatic load balancers in user nodes in order to enhance the voltage profiles as a function of the measured load changes.To solve the phase-balancing problem in three-phase asymmetric networks and reduce overall grid power losses, a methodology was developed in [31] which is based on the hybridization of the fuzzy logic and particle swarm optimization approaches.Excellent numerical results were obtained when this proposal was put into practice in a mediumlow voltage feeder in Bariloche, Argentina.However, the authors do not compare their approach against any other, making it a challenge to assess the efficacy and robustness of their methodology.The Internet of Things was used by the authors of [16] to provide an autonomous system for balancing loads in low-voltage networks.An effective genetic optimization algorithm was used as the load balancer.Using optimization techniques such as the multi-objective version of the genetic algorithm to reconfigure these loads in realtime applications, an experimental validation performed on a tiny three-phase platform demonstrates that it is possible to reduce energy losses and voltage imbalances.The work by [32] proposes the application of the multi-objective gray wolf optimizer to deal with the minimization of power losses and voltage imbalances in three-phase asymmetric networks by reconfiguring single-phase residential users while assuming the existence of smart meters.The main contribution of this research is an automatic grid load reconfiguration that considers the presence of multiple small photovoltaic generators along the feeder.Additional research articles regarding the optimal load reconfiguration problem in three-phase networks include the application of the Birkhoff polytope using group theory [11], artificial neural networks [4], mixed-integer convex approximations based on average powers and currents [3,33], the vortex search algorithm [7], and the sine-cosine algorithm [8], among others.
The main studies on the topic of using optimization approaches to solve the issue of optimal load redistribution in three-phase asymmetric distribution networks are shown in Table 1.
Table 1.Summary of the methodologies used for load redistribution in three-phase unbalanced networks.

Solution Methodology
Objective Function

Ref. Year
Fuzzy logic combined with the Newton-Raphson method Power losses minimization [30] 2009 Fuzzy evolutionary particle swarm optimization Power losses minimization [31] 2010 Quadratic optimization model Voltage imbalance minimization [34] 2017 Multi-objective genetic algorithm Voltage imbalance and switching actions minimization [16] 2019 Multi-objective gray wolf optimizer Voltage imbalance and power losses minimization [32] 2020 Birkhoff polytope Reducing the expected value of active power losses [11] 2021 Artificial neural networks and smart meters Energy losses minimization [4] 2021 Mixed-integer convex optimization approach based on the average current Power losses minimization [33] 2021 Mixed-integer convex optimization approach based on the average power Power losses minimization [3] 2021 Vortex search algorithm Power losses minimization [7] 2021 Sine-cosine algorithm Power losses minimization [8] 2021 The following are the key features of the optimization techniques in Table 1: (i) most of them rely on metaheuristic algorithms coupled with three-phase power flow techniques in a master-slave configuration; (ii) minimizing the anticipated grid power losses under peak load conditions is the main goal of most solution strategies; (iii) although the aim of these techniques is to minimize active losses, they also indirectly enhance the phase balance, minimize reactive losses, and result in smaller voltage drops and harmonics; and (iv) the best load balancing problem in asymmetric networks can have a number of potential solutions thanks to the use of metaheuristic optimization techniques, but the convexity and uniqueness of the solution cannot be ensured when using combinatorial optimizers.The exact method, in contrast, is known to produce a set of equations that are simple to use in any programming language and converge to a unique solution.
This study makes the following contributions based on the state of the art: It approximates the optimal load redistribution problem in three-phase asymmetric distribution networks using a mixed-integer convex (MIQC) optimization method.This approximation is based on the average electrical momentum.By analyzing the best load redistribution solution for a three-phase power flow problem, the suggested MIQC model's ability to reduce grid power losses can be validated.
Note that the objective of this research consists of presenting an approximated MIQC model to redistribute the load consumption per node by applying the concept of electrical momentum.Electrical momentum is formulated in analogy to the momentum generated by the force applied at a certain distance from the rotating point, where the distance is defined as the equivalent resistive effect and the force is equivalent to the total active and reactive power consumption.It is worth mentioning that, after reviewing the state of the art, this is the first time that an approximated optimization model based on the average electrical momentum is proposed to address the optimal load redistribution problem.In addition, this approach solves the optimal-phase balancing problem while considering the worst-case scenario, which corresponds to operation under peak load conditions.However, more research will be required to study networks with residential, industrial, and commercial users while considering daily load profiles, which can be a research opportunity for future work.
The remaining sections are organized as follows.The exact mathematical modeling of the ideal load redistribution problem is described in Section 2, where the proposed model's objective function and set of constraints are described.The solution methodology is described in Section 3, which is based on the solution of the suggested MIQC model in GAMS and MATLAB.Additionally, an iterative numerical approach is used to present the generic three-phase power flow solution.The model is evaluated in Section 4 via the IEEE 8-, 15-, and 25-node test systems.The complete results (including comparative solution reports), as well as their analysis and discussion, are provided in Section 5. Finally, conclusions and future work are presented in Section 6.

Mathematical Model
A mixed integer non-linear programming model can be used to simulate the optimal load redistribution problem in asymmetric distribution systems [33], which is given by the power flow formulation, due to the product that appears between the node voltage magnitudes and the trigonometric functions and the presence of binary variables associated with the set of connections for each load.[7].However, in this research, a mixed-integer convex programming model based on the electrical momentum is proposed (this is a convex approximation that neglects the effect of the voltage profile, i.e., it reduces the complexity of the impact of power balancing limitations).This is calculated by multiplying the equivalent resistance of each line by the active and reactive power accumulated for each node which affects the line upstream, representing the load's distance from the feeder.Below is a description of the simplified mixed-integer optimization model.

Objective Function
At each node over the feeder, the objective function is to reduce the electric momentum of the loads per phase.To this effect, the active and reactive power flow variables per line and phase are combined to generate a quadratic function, as seen in Equation ( 1) [3].
where z is the objective function to be evaluated, P l f and Q l f represent the active and reactive power in line l at phase f , and R ave l is the average value of the resistance of the three-phase line.L and F are the sets containing all the lines and phases in the network, respectively.
Remark 1.The main characteristic of the objective function is that it is a positive definite function (as R ave l is a positive parameter ∀l ∈ L multiplied by a sum of two square variables, i.e., P l f + Q 2 l f ), which implies that it is a strictly convex function that will ensure a global optimal solution with an efficient solution technique [35].Figure 1 illustrates the convex structure of the objective function in (1).
To illustrate the concept of electrical momentum and its relationship with power losses in an asymmetric three-phase imbalanced distribution grid, consider the electrical circuit in Figure 2. In order to keep things simple, it is assumed that the grid is entirely resistive, with loads represented by constant currents.It is assumed that the force applied at each node is defined as the current square demanded by it, and the distance is equivalent to the resistive value.Thus, the electrical momentum at node j can be defined as Considering the aforementioned definition, the total electrical momentum of the three-phase network in Figure 2 is defined below: Now, to observe the effect on an unbalanced three-phase demand current on the total electrical momentum defined above, consider the set of three possible load current consumptions listed in Table 2, where Option 0 is the original load connection, and Options 1 and 2 are two possible load redistributions.Electrical momentum The numerical results in Table 2 show that the electrical momentum for the tested load connections is different, which implies that redistributing the load consumption at each node directly affects the objective function.In addition, for the tested numerical example, it is observed that electrical momentum is equal to the function used to calculate power losses in a resistive grid, which confirms that both variables have an intrinsic relation.This relationship is exploited in this research in order to minimize power losses in three-phase asymmetric grids, as presented in the Section 5.

Set of Constraints
The set of constraints for optimal load redistribution is given by the binary nature of the decision variable x k f g and the accumulated power in each of the lines, among other constraints [3].The complete set of constraints is presented below.
where P d k f and Q d k f are the active and reactive power consumed after load redistribution at node k for phase f ; P d kg and Q d kg correspond to the active and reactive power consumed at node k in phase f before load redistribution; x k f g is the binary variable (matrix) that determines whether the load connected in phase g should be reassigned to phase f at node k, P r l f and Q r l f are the active and reactive power accumulated downstream, which have an impact on line l, and T lk corresponds to the position of the upper-triangular matrix that relates line l with node k (the general structure of the T matrix will be presented in the next section).Note that N is the set containing all nodes in the three-phase asymmetric network under analysis.
Remark 2. The main characteristic of the solution space implied by (2)-( 7) is that it defines a mixed-integer solution space, which means that, for each combination of the binary variables x k f g , it becomes an affine solution space that allows reaching the best possible solution via the branch-and-cut and interior-point methods with a logarithmic barrier.
The possible phase distributions of a three-phase system for the binary variable x k f g are presented in Table 3 (last column).From the first three connection types (1, 2, and 3), it is also evident that the sequence does not change, i.e., it remains positive.However, for the remaining connection types (4, 5, and 6), the sequence becomes negative [8].
Table 3. Possible load connections in a three-phase network.

Type of Connection
Phases Sequence For a system with n nodes, the dimension of the solution space is 6 n−1 , which means that it grows exponentially with the number of nodes.Therefore, efficient methods and models are required to solve this type of problem.When applying these types of connections, it is necessary to consider the sensitivity of the loads connected to the system in the face of operation in both positive and negative sequences.For phase exchange at each of the nodes, distribution systems must be automated.It should be taken into account that each user should have a smart meter that shares information with each local distributor in order to make the right decision at the right time.This, in addition to having switchgear that allows for the remote exchange of each of the phases in each of the transformers.

Solution Methodology
With the purpose of solving the mixed-integer quadratic convex (MIQC) optimization model to define the optimal redistribution of loads in asymmetric three-phase networks, as defined by optimization models (1)-( 7), a programming language capable of addressing convex models is required.For the proposed model, different solvers such as SBB, CPLEX, and XPRESS, which are available in the GAMS mathematical optimization software, were used to validate and solve the proposed optimization model.Afterwards, in order to evaluate the efficiency of the proposed model, the active power losses implied by the new load configuration were measured using the triangular power flow method in the MATLAB software.The main aspects of the proposed optimization methodology are presented below.

Three-Phase Power Flow
The unbalanced three-phase power flow problem is generally a complex numerical problem composed of nonlinear non-convex Equations [36], which includes a large number of unknown variables resulting from the structuring of a larger number of equations, given the representation of the system in the abc phase components.
The triangular power flow method [2] was selected for this study, in which the relationship between the nodal voltages on the distribution system and the three-phase currents through the branch is established.All this, through the implementation of the three-phase upper-triangular matrix, to finally relate the three-phase voltages and the apparent power at each demand node.The process is outlined below for the system shown in Figure 3. Figure 3 shows the current flowing through branch i.Here, J i can be defined as a function of the net current injection at each node k (i.e., I k ).The relationship between these variables is expressed in (8).Note that J i and I k are three-phase variables, i.e., they have ABC sub-indices.
where J 3ϕ ∈ C 3(b×1) is a vector that contains all of the three-phase branch currents in the complex domain (b branches), I 3ϕ ∈ C 3(n−1)×1 represents the vector containing all of the three-phase demanded currents at the nodes (n) of the network except the slack node, and T 3ϕ ∈ R 3b×3b represents the three-phase upper-triangular matrix.
The upper-triangular matrix can also be obtained by transposing the lower-triangular matrix obtained from the voltages.Consider that V k3ϕ represents the voltages at each of the demand nodes, V 1A , V 1B , and V 1C denote the voltage at the slack node, and E i3ϕ is the three-phase voltage drop in each of the branches.By applying Kirchhoff's second law for each path that contains the slack source, for example, V 1A , V 1B , and V 1C and all the nodes in the network, Equation ( 10) is obtained: 1 0 0 0 0 0 0 0 0 0 1 0 0 0 0 0 0 0 0 0 1 0 0 0 0 0 0 1 0 0 1 0 0 0 0 0 0 1 0 0 1 0 0 0 0 0 0 1 0 0 1 0 0 0 1 0 0 1 0 0 1 0 0 0 1 0 0 1 0 0 1 0 0 0 1 0 0 1 0 0 1 which is generally compared as defined below: where E 3ϕ ∈ C 3(b×1) is a vector that contains all of the three-phase voltage drops in each of the branches, V 3ϕ ∈ C 3(n−1)×1 represents the vector that holds the three-phase voltages at the demand nodes, and V 13ϕ is the three-phase voltage at the substation node.Note that 1 3ϕ is a vector of ones with appropriate dimensions.
After applying Ohm's law, the three-phase voltage drop in the branches can be expressed as a function of the three-phase branch currents.This is carried out using the three-phase primitive impedance matrix Z 3ϕ , which can be represented as shown in Equation (12) E 3ϕ = Z 3ϕ J 3ϕ (12) where Z 3ϕ ∈ C 3b×3b is a complex diagonal matrix that contains the three-phase branch impedances.For the formulation in Figure 3, the impedance matrix is represented as Z 3ϕ = diag Z 13ϕ , Z 23ϕ , Z 33ϕ .To obtain a general formula for the three-phase node voltages V 3ϕ and the three-phase node currents I 3ϕ , Equations ( 9) and ( 12) are substituted into (11), which yields (13) [37].
Note that (13) can be applied to a three-phase distribution network where the threephase current injection into a node can be calculated as presented below: where S 3ϕ ∈ C 3(n−1)x1 is a complex vector that contains all of the three-phase constant power demands and * represents the complex conjugate vector.Now, if ( 14) is substituted into (13), the following equation is obtained: Remark 3. Equation (15) relates the three-phase node voltages and the current injections at the demand nodes.However, note that the constant power demands must be expressed as current injections.Therefore, there is a need to obtain a recursive formula to determine the three-phase node voltages in terms of the constant power demands.Furthermore, this formula is only applicable to three-phase networks composed of Y loads.
Note that, depending on the type of star Y or delta ∆ connection of the demand, the demanded current I 3ϕ must be calculated differently [2].
Finally, based on (15), the recursive formula for solving the three-phase power flow via the triangular method is defined as Note that in (16), t represents the iterative counter.The convergence of this recursive formula can be found by applying the Banach fixed-point theorem [38], which proposes the existence of a fixed-point value whose function obtains the same value.Furthermore, this convergence is achieved when max V t+1 3ϕ ] − V t 3ϕ ≤ ε, where is the maximum tolerance error, typically assigned as 1 × 10 −10 .Finally, the following formula is used to determine the total grid power losses: where S 3ϕ is the sum of the apparent power losses for each phase over the total demand nodes in the system.Finally, to illustrate the general implementation of the three-phase power flow solution via the upper-triangular formulation, the pseudo-code in Algorithm 1 is present.
Algorithm 1: Generic solution of the power flow problem in three-phase unbalanced distribution networks using the triangular formulation Data: Determine the three-phase network under analysis 1 Obtain the per-unit equivalent for the system; 2 Generate the three-phase upper-triangular matrix T 3ϕ ; 3 Calculate the primitive branch impedance matrix Z 3ϕ ; 4 Define an auxiliary impedance matrix as Z bus 3ϕ = T T 3ϕ Z 3ϕ T 3ϕ ; Remark 4. Note that, prior to using the formula (16) in Algorithm 1, in the presence of loads, some additional calculations regarding the demanded current will be required [37].

Overall Solution Methodology
To illustrate the general implementation of the proposed solution methodology to solve the optimal load redistribution problem in three-phase unbalanced systems, the pseudo-code in Algorithm 2 is presented.
Algorithm 2: Optimal load redistribution in DS using a mixed-integer convex model based on electrical momentum Data: Define the distribution system under study 1 Solve the initial three-phase power flow problem using Algorithm 1; 2 Report the initial power losses; 3 Program the MIQC model based on the electrical moment with Equations ( 1)- (7) in the GAMS software; 4 Solve the MIQC model using a convex optimization tool; 5 Redistribute all of the demanded loads based on the solution provided by the MIQC solution; 6 Solve the final three-phase power flow problem using Algorithm 1; 7 Report the final power losses; 8 Compare the effect of the optimal load redistribution on the total grid power losses;

Test System
To evaluate the mixed-integer convex model based on electrical momentum for optimal load redistribution in distribution systems, three standard IEEE systems with 8, 15, and 25 buses were selected.The main characteristics of each one of these test feeders are described below.

8-Bus Test System
The unbalanced 8-bus system is a medium-voltage network consisting of eight nodes and seven lines, with a nominal line-to-line voltage of 11 kV.Additionally, the active and reactive powers consumed by phase a are 1005 kW and 485 kvar, respectively, those of phase b are 785 kW and 381 kvar, and those of phase c are 1696 kW and 821 kvar.These accumulated values correspond to node 1.The active power losses without system reconfiguration are 13.9925 kW (these values are obtained after evaluating the power flow for the initial grid conditions) [7].The electrical configuration of this system can be appreciated in Figure 4.All the parametric information regarding three-phase loads and branch impedances of the 8-bus system can be consulted in [7].

15-Bus Test System
The unbalanced 15-bus system is a medium-voltage network consisting of 15 nodes and 14 lines, with a nominal line-to-line voltage of 13.2 kV.Additionally, the active and reactive powers consumed by phase a are 9605 kW and 5226 kvar, those of phase b are 6480 kW and 4940 kvar, and those of phase c are 11977 kW and 8778 kvar.The active power losses without load reconfiguration are 134.2472kW.The electrical configuration of this system can be appreciated in Figure 5 [33].All the parametric information regarding three-phase loads and branch impedances of the 15-bus system can be consulted in [33].

25-Bus Test System
The 25-node unbalanced system is a medium-voltage network consisting of 25 nodes and 24 lines with a nominal line-to-line voltage of 4.16 kV.In addition, the active and reactive powers consumed by phase a are 946 kW and 648 kvar, those of phase b are 573.6 kW and 430.6 kvar, and those of phase c are 771.8kW and 554 kvar.The active power losses without reconfiguring all of the demand nodes in this system are 75.4207kW.The electrical configuration of this system can be appreciated in Figure 6 [2].All the parametric information regarding three-phase loads and branch impedances of the 15-bus system can be consulted in [7].

Results Analysis
The optimization model presented herein for optimal load redistribution in distribution systems, i.e., the MIQC model based on electrical momentum, was solved in the GAMS software using the SBB, CPLEX, and XPRESS solvers [39].In addition, the active power losses were evaluated before and after solving the MIQC model via the triangular power flow method, which was executed in MATLAB.Numerical results for the 8-, 15-, and 25-bus grids are presented below.

8-Bus System
Table 4 shows the active losses obtained through the solution of the electrical moment optimization method with the CPLEX, SBB, and XPRESS solvers, as well as the type of connection made.According to the results in Table 4, the best solution methods are CPLEX and XPRESS, with which a 24.34% loss reduction is obtained regarding the base case.On the other hand, Figure 7 presents the active and reactive power for the benchmark case and those obtained by the CPLEX solver.Note that the total phase imbalance regarding active power in the base case is 30.63%, and the imbalance value is reduced by 5.33% when the solution of the MIQC model is implemented (CPLEX solution).This implies a reduction of about 82.59% with respect to the benchmark case, demonstrating the positive effect of the optimal load redistribution in the expected grid imbalances.Figure 8 shows the magnitudes of the voltages for each phase in the nodes of the test system; Figure 8a shows the base case and Figure 8b shows the CPLEX solver.It can be observed that the behavior of the voltage profile improves after the optimization method is applied, and the imbalance is also reduced.According to Figure 8a, the largest deviation in the benchmark case occurs at node 4 in phase c, with 0.4081%.The deviation for the optimal solution reached with the CPLEX solver for the same node and phase was reduced to 0.0368%.It should be noted that the highest average deviation value for the solution with the CPLEX solver is 0.0781% and is located at node 8.For the same node, it is 0.2409% in the benchmark case.These results show the positive effect of the optimal load redistribution problem on the voltage profile performance of the network.Table 5 presents the numerical solutions regarding the active power losses for the phase balancing problem in the 8-bus system, as well as a comparison with metaheuristic approaches and exact methods found in the literature [33].Table 5 compares different literature reports regarding the optimal load balancing problem in three-phase asymmetric distribution networks for 8-bus systems.The numerical results show that (i) the proposed MIQC model based on the electrical momentum (MIQC-EM) reaches the best numerical solution reported for the 8-bus system, with a final power losses value of about 10.5869 kW-the same numerical solution found by the black-hole optimizer (BHO) and the sine-cosine algorithm (SCA)-and (ii) the previous MIQC model based on the average power (MIQC-AP) is stuck in a locally optimal solution, which is explained by the fact that the objective function did not include the effect of the average resistive value in all of the distribution branches, which is considered in our proposed solution strategy.

15-Bus System
Table 6 reports the results obtained by applying the MIQC-EM to the 15-bus system with three solvers available in the GAMS software: CPLEX, SBB, and XPRESS.According to the results presented in Table 6, the best solution is reached by CPLEX, with a 22.05% loss reduction compared to the base case.In addition, SBB and XPRESS are stuck in local optima.Figure 9 shows the active and reactive powers for the initial case and those obtained by the CPLEX solver.Note that the phase imbalance regarding active power in the benchmark case is 20.48%, which is reduced to 2.30% when the CPLEX solver is used, which corresponds to an improvement of about 88.76%, thus confirming the effectiveness of optimally redistributed loads in asymmetric distribution networks regarding the final load balance in the active and reactive power per phase.Figure 10 illustrates the positive effect of optimal load balancing on the voltage profile performance of the 15-bus grid.Figure 10 presents the magnitudes of the voltages for each phase in the nodes of the test system, i.e., in Figure 10a for the benchmark case and in Figure 10b for the CPLEX solver.It can be seen that the voltage profile behavior improves after the optimization method is applied, and the imbalance is also reduced.The largest deviation in the benchmark case occurs at node 14 in phase c, with 0.3482%.The deviation obtained by the CPLEX solver for this same node and phase was reduced to 0.0100%.It is worth noting that the highest average deviation value for the solution with the CPLEX solver is 0.0379% and is located at node 10.For this same node, it is 0.0445% in the base case, which demonstrates the positive effect of optimal load balancing on the expected phase voltage performance of the network.Table 7 presents the numerical solutions regarding the active power losses for the phase balancing problem in the 15-bus system, as well as a comparison with metaheuristic approaches and exact methods found in the literature [33].The numerical results in Table 7 show that (i) the proposed MIQC-EM finds a highquality solution for the optimal phase-balancing problem in the 15-bus grid, with a final value of 109.2218kW, i.e., a reduction of about 18.64% with respect to the benchmark case.The convex methods reported in [3,33] (MIQC-AP and MIQC-AI) also find adequate solutions regarding the final expected power losses value, which are better than the results reached with the BHO and the SCA.However, all of them are stuck in local optima.

25-Bus System
The solution of the proposed MIQC-EM-based model for the optimal load balancing in the 25-bus grid is reported in Table 8 for all the solvers available in the GAMS software, Version 28.1.2(CPLEX, SBB, and XPRESS).According to the results presented in Table 8, the best solver is XPRESS, which yields a 4.14% loss reduction compared to the base case.Note that, with respect to the solutions reached for the 8-and 15-bus grids, in this case, the XPRESS solver found a better numerical solution than CPLEX.This can be attributed to the increased solution space, which is highly dependent on the number of nodes under analysis.However, both numerical solutions are similar, with expected reductions of about 4.12% and 4.14% regarding the benchmark case, respectively.The active and reactive powers for the benchmark case and those obtained by XPRESS are shown in Figure 11.The phase imbalance regarding active power in the benchmark case is 16.60%, and the imbalance value is reduced by 0.97%, representing a 94.53% reduction when the proposed MIQC-EM model is solved via the XPRESS solver.Figure 12 shows the magnitudes of the voltages for each phase in the nodes of the 25-bus system, i.e., in Figure 12a for the benchmark case and in Figure 12b for the XPRESS solver.Note that for benchmark case (Figure 12a) exhibits undervoltage behavior in phase a when compared to phases b and c.However, this behavior is explained by the load conditions of phase a (it is the phase with the highest total load) in the benchmark case in comparison with the total load of the other two phases.Nevertheless, it can be seen that the behavior of the voltage profile improves after the optimization method is applied, and the imbalance is also reduced. of about 4.14% with respect to the benchmark case, which is second only to the MIQC-AI approach, with a reduction of about 4.16%, and (ii) all of the optimization methods exhibit reductions between 3.99% and 4.16% in the final objective function value with respect to the benchmark case.This implies that, even though the BHO, SCA, and MIQC-AP approaches are stuck in local optima solutions, they can be considered adequate for solving the studied problem.However, more research will be required to determine whether there are other solutions with better performance regarding the final objective function value in the IEEE 25-bus grid.

Conclusions and Future Works
This paper studied the optimal load redistribution problem for radial distribution systems.This problem was formulated using a mixed-integer quadratic convex model based on electrical momentum (MIQC-EM).This formulation uses a 3 × 3 binary matrix that relates the current and future state of the load in its corresponding phase as decision variables.The objective of the MIQC-EM model was to reduce the active power losses after load redistribution in the phases, for which a three-phase asymmetric power flow method was implemented, in order to evaluate the benchmark case and the solution reached with the MIQC-EM approach.Numerical results demonstrated that the method converges to high-quality solutions for test systems with 8, 15, and 25 nodes, which show reductions of 24.34%, 18.64%, and 4.14%, respectively.This article's main contribution was simplifying the optimization model of the optimal phase-balancing problem in three-phase asymmetric networks by considering a lower number of constraints and fewer variables.This was achieved via the MIQC-EM model.In addition, this study considers the system's unbalanced phases at each node and in the transmission lines of the network, as well as their impact on the system and its parameters, leading to more realistic representations of the analyzed system's behavior.Note that one of the main limitations of the proposed model is the need for programming languages that can solve MIQC models, as well as resources regarding time and decision variables required for the convergence of these solvers.
In this research, better results were obtained for the proposed method by using the CPLEX and XPRESS.Although exact methods are stable throughout the iterative process while using the same solver, the differences in their convergences and corresponding times are important aspects of the analysis.Additionally, the MIQC-EM model's effectiveness in solving the load-balancing problem was demonstrated, with the greatest loss reduction being reported in the 15-node system, as well as a global minimum for the 8-node system and a difference of 0.4808% in loss reduction when compared to the other MIQC method used in the IEEE 25-bus grid.
In future works, it may be possible to conduct the following studies: (i) the implementation of compensators or renewable energy sources in one, two, or three phases if required by the network for minimizing grid power losses and energy purchasing costs in the substation bus, and (ii) a three-phase network analysis in distribution systems with nodes with variable behavior over time within a 24-h horizon, along with the implementation of residential demand curves where the greatest number of single-phase and two-phase loads occur, for the optimal hourly distribution of loads.

Figure 1 .
Figure 1.Hyper-paraboloid that represents the objective function for the distribution line 1 in phase a, i.e., z 1 = R ave 1

Figure 2 .
Figure 2. Three-phase network for illustrating the concept of electrical momentum in the minimization of power losses.

4 Figure 3 .
Figure 3. Four-node test feeder to exemplify the relationship between net current injections I and branch currents J.

Figure 4 .
Figure 4. Single-line diagram of the 8-bus IEEE system.

Figure 5 .
Figure 5. Single-line diagram of the 15-bus IEEE system.

Figure 6 .
Figure 6.Single-line diagram of the 25-bus IEEE system.

Figure 7 .
Figure 7.Comparison of power demanded per phase for the benchmark case and the CPLEX solver in the 8-bus IEEE system: (a) active and (b) reactive power.

Figure 8 .
Figure 8. Voltage profiles for the IEEE 8-bus system: (a) benchmark case, and (b) best solution using the CPLEX solver.

Figure 9 .
Figure 9.Comparison of power demanded per phase for the benchmark case and the CPLEX solver in the 15-bus IEEE system: (a) active and (b) reactive power.

Figure 10 .
Figure 10.Voltage profiles for the IEEE 15-bus system: (a) benchmark case, and (b) best solution using the CPLEX solver.

Figure 11 .
Figure 11.Comparison of power demanded per phase for the benchmark case and the XPRESS solver the IEEE 25-bus system: (a) active and (b) reactive power.

Figure 12 .
Figure 12.Voltage profiles for the IEEE 25-bus system: (a) benchmark case and (b) best solution using the XPRESS solver.

Table 2 .
Numerical example for a three-phase unbalanced network.
I ja (A) I jb (A) I jc (A) I ja (A) I jb (A) I jc (A) I ja (A) I jb (A) I jc (A)

5
Select the maximum number of iterations t max and define the number of nodes as n = |N |; 6 Select the convergence error ; 8 Make t = 0; 9 Select the initial voltage as V t 3ϕ = 1 3ϕ V 13ϕ ; 10 k = 1;

11 for t ≤ t max do 12 for k = 2 : n do 13
Compute the demanded current at node k depending on the load type, i.e., Y or ∆;

Table 4 .
Optimization results found using different solvers with the proposed methodology for the 8-bus system.

Table 5 .
Numerical results for the 8-bus system.

Table 6 .
Optimization results found using different solvers with the proposed methodology for the 15-bus system.

Table 7 .
Numerical results for the 15-bus system.

Table 8 .
Optimization results found using different solvers with the proposed methodology in the 25-bus system.