Power Flow Solution in Bipolar DC Networks Considering a Neutral Wire and Unbalanced Loads: A Hyperbolic Approximation

: This paper addresses the problem of the power ﬂow analysis of bipolar direct current (DC) networks considering unbalanced loads and the effect of a neutral wire, which may be solidly grounded or non-grounded. The power ﬂow problem is formulated using the nodal admittance representation of the system and the hyperbolic relations between power loads and voltages in the demand nodes. Using Taylor series expansion with linear terms, a recursive power ﬂow method with quadratic convergence is proposed. The main advantage of the hyperbolic approximation in dealing with power ﬂow problems in DC bipolar networks is that this method can analyze radial and meshed conﬁgurations without any modiﬁcations to the power ﬂow formula. The numerical results in three test feeders composed of 4, 21, and 85 bus systems show the efﬁciency of the proposed power ﬂow method. All of the simulations were conducted in MATLAB for a comparison of the proposed approach with the well-established successive approximation method for power ﬂow studies in distribution networks.


Motivation
Recently, the interest in direct current distribution networks has increased due to their advantages in terms of efficiency and controllability and the advances in semiconductor technologies [1,2]. In addition, several distributed energy resources, such as storage devices, photovoltaic panels, and loads, operate directly in DC, reducing the presence of power electronics [3][4][5]. DC networks are a competitive alternative to AC transmission systems at both high and low voltage levels [6]. Systems at high voltage levels are called HVDC systems (i.e., high-voltage DC), and those at low voltage levels are called LVDC systems (i.e., low-voltage DC) [7]. These grids do not require synchronization and pose economic advantages due to their elimination of the adverse effects caused by reactive power requirements in AC networks [8]. Furthermore, it has been reported that DC networks have low power losses and improved voltage profiles [9].
Models, methods, and algorithms for analyzing conventional AC networks must be extended to DC networks [10,11]. Among these methods, the power flow is perhaps the most important, since it is part of every operation, control, or planning model of electric systems. Despite its apparent simplicity, a power flow for DC networks brings challenges, especially for bipolar networks, namely, networks that include positive, negative, and neutral wires [12]. It has been reported that bipolar LVDCs are more efficient and flexible than

Contributions and the Article's Outline
This article proposes a general model for power flow analysis in a bipolar DC grid. The contributions are twofold: • A generalized power flow approximation based on the Taylor expansion of the nonlinear relation between currents and voltages, which includes a neutral pole that is either solidly grounded or non-grounded. This approximation could be used for different studies, such as planning and operation, among others. • A set of numerical experiments in different benchmark test systems confirms that the hyperbolic approximation is applicable for meshed and radial grids, demonstrating the scalability of the model and the low number of iterations, even for large networks.
The remainder of the paper is organized as follows: Section 2 presents the formulation of the power flow in a bipolar LVDC system, and the approximation for the hyperbolic expressions is discussed to obtain a hyperbolic recursive formulation for solving the PF problem. Section 3 presents the test feeders that will be used in the present paper. Section 4 presents the computational validation of the proposed methodology; finally, Section 5 presents the conclusions and future work.

General Representation
A general representation of a bipolar DC network is given in Figure 1. It consists of a grid with a positive pole p, a negative pole n, and a neutral pole o that may be solidly grounded. In a given node k, the bipolar grid could contain multiple loads connected between positive and negative poles p dk,p−n or between positive and neutral poles p dk,po (analogously, a negative and neutral pole p dk,no ). Naturally, the system may be unbalanced. To deal with multiple constant power loads, we start with the analysis of the solidly grounded case and then extend the formulation to the non-grounded case. The bipolar DC grid is represented by an oriented graph F = {L, E }, where L = {1, 2, . . . , l} is the set of hypernodes and E ⊆ L × L is the set of hyperbranches. Each hypernode and hyperbranch has three components represented by the set B = {p, o, n}; these contain the positive, neutral, and negative poles. In this manuscript, the subscripts k and m are entries of L, and the electrical connections are represented by the nodal conductance matrix G bus . Capital letters represent vector variables or constants, while capital bold letters represent matrix variables or constants, and their entries are lower-case letters. The main variables, constants, and indices are presented in the nomenclature below. In the present work, we call the loads connected between a positive and negative pole a C-type connection, and we call the loads connected between a positive and neutral pole an E-type connection (see Figure 2).

