Stability Analysis and Stability Enhancement Based on Virtual Harmonic Resistance for Meshed DC Distributed Power Systems with Constant Power Loads

Abstract: This paper addresses the stability issue of the meshed DC distributed power systems (DPS) with constant power loads (CPLs) and proposes a stability enhancement method based on virtual harmonic resistance. In previous researches, the network dynamics of the meshed DC DPS are often ignored, which affects the derivation of the equivalent system impendence. In addition, few of them have considered the meshed DC DPS including multiple sources with voltage-controlled converters and CPLs. To tackle the aforementioned challenge, this paper mainly makes the following efforts. The component connection method (CCM) is employed and expanded to derive the stability criterion of the meshed DC DPS with CPLs. This stability criterion can be simplified to relate only with the network node admittance matrix, the output impendences of the sources, and the input admittances of the CPLs. A virtual harmonic resistance through the second-order generalized integrator (SOGI) is added in the source with the voltage-controlled converter to lower the peak of the source output impendence, which can enhance the stability of the meshed DC DPS. The effectiveness of the proposed stability criterion and stability enhancement method are verified by nonlinear dynamic simulations.


Introduction
With the increasing penetration of renewable energy generation into the modern electric girds, meshed DC distributed power systems (DPS) have been widely adopted [1][2][3].The meshed DC DPS physically is composed of smaller power modules/subsystems.Usually, these smaller power modules/subsystems are designed only based on the stability requirement in its stand-alone operation, and thus it can operate well in the stand-alone application.However, the meshed DC DPS may become unstable due to the interaction among the modules/subsystems and the network [4].Furthermore, the negative resistance characteristic of the constant power load (CPL) is an important unstable factor.
Many stability methods in previous researches are given to analyze the stability of the single-source-single-load system or the parallel-source-parallel-load system [4][5][6].However, the network dynamics of the meshed DC DPS are often ignored, which affects the derivation of the equivalent system impendence.In addition, few of them have considered the meshed DC DPS including multiple sources with voltage-controlled converters and CPLs.A general approach for the stability analysis of the meshed DC DPS is to build its whole state-space model, and to identify the eigenvalues of the state matrix [7].However, this method requires the detailed models of loads and network dynamics, and the formulation of the system matrices may be very complex.To overcome this Energies 2017, 10, 69 2 of 15 problem, the component connection method (CCM) is introduced for the stability analysis [8].In the CCM, the dynamics of network and power modules/subsystems are modeled separately by a set of two vector-matrix equations, which is beneficial to reduce the computation burden of formulating the system transfer matrices [9].For the ac power-electronics-based power system, the CCM is employed to derivate the closed-loop response transfer matrix of the overall power which can predict the instability of the power system [10], but it does not separate the source impedance and load admittance with the network node admittance matrix.Thus, it is difficult to obtain the general stability criterion intuitionally.In this paper, the CCM is also employed and expanded to derive the stability criterion for the meshed DC DPS with the CPLs.The stability criterion can be simplified to relate only with the network node admittance matrix, the output impendence of the source, and the input admittance of the CPLs which can be separately obtained by theoretical calculation or impedance analyzer measurement.
Many stability enhancement methods have been given in previous researches.An adaptive active capacitor converter (AACC) is introduced to effectively stabilize the cascaded system [11], but this method needs additional hardware devices.Another approach is to add a virtual adaptive parallel impedance [12] or adaptive series-virtual-impedance [13] in the input of the load converter to improve the stability in the cascaded DC/DC converter system.However, in the meshed DC DPS, it is not very effective due to the existence of an electromagnetic interference (EMI) filter for the CPL with a switching-mode power converter.The input admittance of the CPLs is mainly determined by the EMI filter, and the virtual impedance is difficult to be added by the load converter control part.Another way is to adjust the source output impedance.Modifying the source output impedance to be zero [14] and adding virtual impedances in the source output filters [15] are both feasible strategies by the source converter control part, but both of them will influence the low and high frequency characteristics of the source output impedance.To overcome the problem, this paper proposes an effective and achievable method, which is to add a virtual harmonic resistance through the second-order generalized integrator (SOGI) in the source converter control part to lower the peak of the source output impendence.
The main contributions of this paper are follows.Firstly, the stability criterion for the meshed DC DPS is derived, which has been simplified to relate only with the network node admittance matrix, the output impendence of the source, and the input admittance of the CPLs.It is easy to implement, as the impedance and the admittance can both be separately obtained by the theoretical calculation or impedance analyzer measurement.Secondly, an effective stability enhancement method by adding a virtual harmonic resistance is proposed, with which only the small range middle frequency characteristic of the source output impedance is modified, while the low and high frequency characteristics of the source output impedance are still unmodified.The results of the nonlinear dynamic simulations verify the effectiveness of the proposed stability criterion and stability enhancement method.
The remainder of this paper is organized as follows.The stability criterion for the meshed DC DPS is given in Section 2. The modeling of the source with voltage-controlled converters and the CPLs are derived, and the proposed stability enhancement method and the stability analysis for the meshed DC DPS are given in Section 3.Then, the nonlinear dynamic simulations are conducted in Section 4. At last, conclusions are drawn in Section 5.

