Stability of DC Distribution Systems: An Algebraic Derivation

Instability caused by low inertia and constant power loads is a major challenge of DC distribution grids. Previous research uses oversimplified models or does not provide general rules for stability. Therefore, a comprehensive approach to analyze the stability of DC distribution systems is desired. This paper presents a method to algebraically analyze the stability of any DC distribution system through the eigenvalues of its state-space matrices. Furthermore, using this method, requirements are found for the stability of three example systems. Additionally, a sensitivity analysis is performed, which considers the effect of several system parameters on the stability and disputes some overgeneralized conclusions of previous research. Moreover, various simulations are performed to illustrate the dynamic behavior of a stable and an unstable DC distribution system.

A challenge for future distribution grids is instability, which is caused by the (anticipated) low inertia of future grids and the interaction between power electronic converters. Conventional generation provides (virtual) inertia to the grid by naturally governing the kinetic energy of its angular mass or via the control of its power electronic converter. With the increasing share of renewable (distributed) energy sources, this (virtual) inertia of distribution grids is decreased significantly [14]. Moreover, tightly regulated load converters behave as constant power loads (CPLs). Constant power loads exhibit a negative incremental input impedance, which has a destabilizing effect on distribution systems [15].
Two current methods have been identified to counteract instability caused by CPLs in DC distribution systems. First, passive stabilization methods utilize passive elements (e.g., resistance, filters or storage) to dampen the oscillations in the system [16,17]. However, passive stabilization methods increase the cost, complexity and losses of the system. Second, active stabilization methods utilize advanced control methods to improve stability. Examples of different active stabilization techniques can be found in literature [18,19]. However, it is more appealing to ensure inherent stability (if possible).
To analyze the stability of DC distribution systems, several models of the network and converters have been reported in literature. Most methods average and linearize the power electronic converters in the DC distribution grid. Since these small-signal models are valid for frequencies up to half of the switching frequencies, and the bandwidths of practical converters are typically lower than one-tenth of the switching frequency, these models provide reasonable accuracy around the operating point [18].
Four main approaches to the stability analysis of DC distribution grids can be identified. Firstly, the minor loop gain (the relationship between the complex load and source impedance) can be analyzed. Different limits on the minor loop gain have been proposed to ensure stability. However, this approach assumes unidirectional power flow and measurements are essential for accurate impedance estimations [18,[20][21][22]. Secondly, a root locus analysis of the system can be performed and the locations of the poles can be analyzed [23][24][25][26][27]. This approach does not provide rules for stability and often provides little (mathematical) insight into the origin of the (in)stability. Thirdly, an analysis based on Lyapunov methods can be conducted [28][29][30][31][32]. However, these are not easily applied and require a suitable construction of the Lyapunov storage function. Lastly, the roots of the system can be derived from the eigenvalues of its state-space matrix [19,[33][34][35][36]. This approach, which is also used in this paper, relies on the validity of averaging, linearizing, and simplifying the power electronic converter.
Previous research, however, only analyzes specific systems, uses oversimplified models and/or does not provide general rules for stability. For example, only star type systems with a source at the central node are analyzed [35], the node capacitance is not considered in the equations [37], or no (general) rules for stability are provided [25]. Therefore, any conclusions on the effect of system parameters (e.g., inductance, capacitance and droop coefficients) on the distribution system's stability can not be generalized.
The contribution of this paper is a generalized method to algebraically analyze the stability of DC distribution systems, regardless of configuration. The method derives necessary and sufficient conditions for the stability by determining the system's eigenvalues algebraically. Furthermore, the presented model is used to prove and refute some of the common conclusions on the effect of certain system parameters on stability.
The remainder of this paper is organized as follows: in Section 2, the model which is used for each component is described. Section 3 discusses the stability of the simplest DC distribution systems. A method to analyze the stability of any DC distribution grid is presented in Section 4. Subsequently, in Section 5, the presented method is applied to three different examples. Section 6 illustrates the stable and unstable dynamic behavior of DC distribution grids via various simulations. Finally, conclusions are drawn in Section 7.

DC Distribution System Model
In Figure 1, an example of a DC distribution system is shown. Any DC distribution system consists of n nodes, l distribution lines, o phase conductors, and m sources and loads connected to nodes via power electronic converters. For the sake of simplicity, the DC distribution systems in this paper are assumed to be monopolar (i.e., have a single phase conductor), but the model can be extended. The interconnectivity of the system can then be described by the incidence matrix Γ [38]: where j and i are the indices for each distribution line and node, respectively, and I j is the current flowing in line j. Furthermore, the boldface of a variable indicates that it is a vector or matrix.

Nodes
The voltages at the system's nodes are determined by the set of differential equations where C is the (diagonal) capacitance matrix, U N are the voltages at each node, I L are the currents in each distribution line, and I M are the sums of the currents from the converters connected to each node.

Distribution Lines
The currents in the system's distribution lines are determined by the set of differential equations where L and R are the (diagonal) inductance and resistance matrices, respectively.

Sources
It is assumed that the source and load converters have a large enough bandwidth that they react instantaneously to changes in the system. The equivalent circuits of droop controlled sources and constant power loads are given in Figure 2. The current flowing from a droop controlled source into a node is then where U 0 is the reference voltage, U n is the node voltage, and Z s is the (virtual) impedance of the droop controller.

Loads
The linearized load current of constant power loads, assuming a large enough bandwidth of the converter, is where P l is the load's power, U l is the voltage at which the load is linearized, and Z l is the (negative) impedance of the load [26].

Stability of Simple DC Distribution Systems
The simplest (potentially unstable) DC distribution grid is a voltage source connected to a constant power load via a distribution line, which is shown in Figure 3. The differential equations of this system are where U C is the voltage of the capacitor and I L is the current flowing in the inductor. By substituting Equation (5) into the differential equations, the state-space formulation is derived to be The poles of any state-space system can be found by determining the eigenvalues of the "A" matrix (the left matrix in Equation (9)). In fact, the eigenvalues of this matrix can be found by solving |A − λI| = 0. The characteristic equation (|A − λI|) of this system is Therefore, the poles of this system are Any system is considered stable if all poles have negative real parts. Consequently, this system is stable if and only if A slightly more complex DC distribution grid is given in Figure 4, where the voltage source of the previous example is replaced by a droop controlled source. To generalize droop controlled sources and constant power loads, the current flowing from any converter into the node it is connected to is given by where I i,0 and Z i are positive for droop controlled sources, and negative for constant power loads. The state-space formulation of this DC distribution system is then Once again, the characteristic equation is derived by determining |A − λI|: It is required, but not sufficient, that all coefficients of a polynomial are positive for it to have only poles with negative real parts; therefore, To make the conditions sufficient, in this third order case, one additional condition (that the product of the second and third coefficients is greater than the fourth) is required.
If, for the sake of simplicity, the resistance (R) is neglected and the source is located at node 1, the system presented in Figure 4 is stable if and only if

Stability of Any DC Distribution System
The approach used in the previous section is in this section extended for any DC distribution system. By combining Equations (2), (3) and (14) the state-space formulation for any DC distribution consisting of nodes, distribution lines, droop controlled sources and (constant power) loads is derived to be where Z is the matrix containing source and load impedance and I M,0 is the vector of the constant term currents (I i,0 ) from the converters connected to each node.

The Location of the Equilibrium
Besides the stability, the equilibrium (steady state) of the system can be derived from the state-space equations. In steady state, the time derivatives in the system are 0. The steady state voltages and currents can therefore be calculated by The inverse of the left-hand matrix (A) can be decomposed as By substituting Equation (27) into Equation (26), the steady state node voltages and line currents are derived to be Since the matrix of Equation (28) is a block matrix, the node voltages at the equilibrium are which is equivalent to multiplying the constant term currents by the equivalent impedance of the network.

Stability from the Characteristic Equation
The stability of DC distribution systems can be evaluated by determining the eigenvalues of the A matrix. For the ease of computation, this matrix is represented by where F and G are not necessarily square matrices. The characteristic equation of this matrix (|A − λI|) will always have the form of where N and L are the number of nodes and number of distribution lines, respectively. Therefore, this characteristic equation has N + L coefficients and N + L zeros. The coefficients of the characteristic equation can be determined from the A matrix by utilizing traces of powers or the principal minors of the matrix [39,40]: where Tr is the trace of a matrix, and ∆ k is the k-th order principal minor of A.
Utilizing Equations (33) and (30), the algebraic representation of the first five coefficients of the characteristic equation are determined to be From Equations (35) to (39), it becomes apparent that the coefficients relate to the combinations of the sources', loads' and distribution lines' time coefficients, of increasing order with each subsequent coefficient (without creating loops).
For these systems to be stable, it is required, but not sufficient, that all coefficients are larger than 0 [41]. To make the conditions sufficient (N + L − 2) additional coefficients need to be added [42]. These Routh coefficients can be found in the leftmost column of the Routh array [42]: where x is equal to N + L, and the b i,j elements are recursively determined by For example, the first two Routh coefficients are Any DC distribution system is stable if and only if all the coefficients of the characteristic equation and the Routh coefficients (b 1,j ) are larger than zero.

Stability Analysis of Example Systems
To demonstrate the utility of the presented method, the stability of three DC distribution systems in different configurations is analyzed. For these examples, it is assumed that the source is located at node N 1 of each configuration; however, moving the source to a different node yields similar results. Furthermore, a sensitivity analysis is performed and a few general misconceptions concerning the stability of DC distribution systems are discussed.

Bus Configuration
An example of a DC distribution system in a bus configuration is shown in Figure 5. The idea of such a configuration is one "bus", or set of distribution lines, without branches or meshes. For the shown example, it is assumed that a source is situated at N 1 and a load is situated at both N 2 and N 3 . However, the derivation for other configurations is analogous. For the sake of simplicity, the resistance of the distribution lines is neglected. Using Equations (30) and (33), the coefficients of the characteristic equation are derived to be Now if, for simplicity, we assume that the capacitance in each node (C i ) and inductance in each distribution line (L j ) are equal, we arrive at two simple requirements for the stability of this system where the first requirement is derived from the even coefficients, and the second requirement follows from the odd coefficients. Interestingly, from these requirements, a constraint on the capacitance and inductance follows: which defines a minimum on the ratio between the capacitance and inductance for feasible stability. The derived requirements are necessary but not sufficient for stability. To make the set of requirements sufficient additional constraints can be derived from the Routh coefficients. In general, the requirements derived from the Routh coefficients are more complex, but similar to the requirements derived from the characteristic equation's coefficients (see (23)). For this example, according to Equations (40) and (41), if the ratio between the capacitance and inductance is large enough, it is sufficient if which, agreeably, is a stricter version of Equation (50) derived from the characteristic equation's coefficients.

Ring Configuration
The second example is a DC distribution system in a "ring" configuration shown in Figure 6. In this case, the ring configuration does not have any branches, but it is meshed. In this example, the source is located at N 1 and loads are located at N 2 and N 3 . Analogous to the bus configuration, the requirements for stability of this system can be derived from the characteristic equation. These requirements are derived to be (54)

Star Configuration
The last example is a DC distribution system in a "star" configuration, which is shown in Figure 7. In this case, the star configuration has branches, but no meshes. Since, in most cases, the source is situated in the center node, and it is assumed that the system has a source placed at N 1 and a load is placed at N 2 , N 3 and N 4 . Analogous to the previous configurations, the requirements for stability of this system are derived to be The stability requirements for the three examples show strong congruence with respect to the sensitivity towards, for example, inductance and capacitance. Therefore, it becomes viable to make some general conclusions.

Sensitivity Analysis
Similar results with respect to previous work are found for some parameters. With respect to inductance, the results of this paper are congruent with previous works [26,37], where a decrease in inductance leads to increased stability. Furthermore, from the equations in the previous section, it can be observed that the impedance (slope) of the droop controlled source has an upper and lower bound. The upper bound is related to the constant power loads in the system, while the lower bound is related to the oscillations in the system. These oscillations originate from the interaction of capacitance at different nodes through the lines' inductance. These results are also congruent with previous work [43,44].
Generally, it is thought that an increase in capacitance is beneficial for the stability of the system. From Equations (52), (54) and (55), it is clear that there is a minimum required capacitance and further increase of capacitance leads to an increased stability margin. However, Equations (45) and (49) suggest a negative effect on the overall damping of the system. The second (a 2 ) and last (a N+L ) coefficients are the sum and product of all eigenvalues, respectively. Therefore, when the capacitance is increased the sum and product of all eigenvalues are decreased leading to decreased damping in the system. Conceptually, this can also be understood, a voltage deviation with larger capacitance means a larger disparity in energy, which will take a long(er) time to be damped.
Interestingly, Equations (45) and (49) suggest that decreasing the source capacitance compared to the load capacitance has a positive influence on the damping in the system. However, this not necessarily infinitely valid since the decrease of source impedance also has a negative effect on some of the coefficients (e.g., a 3 ). As validation of this observation, the maximum real part of the eigenvalues as function of the source capacitance is shown in Figure 8, in which it can be seen that a reduction of the source capacitance leads to increased stability, but only up until a certain capacitance. The parameters used for both numerical calculations, unless otherwise specified, are given in Table 1.  Moreover, since the effect of the resistance has been excluded from the mathematical derivations (for the sake of simplicity), another numerical calculation was performed to show the effect of resistance on the damping in the system. Figure 8 shows that some resistance is beneficial for the stability, while larger values deteriorate stability. The positive effect can be explained by the additional damping in the system while the negative effect on stability stems from the current limiting characteristic of the resistance.

Simulation Results
In this section, two simulations are performed for all three example systems to illustrate the results of the previous sections. Firstly, the configurations are simulated with a too high droop impedance to illustrate one of the unstable conditions of DC distribution systems. Secondly, the configurations are simulated with an acceptable droop impedance to illustrate the behavior of stable DC distribution systems and to compare it to the previous simulations.
For the simulations the example configurations, shown in Figures 5-7, are used. For each configuration, the source is situated at node N 1 , while loads are situated at N 2 , N 3 and N 4 . The system parameters, which are used for all simulations, are given in Table 2. In the table, R l , L l and C l are the resistance, inductance and capacitance of the distribution lines, respectively. Furthermore, C c is the (grid-side) output capacitance of the source and load converters. For the first set of simulations, the source reference voltage is initially set at 350 V with a droop impedance of 22 Ω. At t = 5 ms, the reference voltage of the source is changed to 375 V. Subsequently, at t = 10 ms a (constant power) load of 500 W is switched on at N 2 . Lastly, at t = 15 ms, a load of 1000 W is switched on at N 3 .
The results of these simulations are shown in Figure 9. It is seen that the system reacts stably to changes in load and voltage reference until the total load exceeds the stability conditions. Consequently, after t = 15 ms, the current rises exponentially over time and the voltage falls to zero. For the second simulations, the source reference voltage is initially set at 350 V with a droop impedance of 2 Ω. At t = 5 ms, the reference voltage of the source is changed to 375 V. Afterwards, at t = 10, 15, 20 ms loads of 500, 1000 and 750 W are switched on at N 2 , N 3 and N 4 , respectively.
The results of these simulations are shown in Figure 10. Firstly, it is seen that the voltage deviation from the reference voltage is lower compared to the previous simulations. Secondly, due to the lower droop impedance, the systems achieve their steady state faster than in the previous simulations (but this also causes some oscillations). Lastly, the systems are stable as the oscillations damp out and the voltages remain stably at their operating point.
In Figure 9, it was seen that instability occurs when the source impedance is too high. In this case, there does not exist a steady state voltage as a result of the interaction between the droop controlled source and the constant power load. Furthermore, the systems reach a state where they are no longer stable (and have eigenvalue(s) with a positive real part). Figure 10 shows that, since the droop controlled source has a low impedance, the system reaches a steady state faster. The second simulations also show oscillations, which are caused by the stronger reaction of the source to a deviation in voltage. These oscillations originate from the interaction between the capacitance of the nodes through the inductance of the distribution lines. If the source impedance and line resistance are too low, these oscillations might not be sufficiently damped out and a different type of instability can occur. In this case, the steady state voltage might exist, but oscillations around this steady state voltage increase or are not damped (enough). Conceptually, this can also be explained. If the impedance of the droop-controlled source is small, it will approximately behave as a constant voltage source. Consequently, the source will provide no damping in the system at all (see the first example in Section 3). In other words, for DC distribution systems to be called (sufficiently) stable, there should exist a stable steady state and there should be enough damping on the oscillations in the system.

Conclusions
This paper introduced a state-space model that represents any DC distribution system by its nodes, distribution lines, and source and load converters. The model allows for straightforward formulation of the state-space matrices, which relies on the validity of averaging, linearization and simplification of the electronic power converters.
The stability of DC distribution systems is derived from the eigenvalues of their state-space matrices. The stability of two simple examples was derived to illustrate the approach. Furthermore, a generalized method was presented to algebraically analyze the stability of any DC distribution system. To guarantee that the eigenvalues have only negative real parts, and thus guarantee stability, a method to derive the necessary and sufficient conditions was presented.
Subsequently, three examples were used to demonstrate how the presented method can lead to requirements for the stability of DC distribution systems. The three example systems showed strong congruency with respect to their sensitivity towards, for example, inductance and capacitance. Moreover, a sensitivity analysis was performed to verify the relationship between the inductance, droop source impedance, capacitance and resistance to the stability (and damping) of a system. The results showed that increasing inductance has a negative effect on stability, that increasing load capacitance has a positive effect, and that the droop impedance has a upper and lower limit (in agreement with previous research). However, it was also shown that both large and small values of the resistance and source capacitance have a negative effect on the stability of the system. This algebraically appends and verifies conclusions on parameter sensitivity of previous research for a large variety of DC distribution systems.
Additionally, various simulations were performed to show the dynamic behavior of a DC distribution system under stable and unstable conditions. These simulations showed two main modes of instability. Firstly, when the droop controlled sources' impedance is too large. In this case, the constant power loads' current increases as the voltage drops, and, as a result, the voltage of the entire system drops to 0 V. Secondly, if there is not enough damping from the source (and distribution lines' resistance). In the latter case, the insufficient damping results in increasing oscillations (which originate from the interactions between the different node capacitance).
In the future, the developed method could be used to evaluate and improve the stability of DC distribution systems. Furthermore, it could be used by the control system to, for example, assess if loads are allowed to turn on, if loads should switch off, or to change the control settings of sources.