A Comparative Study on Power Flow Methods for Direct-Current Networks Considering Processing Time and Numerical Convergence Errors

: This study analyzes the numerical convergence and processing time required by several classical and new solution methods proposed in the literature to solve the power-ﬂow problem (PF) in direct-current (DC) networks considering radial and mesh topologies. Three classical numerical methods were studied: Gauss–Jacobi, Gauss–Seidel, and Newton–Raphson. In addition, two unconventional methods were selected. They are iterative and allow solving the DC PF in radial and mesh conﬁgurations. The ﬁrst method uses a Taylor series expansion and a set of decoupling equations to linearize around the desired operating point. The second method manipulates the set of non-linear equations of the DC PF to transform it into a conventional ﬁxed-point form. Moreover, this method is used to develop a successive approximation methodology. For the particular case of radial topology, three methods based on triangular matrix formulation, graph theory, and scanning algorithms were analyzed. The main objective of this study was to identify the methods with the best performance in terms of quality of solution (i.e., numerical convergence) and processing time to solve the DC power ﬂow in mesh and radial distribution networks. We aimed at offering to the reader a set of PF methodologies to analyze electrical DC grids. The PF performance of the analyzed solution methods was evaluated through six test feeders; all of them were employed in prior studies for the same application. The simulation results show the adequate performance of the power-ﬂow methods reviewed in this study, and they permit the selection of the best solution method for radial and mesh structures.


General Context
Power-flow analysis is an essential tool in the design and operation of alternating current (AC) and direct-current (DC) networks [1][2][3]. It employs loads that are typically modeled as constant power loads (CPLs), which generate hyperbolic non-convexities on the mathematical model that represents the power-flow problem (PF) in electrical grids [4,5]. The focus of the power-flow analysis is to find the voltage profiles of all the buses that compose the grid under steady-state conditions [6]. These profiles in conjunction with the parameters of the grid allow the identification of all technical indexes that represent the behavior of the electrical network, such as the power loss, voltage stability, power flow on the branches, and power factor [7]. The PF has been widely studied in the past decades on AC networks [8][9][10]. However, concerning power-flow analysis in DC grids, the exploration of this problem started recently due to the paradigm shift regarding energy distribution with renewable energy and batteries, which typically operate with DC technologies [11][12][13][14][15][16]. These technologies are currently under development. In addition, DC grids present different advantages with respect to AC grids. For instance, DC grids do not make use of synchronizing generators, need fewer power converters, and most importantly, present a simpler implementation given that reactive power and frequency analyses are not required [17]. This has promoted the implementation of DC networks in the last few years. Finally, it is important to mention that for planning and management strategies of DC electrical networks, it is necessary to execute hundreds or even thousands of power-flow analysis to select a solution that provides the best results in terms of solution quality and processing time [18]. For this reason, it is absolutely necessary to make an adequate selection of the method to solve the PF for electrical grids with both radial and mesh topologies.

Literature Review
Owing to the mathematical and computational complexity that represents solving the DC PF, specialized software was employed in the past decades, e.g., Digsilent [19], Matpower [20], and ETAP [21]. The problems associated with the specialized software are high-cost acquisition, high implementation complexity, and the processing time required to solve the DC PF. Because of these drawbacks, classical solution methods were adapted in the last decade to solve the PF in DC networks without requiring specialized software. Examples are the Gauss-Jacobi, Gauss-Seidel, and Newton-Raphson methods [22,23]. The results reported in previous studies demonstrated the effectiveness of using classical numerical methods to solve the DC PF. In this study, the convergence of the Gauss-Seidel method when solving the DC PF was demonstrated on the basis of a fixed-point theorem described in [22]. The convergence of the Newton-Raphson method was proved by employing the Kantorovitchs theorem reported in [23]. Power-flow solution methods using semidefinite programming [24] and second-order cone programming [25] were also reported in the literature. These solution methods transform the PF in a convex equivalent problem using positive semidefinite matrices and implementing interior-point methods to solve the problem. The main drawback of these solution methods is that they increase the number of variables of the problem, which affect the computational efforts required [26]. Furthermore, these methods need complex solvers to solve the equivalent mathematical model. In AC grids, this problem has been widely studied during the last decades [27,28]. This is why many solution methods have been reported in previous studies. Given that the PF in DC grids is currently under development, different authors proposed solution methods based on prior studies on AC grids, such as linear approximations [29], successive approximations [30], triangular matrix formulation [31], graph theory, and backward/forward methods [18,32]. They established the improvement of convergence and the processing time required by the solution method as the performance criteria. Furthermore, the authors used different radial and mesh test systems, such as 10-, 21-, 33-, and 69-bus test systems [31,33], with the aim of evaluating the effectiveness and robustness of the proposed solution methods. Finally, it is important to stress that most methods proposed in the last few years to solve the DC PF are based on sequential programming to eliminate the use of specialized software, owing to the aforementioned drawbacks.