Stability Criterion for the Meshed DC DPS
The general form of DC DPS is given in [4], but it ignores the network dynamics.In [4], any converter in a DC DPS can be classified as either a bus voltage controlled converter (BVCC) or a bus current controlled converter (BCCC).In this study, a four node meshed DC DPS as Figure 1 is considered as an example, which is composed by two sources with voltage-controlled converters and two CPLs.The sources with voltage-controlled converters are controlled with constant voltage, hence they can be treated as BVCCs.The CPLs are controlled to absorb constant power from the meshed DC DPS, hence they can be treated as BCCCs.
interactions through the power cables in the meshed DC DPS.This consequently necessitates the use of CCM to analyze the stability of the meshed DC DPS.
The block diagram of the CCM applied for the studied meshed DC DPS is shown in Figure 1, where the CCM decomposes the overall DPS into four subsystems and the network.The two BVCCs are modeled by the Thevenin equivalent circuits [16], and the two BCCCs are modeled by the Norton equivalent circuits [17].In Figure 1 is the node voltage column vector, and [ ] is the node injection current column vector.[ ] is the cable current column vector.[ ] is the inductance column vector of cables, [ ] is the parasitic resistance column vector of cables, and [ ] is the capacitance column vector of cables.[ ] is the voltage reference-to-output transfer function column vector of the source and CPLs.
is the reference voltage column vector of the source and CPLs, respectively.Generally, the source is designed to be a voltage source that is stable when unloaded, and the CPL is designed to be stable when supplied by an ideal voltage source.That is, the unterminated behaviors of inverters [ ] are stable, and [ ] are stable, which means that there are no right-half plane poles.
To facilitate the formulation of the nodal admittance matrix, the Thevenin circuits of the BVCCs are converted to Norton circuits [10].Then, the studied meshed DC DPS is represented as Figure 2.Although all of the sources and CPLs are independently designed well, there will be complex interactions through the power cables in the meshed DC DPS.This consequently necessitates the use of CCM to analyze the stability of the meshed DC DPS.
The block diagram of the CCM applied for the studied meshed DC DPS is shown in Figure 1, where the CCM decomposes the overall DPS into four subsystems and the network.The two BVCCs are modeled by the Thevenin equivalent circuits [16], and the two BCCCs are modeled by the Norton equivalent circuits [17].
In Figure 1, Z v1 and Z v2 are the closed-loop output impedances of the sources, and Y c1 and Y c2 are the closed-loop input admittances of the CPLs.
is the node voltage column vector, and is the node injection current column vector.
is the cable current column vector.R L3 and R L4 are the resistive loads connected to node 3 and node 4, respectively.
is the parasitic resistance column vector of cables, and is the voltage reference-to-output transfer function column vector of the source and CPLs.v v1re f v v2re f v c1re f v c2re f T is the reference voltage column vector of the source and CPLs, respectively.
Generally, the source is designed to be a voltage source that is stable when unloaded, and the CPL is designed to be stable when supplied by an ideal voltage source.That is, the unterminated behaviors of inverters G v1 G v2 G c1 G c2 T are stable, and Z v1 Z v2 Y c1 Y c2 T are stable, which means that there are no right-half plane poles.
To facilitate the formulation of the nodal admittance matrix, the Thevenin circuits of the BVCCs are converted to Norton circuits [10].Then, the studied meshed DC DPS is represented as Figure 2. Applying the node voltage equation to the external dashed box part in Figure 2, it can be derived as where sys Y is the network nodal admittance matrix as =  .Then, the node voltage can be derived from (1) as follows:  Applying the node voltage equation to the external dashed box part in Figure 2, it can be derived as where Y sys is the network nodal admittance matrix as where Y L3 and Y L4 are the admittances for loads connected to node 3 and node 4, respectively.Y cable1 ~Ycable4 are the admittances for cables and Then, the node voltage can be derived from (1) as follows: where Z v is the set of the closed-loop output impedances of the sources, and Y c is the set of the closed-loop input admittances of the CPLs.v is the set of node voltages, and i is the set of the node injection currents.G v is the the set of voltage reference-to-output transfer functions of the source, G c is the the set of voltage reference to output transfer functions of the CPLs, v vre f is the set of voltage references of the source, and v cre f is the set of voltage references of the CPLs.I c is a unit diagonal matrix whose dimension is equal to Y c .I v is a unit diagonal matrix whose dimension is equal to Z v .And they satisfy the following relations: Applying the node voltage equation to the internal dashed box part in Figure 2, i can be derived as: In the above equation, Y sys is the network nodal admittance matrix for a real physical system, which means Y sys is stable and has no right-half plane poles.G v and G c are both stable as in the above analysis.Therefore, the stability criterion of the meshed DC DPS T m is If the stability criterion T m has right-half plane poles, the meshed DC DPS will be unstable.Through the analysis of the stability criterion T m , it can be found that the stability of the meshed DC DPS is only related to the network node admittance matrix Y sys , the output impendence of the sources Z v , and the input admittance of the CPLs Y c , which can be separately obtained by theoretical calculation or impedance analyzer measurement.