Solidly Grounded Case
The voltage v o is zero for solidly grounded systems. Therefore, the general power flow equation for a multinodal bipolar DC can be represented by (1) (in the following equations, h, b, w are entries of the set B): where i h gkE is the current generation in the power source connected at node k with an E-type connection in the pole h; i h dkE is the current demanded at node k for a constantpower terminal with an E-type connection in the pole h; i h dkC is the current demanded at node k for a constant-power terminal with a C-type connection in the pole h; g hb kj is the component (entry) of the conductance matrix that relates nodes k and j and the poles h and b, respectively; v b k is the voltage at node k for the pole b. The definitions of the demanded current for the E-type and C-type connections are presented in Equations (2) and (3).
where p h dk is the power consumed at node k in the pole h with respect to the neutral wire; p h−w dk represents the power consumption at node k between poles h and w. Both (2) and (3) have hyperbolic form.
A linear approximation of (2) and (3) is proposed here. It consists of a Taylor expansion of the hyperbolic representations f (x) = x −1 and g(x, y) = (x − y) −1 around the operation point (x 0 , y 0 ), which is defined in (4) and (5).
Equations (4) and (5) are replaced in (2) and (3), so the equivalent currents for the E-type and C-type connections are found: We separate the general power flow Equation (1) into two terms: The first term is related to the slack node, and the second term is associated with the demand nodes. This separation is presented in Equations (8) and (9), respectively.
where i pn gE ∈ R 2×1 contains the power injections in the slack node for the positive and negative poles (pn); i pn dE ∈ R (2l−2)×1 and i pn dC ∈ R (2l−2)×1 are vectors that represent the demanded currents in the E-type and C-type connections. In this definition of the variables, l represents the number of nodes in the bipolar DC grid; v pn s ∈ R 2×1 and v pn d ∈ R (2l−2)×1 are vectors that contain the voltage output in the slack node and the voltage variables in the demand nodes for poles pn, respectively; g pn ss ∈ R 2×2 , g pn sd ∈ R 2×(2l−2) , g pn sd ∈ R (2l−2)×2 , and g pn dd ∈ R (2l−2)×(2l−2) are submatrices of the general conductance matrix g pn , which relates the slack and demanded nodes.
Let us define the following vectors and matrices: where P p−n d ∈ R (2l−2)×1 is the vector of bipolar constant-power loads, P pn d ∈ R (2l−2)×1 is a vector that contains monopolar constant-power terminals, and J ∈ R (2l−2)×(2l−2) and M ∈ R (2l−2)×(2l−2) are the identity matrix and an auxiliary matrix that allow the calculation of the voltage differences for bipolar constant-power terminals.
Equations (10) and (11) allow the generalization of the currents in (6) and (7), as presented below: where V pn,t d is a linearized vector with the voltage at the iteration t. If Equations (12) and (13) are replaced in Equation (9), then a recursive power flow formula is derived for the bipolar power flow problem with the neutral cable solidly grounded, as presented in Equation (14).
Finally, if we rearrange Equation (14) to obtain a general expression for V pn,t+1 d , the following general power flow formula is reached: where Note that the recursive power flow Formula (15) is evaluated from the initial point t = 0 by assigning the voltage at each node that is equal to the substation bus until the convergence criterion is met: where ε is the acceptable tolerance. Notice that, in this study, the metric for evaluating the convergence process is the voltage variation between iterations; however, some authors make use of the power loss formula as the metric. This requires two additional calculations, the first of which is that of the current across lines and the second of which is that of the power losses in each line; this is less direct than assuming the voltage variation as the metric of convergence [27].