Motivation
A careful analysis of the mathematical representations of the AC and DC PFs [1,23,34] shows that the mathematical representation of the DC PF is completely different from that of the AC PF. Consequently, the solution methods proposed for the AC problem do not guarantee an adequate solution for the DC PF. Because of this, along with the paradigm shift resulting from the integration of DC grids on electrical distribution systems, and the important advances in DC distribution networks and their advantages with respect to AC networks, efficient and novel methods were proposed in the last decade to analyze DC grids in terms of power flow. These methods adapted non-linear solution methods proposed for AC grids given that the DC-grid power flow is, in fact, non-linear and non-convex, which hinders obtaining analytical solutions [23,35]. Currently, different methods allow solving the DC PF in DC grids with radial and mesh structures based on classical and unconventional numeral methods. However, it is necessary to identify and select methods based on sequential programming to eliminate the need for specialized software. To this end, test scenarios allow evaluating the performance of these methods in terms of solution quality and processing time. The objective is to select the methods that achieve the best effectiveness and robustness to solve the DC PF for both radial and mesh DC grids.

Contributions and Scope
The study of the PF in DC grids is critical. It is necessary to identify the solution method that attains the best performance in terms of solution quality and computational effort for both radial and mesh structures. In this study, we reviewed, analyzed, and evaluated eight solution methods based on three classical methods (Newton-Raphson, Gauss-Jacobi, and Gauss-Seidel) [23], and five unconventional methods based on the Taylor series expansion and successive approximations [30], triangular matrix formulation [31], graph theory, and backward/forward methods [18,32], to solve the PF in DC grids. Six different test systems were used. The main objective was to select the best solution methods in terms of convergence and processing time for mesh and radial electrical structures. Furthermore, it is important to stress that all the solution methods selected are based on sequential programming, thereby eliminating the need for specialized software and reducing the complexity and implementation cost. The main contributions of this study can be summarized as follows: Identify and select solution methods proposed in the literature based on sequential programming, thereby eliminating the need for specialized software, such as Digsilent, Matpower, and ETAP. This reduces the cost and complexity of solving the DC PF. Presentation of each power flow approach in a tutorial manner. This allows studying DC networks from beginners (engineering students) to experts (power-system researchers) regarding numerical methods and algorithms. Identify through simulation results the solution methods with the best PF performance in DC grids with radial and mesh structures, which are essential in terms of numerical convergence and computational efforts for planning and management strategies of DC grids.

Paper Organization
This paper is organized as follows: Section 2 presents the mathematical formulation of the PF considering both resistive and constant power loads. Section 3 presents the solution methods selected in this study to solve PFs in DC networks with radial and mesh structures other than the classical Gauss-Jacobi, Gauss-Seidel, and Newton-Raphson methods. Section 4 describes the test feeders, comparison techniques, and parameters used to analyze the results. Section 5 presents numerical results. Section 6 discusses the conclusions and possible future work lines resulting from this study.