Modeling of the CPL
The modeling of the CPL connected to node 3 is similar to the CPL connected to node 4. Without loss of generality, taking the CPL connected to node 4 as an example.The CPL usually exploits a switching-mode power converter as the interface with the node, and the compliance with EMI standards usually requires the insertion of an EMI filter between the switching-mode power converter and the node [18][19][20].The close-loop circuit of the CPL with an EMI filter is shown as Figure 3.
Energies 2017, 10, 69 5 of 16 [ ] Applying the node voltage equation to the internal dashed box part in Figure 2, i can be derived as: In the above equation, sys Y is the network nodal admittance matrix for a real physical system, which means sys Y is stable and has no right-half plane poles.v G and c G are both stable as in the above analysis.Therefore, the stability criterion of the meshed DC DPS m T is If the stability criterion m T has right-half plane poles, the meshed DC DPS will be unstable.
Through the analysis of the stability criterion m T , it can be found that the stability of the meshed DC DPS is only related to the network node admittance matrix sys Y , the output impendence of the sources v Z , and the input admittance of the CPLs c Y , which can be separately obtained by theoretical calculation or impedance analyzer measurement.

Modeling of the CPL
The modeling of the CPL connected to node 3 is similar to the CPL connected to node 4. Without loss of generality, taking the CPL connected to node 4 as an example.The CPL usually exploits a switching-mode power converter as the interface with the node, and the compliance with EMI standards usually requires the insertion of an EMI filter between the switching-mode power converter and the node [18][19][20].The close-loop circuit of the CPL with an EMI filter is shown as Figure 3. Applying switch period average, the small-signal model of the CPL with EMI filter can be achieved as Figure 4. Applying switch period average, the small-signal model of the CPL with EMI filter can be achieved as Figure 4.
Energies 2017, 10, 69 6 of 16 It can be derived from Figure 4 that the small-signal disturbance of the output current 4 î can be expressed as: (10)