Non-Grounded Neutral Wire
This subsection explores the extension of the power flow formulation for bipolar DC networks where the neutral wire is uniquely grounded in the substation bus. The main difference in this case is that the current per pole depends on the difference between the pole and load voltages. We define the vector of demanded current per node as i pon dk : where i p dk , i o dk , and i n dk are the net current injections at node k in the positive, neutral, and negative poles, respectively. These currents are calculated as follows: Equations (18)- (20) have the same hyperbolic relation defined for the bipolar current i h dkC in Equation (3). This implies that the currents for the positive, neutral, and negative poles can be approximated as presented in (7) by using the Taylor series expansion with its linear term.
To obtain the power flow formulation for the non-grounded neutral bipolar DC case, we rewrite the matrix Equation (9), as presented in Equation (21).
where i are the vectors of demanded currents and voltages in the demand nodes. V pon s ∈ R 3×1 is the voltage output in the slack node. G pon ds ∈ R 3×(3l−3) and G pon dd ∈ R (3l−3)×(3l−3) are submatrices associated with the bipolar conductance matrix considering the neutral wire. To define the general power flow formula for a neutral non-solidly grounded wire, the following vectors are defined.
in addition, we define the following auxiliary matrices: , were constructed by the authors in order to vectorize the recursive formula, and they are defined as follows: Now, the following approximation for the demanded current i pon d is found: With the approximation of the currents in Equation (24), it is possible to reach a general formula for bipolar DC grids by combining it with Equation (21), as presented in Equation (25).
Note that the recursive power flow Formula (25) is evaluated from the initial point t = 0 by assigning the voltage at each node that is equal to the substation bus until the convergence criterion is met, i.e.,

Test Feeders
This section presents three bipolar DC grids that were used in the numerical validations of the proposed hyperbolic power flow approximation. These systems were composed of 4, 21, and 85 buses. All of the information of these test feeders is presented below.

The Four-Bus System
The four-bus system is a radial bipolar DC network composed of four nodes and three lines; the slack bus is located at node 1, and it is operated at ±220 V in the positive and negative poles and at 0 V in the neutral pole. The demand information per pole in this system is presented in Table 1,

The 21-Bus System
The 21-bus grid is an adaptation of the monopolar DC network reported in [26]. The electrical configuration of this network is presented in Figure 3. The 21-bus grid has a radial configuration with a voltage output in the slack node of ±1 kV and with 0 V in the neutral point. The complete parametric information for this test feeder is listed in Table 2.

The 85-Bus System
This 85-bus grid with bipolar configuration corresponds to an adaptation of the IEEE 85-bus grid presented in [31]. For this test feeder, the substation bus is located at node 1 with a voltage output of ±11 kV for the positive and negative poles and 0 V in the neutral point. The electrical configuration of the 85-bus grid is depicted in Figure 4.
Information regarding the branch and load parameters for this test feeder is reported in Table 3.   Figure 4. Grid topology of the IEEE 85-bus system for bipolar power flow applications.

Computational Validation
The implementation of the proposed hyperbolic approximation for the power flow solution in DC bipolar grids with multiple unbalanced loads was executed in MATLAB 2021b on a PC with an AMD Ryzen 7 3700 2.3-GHz processor and 16.0 GB of RAM, running on a 64-bit version of Microsoft Windows 10 Single Language. The methodology was compared to the successive approximation method proposed in [26] with its corresponding bipolar adaptation. Furthermore, the authors of these articles reported that successive approximation is completely equivalent to the classic backward/forward algorithm in the calculation of the power flow model.

Results in the Four-Bus Grid
For this test feeder, we considered two simulation cases, i.e., the neutral wire solidly grounded in all of the buses and the neutral wire grounded only at the substation bus. Table 4 presents the voltage in the positive, neutral, and negative poles. Both cases were solved with the proposed hyperbolic power flow approximation; Table 5 presents the solutions obtained with the successive approximation power flow method.  The numerical results in Table 5 show that: (i) The successive approximation power flow approach and the hyperbolic power flow proposal reached the same numerical voltages with a tolerance assigned to 1 × 10 −10 and the same power losses for both simulation cases; however, the number of iterations of the hyperbolic approximation method was three for the solidly grounded case and eight for the non-solidly grounded scenario, while the successive approximation power flow method took seven iterations for the solidly grounded case and eight iterations for the non-grounded case. This behavior was attributed to the linear convergence of the successive approximation power flow method, while the hyperbolic approximation exhibited quadratic convergence, as demonstrated in [16]. (ii) As expected, the non-grounded neutral wire presented high power losses compared to those in the solidly grounded case. This situation occurred due to the currents that flowed through the neutral wire because of the presence of monopolar constant-power DC terminals. (iii) The voltage regulation for the solidly grounded neutral wire case was 3.6142%, and in the non-grounded scenario, it was 3.6210%. In addition, it was noted that, due to the load imbalances in the monopolar loads, the positive and negative voltage magnitudes per node differed among them.