DC Power-Flow Mathematical Formulation
The mathematical formulation associated with the DC PF is composed of a set of non-linear and non-convex equations that allow finding all electrical variables (voltage and currents) that represent the behavior of an electrical grid in steady-state conditions [36][37][38]. The equations that compose the DC power-flow formulation are obtained by using Kirchhoff's laws and Tellegen's theorem in electrical networks that satisfy the following conditions [22]: (i) the DC grid contains at least one constant voltage bus (slack bus); (ii) isolated buses do not exist across the grid; (iii) all voltage buses are contained inside the voltage bounds V min ≤ v ≤ V max , being V min > 0; and (iv) external perturbations are not considered inside the DC grid operation and steady-state conditions are applied.
According to these conditions and the aforementioned mathematical laws and theorems, the DC PF can be mathematically described as follows [39]: where v ∈ R n , i ∈ R n , and p ∈ R n represent the bus voltages, current, and power injections in all the buses of the network; D (v) ∈ R n×n is a diagonal positive definite matrix such that D ii = v i ; i = 1, 2, ..., n and D ij = 0; i = j; finally, G ∈ R n×n is the nodal admittance matrix. This matrix results from the addition of the conductance matrix, which represents the conductivity effects of each branch that interconnects the different buses of the electrical system (G B ∈ R n×n ), and the conductance matrix associated with the resistive loads connected to the DC grid (G n ∈ R n×n ). In addition, G can be divided in four components, as shown in Equation (2), where G gg and G dd represent the conductance matrix associated with the generators and load connections, respectively. Note that G gd = G dg T denotes the conductance matrix that relates the generators and loads. Finally, in this mathematical formulation, n and b are the number of buses and branches, respectively. The set of equations presented in Equation (1) can be rewritten as follows [39]: In the previous equations, i g ∈ R s , i d ∈ R n−s , v g ∈ R s , and v d ∈ R n−s are associated with the currents and voltages at the ideal generators (known voltage buses) and demand buses (unknown voltages), respectively; p g ∈ R s denotes the power generated by the ideal voltage buses; and p d ∈ R n−s represents the power generated or demanded by the demand buses. This occurs because the demand buses are composed of loads and distributed generators, different from the slack buses (ideal generators). Note that D v g ∈ R s×s and D (v d ) ∈ R (n−s)×(n−s) can be interpreted as D (v). Finally, s corresponds to the number of slack buses; in this way, s ≥ 1, and the power generated by the distributed generators and the energy storage systems, located in demand buses, are considered as negative loads in the mathematical model.
Note that Equation (5) corresponds to a set of linear expressions given that the variables of the problem are p g and v d , where p g is a free variable that guarantees global power balance. Notwithstanding, Equation (6) corresponds to a non-linear expression resulting from the product between variables (v d ). In this way, Equation (6) represents the PF for DC power grids in both mesh and radial structures [40], and Equation (5) is not required. Note that the resolution of the PF through numerical methods is carried out only when the grid has at least one constant power terminal, that is, constant power load/consumption. This condition transforms the mathematical formulation of the PF into a non-linear and non-convex problem with hyperbolic structure, which makes numerical methods necessary [1,41].

Solution Methods
This section describes five solution methods reported in the last few years to solve the PF in DC grids [23]. Two methods are based on the Taylor series expansion and successive approximations [30], which solve the PF in radial or mesh structures by leveraging the mathematical formulation in Equation (6). Furthermore, three additional methods were selected to solve the PF in DC grids with radial topologies. These methods are based on triangular matrix formulation [31], graph theory, and backward/forward methods [18,32], respectively. All these solution methods were published in a period no longer than two years, and they were selected first because of their excellent results in terms of convergence and processing time reported by the authors, and also because they are based on sequential programming, thereby eliminating the need of specialized software.
To validate the performance of all the solution methods analyzed in this study in terms of both solution quality and computational effort, three classical numerical methods proposed in the literature were used for comparison: the first is the Newton-Raphson method (NR) [42], which finds the voltage profiles by using an iterative process based on a Jacobian matrix. The second and third numerical methods are the Gauss-Jacobi (GJ) and Gauss-Seidel (GS) approaches [43], which follow a simple iterative procedure to find the voltage profiles. These solution methods are not explained in this paper because they are classical methods typically described in textbooks of power analysis [43].

Power-Flow Solution Methods for DC Networks with Radial or Mesh Structures
This subsection addresses the PF for DC grids with radial and mesh structures through two different iterative techniques. The first technique applies a linearization based on a series expansion with the Taylor series [29], whereas the second employs the original set of non-linear equations that represent the problem [30]. Both solution methods are explained in the following subsections.