Modeling of the Source with Voltage-Controlled Converters
The modeling of the source with voltage-controlled converters connected to node 1 is similar to the source connected to node 2. Without loss of generality, taking the source connected to node 1 as an example, the close-loop circuit of the source with voltage-controlled converters is shown as Figure 5.It can be derived from Figure 4 that the small-signal disturbance of the output current î4 can be expressed as: where G c2v is the voltage controller, Y c2 is the closed-loop input admittance of the CPL, and G c2 is the voltage reference to output transfer function of the CPL. a 0 ∼ a 4 , b 0 ∼b 5 , and c 0 ∼c 3 are given in Appendix A.

Modeling of the Source with Voltage-Controlled Converters
The modeling of the source with voltage-controlled converters connected to node 1 is similar to the source connected to node 2. Without loss of generality, taking the source connected to node 1 as an example, the close-loop circuit of the source with voltage-controlled converters is shown as Figure 5. Figure 6 shows the small-signal block diagram of the multi-loop control system of the source with voltage-controlled converters.
The small-signal block diagram of the multi-loop control system of the source.
Thus, the dynamic behavior of the close-loop system can be given by where where is the closed-loop output impedance of the source with voltage-controlled converters, and is the voltage reference-to-output transfer function of the source with voltage-controlled converters.

Stability Enhancement Method
The source output impedance has a peak at the resonant frequency.Based on the Nyquist criterion, instability between the source and the network may occur if the network input impedance becomes lower than the source output impedance at this frequency.Therefore, to ensure the system stability and to minimize the potential for inadvertent interactions with the source, it is important to lower the peak of the source output impedance.For this purpose, this paper proposes a method by adding a virtual harmonic resistance through the SOGI to be in parallel with the capacitor.
The transfer function f G s of SOGI can be expressed as Figure 6 shows the small-signal block diagram of the multi-loop control system of the source with voltage-controlled converters.
Energies 2017, 10, 69 7 of 16 Figure 6 shows the small-signal block diagram of the multi-loop control system of the source with voltage-controlled converters.
The small-signal block diagram of the multi-loop control system of the source.
Thus, the dynamic behavior of the close-loop system can be given by where where is the closed-loop output impedance of the source with voltage-controlled converters, and is the voltage reference-to-output transfer function of the source with voltage-controlled converters.

Stability Enhancement Method
The source output impedance has a peak at the resonant frequency.Based on the Nyquist criterion, instability between the source and the network may occur if the network input impedance becomes lower than the source output impedance at this frequency.Therefore, to ensure the system stability and to minimize the potential for inadvertent interactions with the source, it is important to lower the peak of the source output impedance.For this purpose, this paper proposes a method by adding a virtual harmonic resistance through the SOGI to be in parallel with the capacitor.
The transfer function

( )
f G s of SOGI can be expressed as Thus, the dynamic behavior of the close-loop system can be given by where where Z v1 is the closed-loop output impedance of the source with voltage-controlled converters, and G v1 is the voltage reference-to-output transfer function of the source with voltage-controlled converters.

Stability Enhancement Method
The source output impedance has a peak at the resonant frequency.Based on the Nyquist criterion, instability between the source and the network may occur if the network input impedance becomes lower than the source output impedance at this frequency.Therefore, to ensure the system stability and to minimize the potential for inadvertent interactions with the source, it is important to lower the peak of the source output impedance.For this purpose, this paper proposes a method by adding a virtual harmonic resistance through the SOGI to be in parallel with the capacitor.The transfer function G f (s) of SOGI can be expressed as where ω is the resonant frequency, and k is the frequency coefficient.The frequency-domain characteristic of G f (s) is shown as Figure 7. Around the resonant frequency, the amplitude of G f (s) is very small, and the amplitude of G f (s) in other frequency ranges is 0 dB.
Energies 2017, 10, 69 8 of 16 where ω is the resonant frequency, and k is the frequency coefficient.The frequency-domain characteristic of ( ) f G s is shown as Figure 7. Around the resonant frequency, the amplitude of is very small, and the amplitude of G s in other frequency ranges is 0 dB. Figure 8 shows the small-signal control block of the source 1.If a virtual harmonic resistance vh R is required to be added in parallel with the capacitor, it is intuitively to introduce vh R to the block as a dot-dashed part in Figure 8, which is the basic idea of this control method.However, this method cannot be realized by control directly.Therefore, vh R is moved to the output of 1 ( ) G s as the dashed part in Figure 8.