Results in the 21-Bus Grid
Considering that the solution reached by the hyperbolic approximation took fewer iterations and ensured the same numerical results as those of the successive approximation power flow method with fewer iterations, here, we only present the results reached by our proposal. Table 6 presents the voltage profile performance for the 21-bus grid in the solidly grounded and non-grounded simulation scenarios.  The numerical results in Table 6 show that: (i) The numbers of iterations with the hyperbolic approximation were four for the solidly grounded case and 13 for the non-grounded scenario, while the successive approximation power flow method took 10 iterations in the solidly grounded case and 13 in the non-grounded scenario. (ii) The effect of the neutral wire grounded in all of the nodes of the network was very positive in terms of the amount power loss, since for this scenario, this was 91.2701 kW, while in the non-grounded case, the power loss increased until 95.4237 kW, i.e., 4.2136 kW of additional losses. (iii) The voltage regulation when the neutral wire was solidly grounded was defined by the positive pole at node 17, since it had a low voltage profile with a magnitude of 890.1027 pu, which defined a regulation of 10.9897% for the 21-bus system in this simulation case. In addition, when the neutral wire was non-grounded, the minimum voltage was maintained at node 17 with a magnitude of 888.2594 V for the positive pole, which produced a regulation in the whole system of 11.1741%, which was worse than that in the solidly grounded case; however, this situation was explained by the additional voltage droops in the neutral wire due to the absence of the grounded system at each node. (iv) In the non-grounded neutral wire scenario, the voltage at this conductor reached about 24.3408 V at node 17; even though this is a low voltage value, it can cause the misoperation of electronic devices in areas near this node, since power electronic devices can be sensitive to variations in the voltage magnitudes of their local references.
To demonstrate the ability of the proposed hyperbolic approximation to deal with the power flow problem in meshed DC bipolar configurations, here, we add two lines to the 21-bus system that allow the radial configuration to become a meshed one. The added distribution lines connect node 7 with node 19 and node 11 with node 16, with resistances per pole of 0.082 and 0.037 Ω, respectively.
The numerical simulations of the meshed configuration for the 21-bus system showed that: (i) The total power losses in this system with a solidly grounded neutral wire were 75.1832 kW, while for the non-grounded scenario, these were 78.7372 kW, i.e., as expected, there was an increment in the power losses of 3.5540 kW due to the currents that flowed through the neutral conductor. (ii) The performances of the voltage profiles in both simulation cases were very similar, as can be seen in Figure 5; however, there were some differences when comparing the voltage magnitudes between the positive and negative poles. For the solidly grounded simulation case, the magnitude of the voltage at node 17 for the positive pole decreased until 0.9252 pu, while for the negative pole, the minimum voltage occurred at node 18 with a magnitude of 0.9392 pu. It is important to mention that these differences between poles for this particular neutral wire connection are attributable to the imbalances caused by monopolar loads due to the positive pole having a total monopolar load of 554 kW, while the total load in the negative pole was 445 kW. For the meshed configuration, the total number of iterations of the proposed hyperbolic approximation method was four in both simulation cases; the successive approximation method took nine iterations in the solidly grounded case and 11 in the non-grounded neutral wire simulation scenario.