Taylor-Series-Based Approximation (TBM)
To obtain a linear approximation of the mathematical model that represents the problem discussed here, an equivalent linear expression for Equation (6) is required. Given that the inverse value of the voltage profile at the kth bus is a non-linear term, it is used as the alternative form of Equation (1) expressed as follows: The linear approximation of the term 1 v k presented in Equation (7) can be obtained from the first-order Taylor series expansion around the operating point v 0 k [44]; v k is limited, as described in the mathematical formulation. The mathematical expression that represents the Taylor series expansion used here is presented in Equation (8).
Then, applying Equation (8) to the non-linear term 1 v k around v 0 k , it is possible to obtain its equivalent linear form: By substituting Equation (9) into Equation (7), we obtain the linear approximation of the PF in DC grids by using the Taylor series expansion: Given that the main objective of the PF is to find all bus voltage profiles throughout the electrical network, the linear approximation of the PF does not consider the voltages associated with the generators with voltage control capacity (slack buses). This is because these voltage profiles are known [40]. In light of the above analysis, the net injected current, power, and voltage profile in the kth bus are equal to the current (i k = i d ), voltage (v k = v d ), and power (p k = −p d ) on the demand buses. In this manner, the linear approximation of Equation (10) is rewritten as follows: where p d is a vector that contains all constant power loads, and v 0 d denotes the initial voltage profiles composed of a vector filled by those with dimensions b × 1 (plane voltage profiles). Furthermore, by using the diagonal operator property, Equation (11) can be rewritten as follows: By replacing Equation (4) into Equation (12), we obtain the linear Equation (13), which allows the calculation of vd (voltage profile in the demand buses). Note that this equation solves the PF in DC grids with radial and mesh structures.
Finally, aiming at reducing the convergence error associated with the voltage profiles, we propose an iterative version of Equation (13) in Equation (14), which is obtained by replacing v 0 d with v t d (i.e., the value of v d in iteration t) and adding some terms.
3.1.2. Successive Approximation (SA) For the mathematical formulation of the SA, we use Equation (6) and a simple mathematical process to obtain the following expression [22]: Given that G dd is a positive definite symmetric matrix, its inverse always exists. Based on this condition, the following expression allows solving the PF using fixed-point theorems: Finally, an iterative procedure to calculate v d is proposed in the following equation, for which the iterative counter t was added to Equation (16): Finally, for both solution methods previously described, in this study we propose an iterative process that is executed until an acceptable convergence error or a maximum number of iterations (t max ) is reached (see Algorithm 1). To this end, we use the variation in the voltage profile as convergence error during two consecutive iterations: where corresponds to the minimum convergence rate fixed.
Algorithm 1: Iterative power-flow method for direct-current (DC) resistive networks with radial or mesh structures.
Data: Define DC radial test system Construct the conductance matrix G; Define v t d (with t = 0), , and t max ; for t = 0 : t max do Evaluate v t+1 d (i.e., apply Equation (13) or (17));

Power-Flow Solution Methods for DC Grids with Radial Structure
To describe a DC grid with a radial structure mathematically, it is necessary to consider two additional assumptions: A DC grid contains a radial structure if it has n buses and b branches, and b = n − 1. This guarantees a unique route between each pair of buses, as shown in Figure 1. In radial grids, a unique ideal generator with voltage control capability is assumed. This enables global power balance and voltage control on the DC grid.
An example of electrical network with a radial structure is shown in Figure 1. This network presents a typical radial DC grid composed by a unique main generator (an ESS formed by multiples batteries in this example), branches, resistive loads, and CPLs. The electrical grid or distributed generations (DGs) can replace the batteries, provided that these devices have sufficient power to achieve the global power balance on the DC network. Next, we present the three numerical methods used for solving the PF in DC grids with radial structure:

Power-Flow Method Based on Triangular Matrix (TM) Formulation
This subsection describes an iterative sweep method based on a TM formulation [31]. It uses a primitive impedance matrix and the relationship between the bus and branch currents of the DC grid. The main advantage of this method is that it does not require the inversion of non-diagonal matrices. This reduces the mathematical and computational efforts.
For a general DC grid, the primitive resistance matrix is a positive definite matrix (R p ∈ R b×b ), which is defined in Equation (18), where R ij denotes the resistive parameter of the branch that connects buses i and j, respectively. The branch currents i B ∈ R b×1 are associated with the bus-injected currents i d ∈ R (b)×1 through a triangular incidence matrix T ∈ R b×(n−1) , as shown in Equation (19). In this equation, i d denotes the net current demanded at the load buses, which is a function of the power demanded (loads) or generated (distributed generators) by CPLs and resistive loads; i d is calculated through Equation (20).
(3) (4) To carry out the construction of the triangular incidence matrix T, it is first necessary to order the nodes of the DC grids, as reported in [45]. Then, to calculate T, we use the following conditions:

•
If branch k is located upstream bus l, the component T(k, l) takes a value of 1.

•
If branch k is located downstream bus l, the component T(k, l) takes a value of 0.
As an illustrative example of the structure of T, consider the radial DC electrical network depicted in Figure 1. This network contains six demand buses, a voltage-controlled bus 0, and six branches represented with letters from a to f , respectively. Then, by using Equation (19) and the conditions described in the previous paragraph, T for the DC grid in Figure 1 takes the following form: After obtaining T, it is possible to calculate the voltage drop associated with each branch that composes the electrical grid. For this, we use the relation between the voltage and currents in the branches reported in Equation (22), where ∆v B ∈ R b×1 represents the vector that contains the voltage drops in all branches of the DC network.
Subsequently, by applying second Kirchhoff's law for each closed loop between the slack bus (bus 0 in Figure 1) and the demand buses, it is possible to obtain Equation (23) to calculate V d , where v 0 denotes the voltage at the slack bus and 1 ∈ R (n−1)×1 is formed by a vector filled by ones.
Finally, by combining Equations (19), (20), and (22) with Equation (23), we obtain the following expression: where the matrix operation produces a constant value that depends on the bus connections and branch parameters: T T R p T. This term can be defined as the equivalent impedance matrix (R bus ) and can be integrated into Equation (24). The mathematical formulation that represents the power-flow method based on TM formulation to solve the DC PF by considering radial structures is expressed as follows: Given that this equation corresponds to a non-linear and non-convex formulation, it is necessary to add an iterative counter t to Equation (25); see the resulting Equation (26). To solve the DC PF, the iterative procedure presented in Algorithm 2 must be applied.
Algorithm 2: Proposed iterative power-flow method for DC resistive networks based on triangular matrix (TM) formulation. Data: Define DC radial test system Apply the nodal ordering method [45]; Construct the triangular incidence matrix T; Construct the conductance matrix G L ; Construct the primitive resistance matrix R p ; Build and store the R bus matrix; Define v 0 , v t d (with t = 0), , and t max ; for t = 0 : t max do Evaluate v t+1 d with (13) ;

Sweep Method Based on Graph Theory (SMBGT)
The Sweep Method Based on Graph Theory (SMBGT) correlates the branch and nodal currents through an incidence matrix by using an iterative procedure that provides simplicity in its mathematical formulation and shorter convergence time. The basic concept of the branch-to-bus incidence matrix is described below.
The incidence matrix A for radial electrical networks is always square and invertible. This matrix is constructed as follows: • A(k, l) = 1 if branch k is connected to the bus l and the current leaves bus l. • A(k, l) = −1 if branch k is connected to bus l and the current is arriving to bus l. • A(k, l) = 0 if branch k is not connected to bus l.
Moreover, in radial distribution networks, the current can be assumed to flow from the substation to the load buses, that is, the direction of the current through branch k, connected between buses i and j, is from i to j if bus i is upstream of bus j; otherwise, the current flows from bus j to bus i.
To provide a numerical example of the branch-to-bus incidence matrix, consider the radial DC grid depicted in Figure 1, which is represented by incidence matrix A in the following expression: The rows of A are branches (a,b,c,d,e,f), and the columns represent buses (1, 2, 3, 4, 5, 6). The slack bus (bus 0 in Figure 1) is not contained in the incidence matrix; instead, it is represented in an additional matrix A 0 . The voltage drop in each branch of the network is calculated by Equation (28), where v s represents the voltage on the main generator (voltage-controlled bus).
The relation between branch currents is also established: Moreover, the relation between the voltage drops and currents in all the branches is obtained by applying Ohm's law: where R p is the primitive resistance matrix described in the previous subsection. In the example depicted in Figure 1, R p takes the following form: Combining Equations (28) to (30) and rearranging some terms, we obtain the following expression: The formula defined in Equation (31) relates the bus voltages and currents; nevertheless, an iterative procedure is needed to solve this set of equations, given that i d depends on bus voltages v d , as shown in Equation (20). Finally, by replacing Equation (20) in Equation (31), and adding an iterative counter t, we obtain Equation (32), which allows solving the PF in DC grids by using the SMBGT.
The formulation of the power-flow method based on SMBGT requires a recursive actualization of the variables; this process is described in Algorithm 3.