Figure 8. Adding a virtual harmonic resistance through the second-order generalized integrator (SOGI).
( ) R G s is expressed as With the parallel virtual harmonic resistance , the frequency-domain characteristic of the source output impedance is shown as Figure 9. Figure 8 shows the small-signal control block of the source 1.If a virtual harmonic resistance R vh is required to be added in parallel with the capacitor, it is intuitively to introduce R vh to the block as a dot-dashed part in Figure 8, which is the basic idea of this control method.However, this method cannot be realized by control directly.Therefore, R vh is moved to the output of G v1v (s) as the dashed part in Figure 8.
Energies 2017, 10, 69 8 of 16 where ω is the resonant frequency, and k is the frequency coefficient.The frequency-domain characteristic of ( ) f G s is shown as Figure 7. Around the resonant frequency, the amplitude of is very small, and the amplitude of G s in other frequency ranges is 0 dB. Figure 8 shows the small-signal control block of the source 1.If a virtual harmonic resistance vh R is required to be added in parallel with the capacitor, it is intuitively to introduce vh R to the block as a dot-dashed part in Figure 8, which is the basic idea of this control method.However, this method cannot be realized by control directly.Therefore, vh R is moved to the output of 1 ( ) G s as the dashed part in Figure 8.

Figure 8. Adding a virtual harmonic resistance through the second-order generalized integrator (SOGI).
( ) R G s is expressed as With the parallel virtual harmonic resistance , the frequency-domain characteristic of the source output impedance is shown as Figure 9. G R (s) is expressed as With the parallel virtual harmonic resistance R vh = 1 Ω, the frequency-domain characteristic of the source output impedance is shown as Figure 9.As can be seen, the peak of the source output impedance has been lower, and only the small range middle frequency characteristic of the source output impedance is modified, while the low and high frequency characteristics of the source output impedance are still unmodified.

Stability Analysis
All the main electrical parameters of the studied DC meshed DPS as shown in Figure 1 are listed in Appendix B.
The closed-loop output impedances of the sources v Z can be represented by using impedance calculated in Equations ( 5) and ( 14), the closed-loop input admittances of the CPLs c Y can be represented by using admittance calculated in Equations ( 5) and (10), and the network nodal admittance matrix sys Y can be represented by using admittance calculated in Equation (2).Then, the the stability of the meshed DC DPS can be analyzed with the stability criterion m T expressed in Equation (8).If the stability criterion m T has right-half plane poles, the meshed DC DPS will be unstable.
At first, it is necessry to compare the proposed stability criterion m T in this paper and the previous stability criterion 2 m T given in [4] as follows: The dominant poles of m T ans 2 As can be seen, the peak of the source output impedance has been lower, and only the small range middle frequency characteristic of the source output impedance is modified, while the low and high frequency of the source output impedance are still unmodified.