Results in the 85-Bus Grid
When the proposed hyperbolic power flow approximation was applied to the 85-bus system, the following results were obtained: (i) When the neutral wire was solidly grounded, the total power losses for this system were 452.2981 kW, i.e., 6.76% of the total power consumption; however, when the neutral wire was non-grounded, these power losses increased to 489.5759 kW, i.e., 7.32% of the net power consumption. Hence, there was an additional 37.2778 kW of power losses. (ii) The total numbers of iterations were about four for the solidly grounded case and 13 for the non-grounded case.
The voltage profiles in the 85-bus system for the positive and negative poles considering the solidly grounded and non-grounded operation cases are presented in Figure 6.  The voltage profiles in Figure 6 show that: (i) For the solidly grounded case, the operative costs of the minimum voltage in the positive pole were 0.9189 pu at node 54, and for the negative pole, they were 0.8950 pu at the same node. (ii) In the non-grounded scenario, these values were 0.9204 pu for the positive pole and 0.8925 pu for the negative pole.
Note that when the positive and negative voltage profiles are compared (see Figure 6a,b), it can be observed that the positive pole had a similar behavior to that of the negative pole, which was very similar to the reflection of the voltage profile with respect to the axis of the nodes; even if the system were to be perfectly balanced, then v p must be equal to the absolute value of v n .

Simulations Using the Quasi-Newton Formulation
To demonstrate the effectiveness of the proposed bipolar power flow formulation based on the hyperbolic approximation of the relation between the voltage and power relations in the power balance equation, here, we use an approximated solution method known as the quasi-Newton approach, which consists of calculating the Jacobian matrix only one time and keeping it constant during the entire iteration process. To implement the quasi-Newton method, the power flow Formulas (15) and (26) take the following structures: where The main advantage of keeping the Jacobian matrix constant for the solidly and nonsolidly grounded cases corresponds to the possibility of calculating it and storing it to be used in the recursive power flow formula without making inversions each time, which makes the processing times required to solve the power flow problem minor in comparison with the time taken by the original power flow formulation. One of the main advantages of keeping the Jacobian matrix constant for the cases studied here (solidly and non-solidly grounded cases) corresponds to the possibility of calculating it one time and storing it as a constant in order to be used in the recursive power flow formula without making inversions at each iteration; this approach reduces the processing times required to solve the power flow problem, but it increases the number of iterations required for the convergence of the algorithm.

Conclusions
A power flow model based on the Taylor expansion of the hyperbolic equation of the bus injection model in bipolar DC networks was addressed in this work. This formulation applies to radial and meshed grids. The method was compared to the well-known successive approximation algorithm, demonstrating fast convergence in terms of numbers of iterations and leading to the same point of operation (with an error of 1 × 10 −10 ); furthermore, the analysis of three benchmarks allowed a conclusion as to the scalability of our proposal, since it applies to systems with different numbers of nodes with reasonable execution times (and a low number of iterations). The numerical results demonstrated that it is desired to have a neutral wire that is solidly grounded in all the nodes, since the losses (and voltage imbalances) increase considerably due to the current in the neutral pole. However, this would be very expensive and impractical in real-life situations. Finally, it is essential to remark that the power flow is the basis for the development of algorithms and the study of other subjects, such as planning (which refers to studies such as the reconfiguration and expansion of a grid, the optimal selection, integration, and placement of distributed energy resources such as photovoltaic units or batteries, and other studies in the long term), operation (which is related to the optimal control of a grid), and protection (which is related to keeping a grid safe when there is a large disturbance on it) of bipolar DC networks.
In future work, it will be possible to develop the following: (i) the proposal of a sequential quadratic optimal power flow for bipolar DC grids with unbalanced loads based on the hyperbolic approximation presented in this research; (ii) the development of a mixed-integer conic formulation for solving the optimal pole-swapping problem in unbalanced bipolar DC networks; (iii) the practical implementation of a bipolar DC grid in order to test the proposed methodology for the power flow calculation; (iv) the detailed exploration of the effect of a non-solidly grounded neutral wire in the expected convergence of the proposed hyperbolic power flow formulation regarding the large variations in the number of iterations with respect to the solidly grounded scenario for the neutral wire.

Data Availability Statement:
The data used in this study are reported in the paper's figures and tables.

Conflicts of Interest:
The authors declare no conflict of interest.