Algorithm 3:
Iterative procedure for power-flow analysis in radial DC networks through the Sweep Method Based on Graph Theory (SMBGT).
Data: Define DC radial test system Apply the nodal ordering method; Construct A and A 0 ; Obtain the inverse of A; Construct the primitive resistance matrix R p ; Define v s , v 0 d , , and t max ; for t = 0 : t max do Evaluate v t+1 d in (32) ;

Backward/Forward Sweep Method (BF)
The mathematical formulation that represents the BF is enabled by the relation between the bus voltage and branch currents [18]. This method uses an iterative procedure based on Equations (29) and (30) to solve the PF in DC grids with a radial structure.
The BF starts by defining the voltages in all the nodes as follows: In Equation (33), t represents the iterative counter, and the condition t = 0 corresponds to the initial solution, denoted by v 0 d : plane voltages in per-unit representation. Then, by applying Kirchhoff's laws for each closed loop between the slack bus and other buses, and by using Equations (29) and (30), it is possible to find the recursive equation that allows calculating the voltage profiles in the buses other than the voltage-controlled bus: Finally, Algorithm 4 describes the iterative process required by the BF, starting from the initial solution, where v 0 d = 1 → is the vector filled by ones with dimensions b × 1. Define v 0 , v 0 d , , and t max ; for t = 0 : t max do Calculate the nodal currents (i d ) by using Equation (20); Calculate the branch currents (i B ) by using Equation (29)

Computational Analysis
The next subsection presents the test feeders used to validate the effectiveness in terms of convergence and processing time of the solution methods employed in this review paper, and the parameters applied to carry out a fair analysis of the results obtained.

Test Systems
In this study, six test systems were used to evaluate the performance of the solution methods selected, including four radial test systems with 10, 21, 33, and 69 buses. The test systems composed of 10 and 69 buses were modified by adding multiple branches to obtain its equivalent mesh structure. The test systems with 10 and 21 buses were taken from the literature [22,33]. The test systems with 33 and 69 buses are modifications of test systems used in AC networks [46][47][48]. The reactance in the lines and the reactive power consumed by the loads were eliminated to obtain the DC equivalent systems. The electrical configuration and data corresponding to the test systems are described below.

10-bus Test System
The traditional electrical configuration associated with the 10-bus test system is presented in Figure 2. This electrical system is formed by 10 buses and nine branches [22], with a unique generator and radial structure, and including resistive and constant power loads. The technical characteristics of the 10 bus-test system with a radial structure are presented in [22]. In this study, a mesh topology of this test system is proposed by adding two lines to the radial version following the connections and parameters given in Table 1. This modified test system is identified in this study as a 10-bus mesh test system. In Figure 2, the branches added to obtain the mesh test systems are in red. Table 1. Proposed connections for testing the mesh grid case in the 10-bus test system. from to R (p.u) from to R (p.u) 5 10 0.0035 8 10 0.0025

21-Bus Test System
The 21-bus test system was reported in [33]. This electrical DC network is formed by 21 buses and 20 branches. Its main characteristic is that it only considers constant power loads in the system. The electrical configuration of this test system is shown in Figure 3. Its base voltage and base power are 1 kV and 100 kW, respectively. The voltage at the slack node (Node 1) is considered flat, that is, 1 p.u. The last condition was applied to all test systems used in this study.

33-Bus Test System
This test system is formed by 33 buses and 32 branches, presenting a unique slack bus and multiple buses with constant power loads. The 33-bus AC electrical system was proposed in [49], and its electrical configuration is shown in Figure 4. This section neglects the values corresponding to the reactance and reactive power demanded in the AC system to obtain the DC equivalent network with 33 buses. For this test system, base values of voltage (12.66 kV) and power (100 kVA) were used to obtain p.u. values.