Stability Analysis
All the main electrical parameters of the studied DC meshed DPS as shown in Figure 1 are listed in Appendix B.
The closed-loop output impedances of the sources Z v can be represented by using impedance calculated in Equations ( 5) and ( 14), the closed-loop input admittances of the CPLs Y c can be represented by using admittance calculated in Equations ( 5) and (10), and the network nodal admittance matrix Y sys can be represented by using admittance calculated in Equation (2).Then, the the stability of the meshed DC DPS can be analyzed with the stability criterion T m expressed in Equation (8).If the stability criterion T m has right-half plane poles, the meshed DC DPS will be unstable.
At first, it is necessry to compare the proposed stability criterion T m in this paper and the previous stability criterion T m2 given in [4] as follows: The dominant poles of T m ans T m2 are displayed as Figure 10, in which the cables are long R cablei = 0.05 Ω , L cablei = 0.5 mH , C cablei = 100 µF, and i = 1 ∼4.As shown in Figure 10, all the dominant poles of T m2 are both in the left-half plane, which indicates that the meshed DC DPS is stable.However, there are two right-half plane poles 0.15 ± 94.5i in the dominant poles of T m , which indicates the meshed DC DPS is unstable and have an oscillation with approximately period 0.066 s.The results of the proposed stability criterion T m and the previous stability criterion T m2 are not consistent, and the results of proposed stability criterion T m will be confirmed to be right with the nonlinear dynamic simulations in Section 4.Moreover, many poles related to the network dynamics are also lost in T m2 compared to T m .and the dominant poles of m T are displayed as Figure 13.
As can be seen from Figure 11, all of the poles of m T are in the left half plane, which indicates that the meshed DC DPS is stable.But there are two right half plane poles 0.48 97i ± in the Figure 12, which indicates the meshed DC DPS is unstable and have an oscillation with approximately period 0.0647 s.Furthermore, with the power reducation of CPL2, the meshed DC DPS becomes stable as shown in Figure 13.With the electrical parameters listed in Table B1, when the power of resistive load 4 P L4 is 50 kW, which means R L4 = 5 Ω, the dominant poles of T m are displayed as Figure 11.When the power of resistive load 4 P L4 reduces from 50 kW to 10 kW, which means R L4 = 25 Ω, the dominant poles of T m are displayed as Figure 12.Furthermore, when the power of resistive load 4 P L4 is 10 kW and the absorbed power of CPL2 P c2 reduces from 100 kW to 50 kW, which means R L4 = 25 Ω and R c2 = 1.25 Ω, and the dominant poles of T m are displayed as Figure 13.and the dominant poles of m T are displayed as Figure 13.
As can be seen from Figure 11, all of the poles of m T are in the left half plane, which indicates that the meshed DC DPS is stable.But there are two right half plane poles 0.48 97i ± in the Figure 12, which indicates the meshed DC DPS is unstable and have an oscillation with approximately period 0.0647 s.Furthermore, with the power reducation of CPL2, the meshed DC DPS becomes stable as shown in Figure 13.As can be seen from Figure 11, all of the poles of T m are in the left half plane, which indicates that the meshed DC DPS is stable.But there are two right half plane poles 0.48 ± 97i in the Figure 12, which indicates the meshed DC DPS is unstable and have an oscillation with approximately period 0.0647 s.Furthermore, with the power reducation of CPL2, the meshed DC DPS becomes stable as shown in Figure 13.
Applying the proposed stability enhancement method in the sources connected with node 1 and node 2 of the meshed DC DPS when P L4 is 10 kW and P c2 is 100 kW, the dominant poles of T m are displayed as Figure 14.There are no right half plane poles, which means that the meshed DC DPS is stable.Compared to the poles in Figure 12, the stability of the meshed DC DPS is enhanced by the proposed stability enhancement method.