69-Bus Test System
This test system is obtained in the same manner as the 33-bus test system previously described, i.e., making some modifications in the electrical configuration, which is shown in Figure 5. The original 69-bus AC test system was proposed in [50] The base values used in this test system to obtain the per-unit values were 12.66 kV and 100 kVA, respectively. Four branches were added to the traditional 69-bus test system [32] to obtain an equivalent mesh system. The data of these branches are described in Table 2 and presented in red in Figure 5.

Simulation Results
This section compares the results obtained with the numerical methods selected. The first scenario compares all the solution methods that allow the PF to be solved in DC networks with a mesh structure. The second scenario evaluates the same problem in grids with a radial structure. The main objective of these scenarios is to find the methods with the best performance to solve the DC PF in each type of structure. The classical techniques (GJ, GS, and NR) together with Taylor-Series-Based Approximation (TBM) and SA were considered in both scenarios because these methods can solve the DC PF in mesh and radial networks.
For the sake of a fair comparison between different solution methods, in all scenarios we set = 1 × 10 −10 , t max = 1 × 10 3 iterations, and 1 × 10 5 consecutive executions to obtain the average time required. These values were obtained in a heuristic form. The main motivation to carry out 1 × 10 5 consecutive executions for each solution method is to eliminate the effect that any other process can cause to the simulation processing time. Furthermore, in this study, we used the NR as the base case to analyze the performance of the solution methods in terms of quality solution and processing time. In [23], this method converges to the solution of the problem discussed here. The adopted criteria to compare all the solution methods are the average voltage error (VoltageError), power loss error (P loss Error), and average processing time. Finally, all simulations of this chapter were carried out on a Dell Precision T7600 Workstation with 32 GB of RAM and an Intel(R) Xeon(R) CPU ES-2670 ()at 2.50 GHz, using MATLAB software.

Simulation Results for DC MGs with Mesh Structure
In this subsection, we analyze the results obtained by the NR, GJ, GS, TBM, and SA methods when they were applied to the mesh test systems. Table 3 shows the numerical results; from left to right: method, VoltageError, P loss Error, and average processing time required by each method. With respect to NR, it is appreciated from the analysis of VoltageError and P loss Error for the solution methods that these error values can be neglected for both mesh test systems given their small values. Concerning VoltageError, the worst case was presented by GJ, with an error of (1.93 × 10 −7 ), whereas the lower error was presented by TBM (2.38 × 10 −13 ). Furthermore, regarding the power loss error, the maximum error was achieved by GJ (1.59 × 10 −5 ), whereas the minimum error was obtained by TBM (1.89 × 10 −12 ). Based on these values, it is possible to conclude that all the solution methods evaluated in this review paper are adequate in terms of convergence to solve the DC PF in grids with a mesh structure.
According to the analysis presented in the previous paragraph, the selection of the method with the best performance solving the DC PF in grids with mesh structure is based only on the average processing time. Figure 6 depicts the processing time required by each method with regard to NR, such that the time needed for the NR method corresponds to 100%. Note that the worst methods in terms of processing time are GJ and GS, which presented an average increase of 51,247% and 28,154% with respect to NR, respectively. The other two methods provide a considerable reduction in the average processing time, with TMB requiring 43% and SA requiring only 18% of the time used by the NR method. SA exhibits a shorter processing time in all test systems; hence, SA is the solution method, among the studied ones, with the best performance solving the PF for DC grids with a mesh structure. Finally, note that the variation in processing time for each method was less than (1 × 10 −5 ) in the different test scenarios proposed. This demonstrates that 1 × 10 −5 consecutive executions eliminated the effect of other processes in terms of processing time.

Simulation Results for DC MGs with Radial Structure
In this subsection, we analyze the results provided by the methods used to solve the PF in DC grids with radial structure; these results are presented in Table 4. Note from the values of VoltageError and P loss Error in the four radial test systems that the same situation as in the case of mesh structure occurs: the errors are small enough to be neglected. The minimum and maximum values of VoltageError were 3.88 × 10 −15 and 8.86 × 10 −8 , respectively, and the minimum and maximum values of P loss Error were 1.77 × 10 −12 and 2.12 × 10 −5 , respectively. Therefore, all the considered methods provide solutions with similar quality.  Figure 7 presents the processing time of the methods in percentage with respect to NR. Note that GJ and GS provide an average increase of 40,933% and 32,804% with regard to NR, respectively. However, the other methods provide an average reduction with respect to such base case: 82% (TBM), 22% (TM), 32% (SA), 24% (SMBGT), and 11% (BF). By analyzing these results, we conclude that BF achieves the best results in terms of processing time, with an average reduction of 84.65%. Finally, note that the variation in processing time for each method was less than (1 × 10 −5 ) for each solution method used in the different test scenarios proposed for DC grids with radial structures.