Nonlinear Dynamic Simulations
To validate the effectiveness of the proposed stability criterion and stability enhancement method, the meshed DC DPS in Figure 1 is built in the nonlinear time-domain simulation by using Matlab/Simulink.The electrical and controller parameters are all given in Appendix B.
To compare the proposed stability criterion T m and the previous stability criterion T m2 , the cables are chosen to be long cables with parameters R cablei = 0.05 Ω , L cablei = 0.5 mH , C cablei = 100 µF, i = 1∼4, the nonlinear dynamic simulation results of the node 1 voltage v 1 (t) and the cable 1 current i cable1 (t) are shown in Figure 15.As can be seen in Figure 15, before 2.5 s, the meshed DC DPS gradually becomes to a stable state with P L4 = 50 kW, and once the power of resistive load 4 P L4 reduces to 15 kW at 2.5 s, the meshed DC DPS becomes to an oscillation state.The oscillation amplitude gradually grows, and the oscillation period is about 0.066 s, which satisfy the theoretical analysis in Figure 10.The results verify the effectiveness of the proposed stability criterion T m compared to the previous stability criterion T m2 which does not take the network dynamics into account., which confirms the theoretical analysis in Figure 11.Then, the power of resistive load 4 4 L P reduces to 10 kW at 1.5 s, the meshed DC DPS becomes to the oscillation state with about period 0.065 s, which confirms the theoretical analysis in Figure 12.The power of CPL2 reduces from 100 kW to 50 kW at 2 s, and the meshed DC DPS becomes to a stable state again, which confirms the theoretical analysis in Figure 13.To verify the effectiveness of the proposed stability criterion T m with power changes of the resistive load and CPL, and the proposed stability enhancement method, the nonlinear dynamic simulation results are shown as the following Figures 16-18.
As can be seen in Figure 16, before 1.5 s, the meshed DC DPS gradually becomes to stable state with P L4 = 50 kW and P c2 = 100 kW, which confirms the theoretical analysis in Figure 11.Then, the power of resistive load 4 P L4 reduces to 10 kW at 1.5 s, the meshed DC DPS becomes to the oscillation state with about period 0.065 s, which confirms the theoretical analysis in Figure 12.The power of CPL2 reduces from 100 kW to 50 kW at 2 s, and the meshed DC DPS becomes to a stable state again, which confirms the theoretical analysis in Figure 13.As can be seen in Figure 17, when the power of CPL2 increases from 50 kW to 100 kW at 3.5 s, the meshed DC DPS becomes to the oscillation state again between 3.5-4.5 s.Then, the vitural harmonic resistances are added in the controllers of source 1 and source 2 at 4.5 s.With the vitural harmonic resistance, the meshed DC DPS becomes to the stable state again, and the node 1 voltage v 1 converges again to the steady value 500 V.The simulation results between 4.5-5.5 s confirm the theoretical analysis results in Figure 14.The power of resistive load 4 P L4 increases to 50 kW at 5.5 s, and the meshed DC DPS is maintained at a steady state.
The currents of cable 1 i cable1 and cable 4 i cable4 are displayed as Figure 18.The steady state and unstable state are similar to the node 1 voltage v 1 as shown in Figures 16 and 17.It is worth noting that the steady state current values in 0-1.5 s and 5.5-6 s are the same, which means that with or without the proposed stability enhancement method, the current distributions in the meshed DC DPS are not affected, which verifies that the low frequency characteristic of the source output impedance has not been modified.

Conclusions
This paper discusses the modeling and stability analysis of the meshed DC DPS including multiple sources with voltage-controlled converters and CPLs, and proposes a stability enhancement method.The stability criterion is found to relate only with the network node admittance matrix, the output impendence of the source, and the input admittance of the CPLs, and the three parts can be obtained separately.Virtual harmonic resistance through the second-order generalized integrator (SOGI) for the source proposed in this paper is a decoupling design which is only based on the characteristic of the source output impedance.Nonlinear dynamic time-domain simulation results verify the effectiveness of the proposed stability criterion and stability enhancement method.Furthermore, the current control sources or loads can be treated as BCCCs, and the input admittances of current control sources or loads are related to the stability analysis.Multiple sources with droop control can be treated as BVCCs, and the output impendence of the droop control sources are related to the stability analysis.So the proposed stability criterion can also deal with the current control sources or loads and multiple sources with droop control.Related analyses are in the future researches.

Figure 1 .
Figure 1.The block diagram of the component connection method (CCM) applied for the studied meshed DC distributed power systems (DPS).
closed-loop output impedances of the sources, and 1 c Y and 2 c Y are the closed-loop input admittances of the CPLs.[ ]

3 L R and 4 LR
are the resistive loads connected to node 3 and node 4, respectively.

Figure 1 .
Figure 1.The block diagram of the component connection method (CCM) applied for the studied meshed DC distributed power systems (DPS).

Figure 2 .
Figure 2. The equivalent Norton circuit of the studied meshed DC DPS.
) where v Z is the set of the closed-loop output impedances of the sources, and c Y is the set of the closed-loop input admittances of the CPLs.v is the set of node voltages, and i is the set of the node injection currents.v G is the the set of voltage reference-to-output transfer functions of the source, c G is the the set of voltage reference to output transfer functions of the CPLs, vref v is the set of voltage references of the source, and cref v is the set of voltage references of the CPLs.c I is a unit diagonal matrix whose dimension is equal to c Y .v I is a unit diagonal matrix whose dimension is equal to v Z .And they satisfy the following relations:

Figure 2 .
Figure 2. The equivalent Norton circuit of the studied meshed DC DPS.

4 iFigure 3 .
Figure 3.The close-loop circuit of the constant power load (CPL) with an electromagnetic interference (EMI) filter.

Figure 3 .
Figure 3.The close-loop circuit of the constant power load (CPL) with an electromagnetic interference (EMI) filter.

Figure 4 .
Figure 4.The small-signal model of the CPL with an EMI filter.

2 cY 5 ~
is the closed-loop input admittance of the CPL, and 2 c G is the voltage reference to output transfer function of the CPL.0 4 0

Figure 4 .
Figure 4.The small-signal model of the CPL with an EMI filter.

Figure 5 .
Figure 5.The close-loop circuit of the source with voltage-controlled converters (where 1 v v G and

Figure 5 .
Figure 5.The close-loop circuit of the source with voltage-controlled converters (where G v1v and G v1i are the voltage and current controllers).

Figure 5 .
Figure 5.The close-loop circuit of the source with voltage-controlled converters (where 1 v v G and

Figure 6 .
Figure 6.The small-signal block diagram of the multi-loop control system of the source.

Figure 7 .
Figure 7.The frequency-domain characteristic of

Figure 7 .
Figure 7.The frequency-domain characteristic of G f (s).

Figure 7 .
Figure 7.The frequency-domain characteristic of

Figure 8 .
Figure 8. Adding a virtual harmonic resistance through the second-order generalized integrator (SOGI).

Figure 9 .
Figure 9.The source output impedance with/without parallel virtual harmonic resistance vh R .

.Figure 9 .
Figure 9.The source output impedance with/without parallel virtual harmonic resistance R vh .

Figure 10 . 4 LP 4 LPLR = Ω and 2
Figure 10.The dominant poles of m T ans 2 m T when the cables are long.

Figure 11 .Figure 10 .
Figure 11.The dominant poles of m T when 4

Figure 10 . 4 L P is 50 kW, which means 4 5 LR 4 LP 4 LP
Figure 10.The dominant poles of m T ans 2 m T when the cables are long.

Figure 14 .
Figure 14.The dominant poles of m T with the proposed stability enhancement method.

Figure 12 .Figure 13 .
Figure 12.The dominant poles of m T when 4

Figure 14 .
Figure 14.The dominant poles of m T with the proposed stability enhancement method.Figure 14.The dominant poles of T m with the proposed stability enhancement method.

Figure 14 .
Figure 14.The dominant poles of m T with the proposed stability enhancement method.Figure 14.The dominant poles of T m with the proposed stability enhancement method.

15 . 4 LP 2 mT
As can be seen in Figure15, before 2.5 s, the meshed DC DPS gradually becomes to a stable state with 4 50 kW L P = , and once the power of resistive load 4 reduces to 15 kW at 2.5 s, the meshed DC DPS becomes to an oscillation state.The oscillation amplitude gradually grows, and the oscillation period is about 0.066 s, which satisfy the theoretical analysis in Figure10.The results verify the effectiveness of the proposed stability criterion m T compared to the previous stability criterion which does not take the network dynamics into account.

Figure 15 .
Figure 15.Comparison between m T and 2 m T .

Figure 15 .
Figure 15.Comparison between T m and T m2 .

Figure 18 .
Figure 18.The currents of cable 1 i cable1 and cable 4 i cable4 .