Conclusions
This study addressed the calculation of the power flow in DC grids considering both radial and mesh structures to select the solution methods with the best performance in terms of numerical convergence and processing time for mesh and radial DC grids. To understand this problem, Section 2 presents the mathematical formulation, where the assumptions and equations that describe the problem are completely discussed.
The review of the state-of-the-art presented in this study stemmed from the need to identify efficient methods to solve the DC PF providing adequate performance in terms of solution quality and processing time. These methods are based on sequential programming to eliminate the need for specialized software, thereby reducing the implementation complexity and costs. For these reasons, three classical and five unconventional methods were selected (they are explained in Section 3): the classical solution methods Gauss-Jacobi, Gauss-Seidel, and Newton-Raphson; two methods based on Taylor series expansion and successive approximations, both applicable on radial and mesh grids with multiple ideal buses; and three methods based on triangular matrix formulation, graph theory, and backward/forward methods for DC grids with radial structures of any size.
To verify the performance of the selected methods, four radial test systems composed of 10, 21, 33, and 69 buses were employed. Test systems with 10 and 21 buses were taken from [22,33]. The 33and 69-bus test systems are modified versions of systems used in AC networks [46,48]. Furthermore, two mesh adaptations of the 10-and 69-bus test systems were defined by adding branches to the radial versions. The performance of the methods was evaluated in terms of solution quality and processing time by using the three aforementioned classical methods reported in the literature as comparison methods [22,43]. The comparison was based on the average voltage error, power loss error, and average processing time. Note that, to calculate the average processing time, 100.000 executions were run. This allowed obtaining differences in processing time below (1 × 10 −5 ) for each execution of the solution methods on the proposed test scenarios. This demonstrates that it is possible to reduce the effect of other computational processes in terms of processing time. In addition, NR was selected as the base case for comparison because its convergence was demonstrated in [23]. This entails that it converges in DC grids with any structure.
In addition, this study compared two test scenarios: the first scenario adopts grids with mesh structure only by using NR, GJ, GS, TBM, and SA. The second scenario evaluates the performance of all the methods in networks with a radial structure. It is essential to stress that the average voltage and power loss errors were neglected in both scenarios owing to the lower values achieved by the solution methods: a maximum average voltage and power loss error of 1.93 × 10 −7 and 2.2 × 10 −5 , respectively. For this reason, we concluded that all the studied methods are adequate to solve the DC PF in grids with mesh and radial structure in terms of solution quality (convergence).
According to the previous results, the selection of the methods with the best performance was based on the average processing time. The simulation results obtained for the mesh test systems demonstrated the robustness and effectiveness of the proposed methods in comparison with the other solutions reported in the literature. From these results, it was possible to identify SA as the power flow method with the best performance for mesh DC grids: it achieved a reduction in processing time equal to 82% with respect to the base case (NR). For the second scenario, (DC grids with radial structures), the sweep iterative methods provided the best results as the size of the DC network increased. BF was the solution method that achieved the best performance, with a reduction in processing time equal to 89% with respect to the NR. Overall, according to the study conducted, the implementation of BF is recommended to solve the PF in DC grids with radial structure, especially when they are analyzed in the context of large systems. Furthermore, the SA performed the first in mesh grids (82%) and the second in radial networks (76%) given that it is not possible to evaluate the performance of the BF in mesh grids. Thus, the best solution method to solve the DC PF for any network size and for mesh and radial grids is the SA method.
As future work, it is possible to improve the efficiency of the mathematical methods for those approaches requiring matrix inversions. In addition, parallel processing tools can be used to enhance the effectiveness of the solution methods in terms of processing time. Furthermore, theoretical complexity analysis for the power-flow solution methods could be considered and then compared with simulation to improve the computational analysis of these methods. Finally, the application of SA and BF for planning and control strategies of DC grids could be studied. The objective would be to reduce the processing time, thereby improving the exploration of these algorithms and their performance in terms of quality of solutions.