A Procedure to Determine the Droop Constants of Voltage Controllers Coping with Multiple DG Interactions in Active Distribution Systems

: In modern distribution systems, the presence of an increasing amount of Distributed Generation (DG) systems causes over-/under-voltage problems, due to the reverse power ﬂows. To face these problems, the voltage-reactive power droop controllers of DG systems are commonly used for their simplicity and are required by international standards. On the other hand, the interaction among voltage droop controllers of different DG systems may introduce instability. The paper presents an effective procedure to determine the droop constants of voltage-reactive power controllers for multiple DG systems. Firstly, a multi-input multi-output model of the distribution system is introduced. Then, using the concept of the interaction measure under decentralized control, a simple constraint is added to the single-input single-output design of each droop controller. Such a constraint guarantees stability with respect to the interaction among the voltage droop controllers of all the DG systems. Eventually, the proposed procedure is applied to an LV test system with 24 nodes and six photovoltaic systems; the results of numerical simulations are presented, giving evidence of the effectiveness of the proposed procedure in various operating conditions of the distribution system.


Introduction
Distributed Generation (DG) systems are widely spreading in electric distribution networks, and this phenomenon highly impacts the grid management and control. In fact, distribution systems were traditionally designed to distribute electric energy from a supplying substation, downstream along the feeders, to the loads. The presence of DG systems causes the inversion of the power flows and, consequently, a radical change of the voltage profiles along the feeders. On the other hand, DG systems are also capable of absorbing/injecting reactive power from/into the distribution grid, and such an ability can compensate the negative effects of active power injections [1][2][3]. DG reactive power is the first resource to be used for voltage support and overvoltage containment. In fact, as an alternative, active power curtailment [4][5][6] can be used, which has a significant impact on the voltage profile due to the high R/X ratio of distribution lines [7]. However, it has negative impacts on the DG economic revenues and on the environment, reducing the exploitation of renewable energy sources. Eventually, both active and reactive power controls could contribute to system voltage regulation (In grid-connected operation, small-sized DG systems do not directly partake in the frequency regulation, which is assured by the bulk transmission system. This is not true in islanded operation of microgrids, which is out of the scope of this paper.) if DG were equipped with energy storage systems [8], but this would significantly increase the DG installation costs.
In this frame, international standards are being revised to ask DG systems for adequate reactive power capabilities so as to support the Distribution System Operator (DSO) in voltage regulation [9,10]. To this aim, the DG systems must be equipped with a reactive power control loop, whose set-point is determined according to the chosen operation mode. Three basic types of operation mode are possible, namely: (1) Active power-reactive power mode (2) Constant reactive power mode (3) Voltage-reactive power mode In the first DG operation mode (active power-reactive power mode), the reactive power set-point Q ci of the ith DG system is determined as a function of the active power injected by the DG. Typical functions are chosen so that the reactive power absorption counteracts the effects of active power injection, especially overvoltages; also, constant power factor functions are used. However, this operation mode does not directly involve the DG in voltage regulation, and no general rule for the choice of the function can be established that guarantees an improvement of network voltage profiles in all operating conditions.
In the second DG operation mode (constant reactive power mode), the reactive power set-point Q ci is constant. This is the case of Q ci = 0, forcing the unitary power factor, which is typically adopted by DSOs in the case of the absence of any kind of real-time voltage control system (as is often the case in LV distribution networks). However, the fit-and-forget strategy is not a sustainable option in view of a large penetration of DG [11]. Indeed, the constant reactive power mode is also adopted in the case of real-time voltage control systems with a centralized architecture. In fact, the international standards impose that the set-point Q ci can be received as the input signal sent by the DSO. In such a centralized control scheme, the DSO monitors and estimates the system operating state and consequently evaluates and sends the reactive power set-points to the DG systems. The main benefit of centralized voltage control is that the voltage profiles along the feeders of the distribution network are optimized. The main drawback is related to the complexity of the measurement and communication infrastructure that is needed in the centralized control scheme: such a complexity involves large investment costs and reduces the reliability of the voltage control system [11].
To overcome the drawbacks of centralized voltage control schemes, decentralized approaches are adopted, which require limited or no data exchange [12]. Decentralized schemes assume that DG reactive power control is operated in the third mode (voltage-reactive power mode). In this mode, the set-point Q ci is determined locally by the DG system. In practice, the voltage amplitude V m i is measured at the node m i of the distribution network at which the ith DG system is connected, namely the Point of Common Coupling (PCC). Using the voltage measure, the DG voltage controller evaluates the reactive power set-point Q ci so as to partake in the real-time voltage control of the distribution system. Various structures of the DG voltage controllers have been proposed and designed. The simplest one is a proportional controller, or voltage droop controller, which implements a linear function Q ci = k ci V m i [13,14]. The voltage droop controller is the most commonly used for its simplicity, and its implementation on DSO's request is imposed by standards [9,10]. The voltage droop control can be enriched by nonlinear functions implementing dead-zones, and if the DG system is composed of multiple generators, as in a wind farm, the droop of each generator can be chosen so as to dispatch the reactive power among the generators. In this way, the different operating conditions of the generators can be accounted for [15], while assuring that the overall response of the DG system at the PCC still results in following the voltage droop characteristics.
Unfortunately, due to the distribution network, all the voltage amplitudes at the PCCs are affected by the reactive power absorbed/injected by all the DG systems. Then, the issue to be tackled in a decentralized voltage control scheme is related to the dynamic interactions among a large number of DG voltage control systems, which may generate unwanted oscillatory (hunting) behavior, causing instability [13,[16][17][18]. This problem is even more complex to solve in the absence of communication among DG systems [19]. Various solutions have been proposed to cope with mutual interaction by adopting more complex structures for the DG voltage controllers, as well as design methodologies based on Multi-Input Multi-Output (MIMO) models, so as to account for the dynamic interactions among DG voltage control systems [20][21][22][23]. However, such solutions lose the simplicity of the voltage droop controller, which is already provided in commercial DG systems according to the standards.
Motivated by the above considerations, this paper proposes a simple procedure to determine the constants of the voltage droop controllers. According to a decentralized approach, each controller uses only the local measure of the voltage amplitude at the PCC and does not require any real-time data exchange with the other DG systems and the DSO. Each droop constant can be designed using a Single-Input Single-Output (SISO) technique in order to stabilize the single loop. However, each loop interacts with the other ones, and consequently, when the all the droop controllers are simultaneously connected to the distribution system, instability phenomena may appear. To solve this problem, this paper, using the concept of the interaction measure under decentralized control presented in [24], proposes a simple constraint that each droop constant must meet at the design step.
The paper is structured as follows. After introducing the MIMO model of the system in Section 2, the concept of the interaction measure is used in Section 3 to obtain a simple procedure for the choice of the droop constants of the voltage-reactive power controllers of the DG systems. The application of the proposed procedure to a case study is presented in Section 4 together with the results of the numerical simulations.

System Modeling
Let a distribution system be considered as shown in Figure 1. The main system is a higher voltage system with imposed voltage frequency and amplitude (slack bus). Through the substation and its transformer, the main system supplies the distribution grid, which is generally operated in a radial configuration. It is assumed that DG systems are connected at various nodes of the distribution grid. Let m i , for i = 1, . . . , , represent the index of the network node to which the ith DG is connected.  In this section, a MIMO model of the distribution system in Figure 1 is introduced, which is suitable for voltage control design. The model is composed of two parts: the DG system and the distribution network (including the substation and the slack bus). Concerning the DG system, it is usually interfaced with the grid by a power converter. In the most typical configuration, the power converter is a Voltage Source Inverter (VSI) equipped with the Phase-Locked-Loop (PLL) and current control loops. The set-points to the inverter are the current components along the d − q axis, which are usually determined by the active and reactive power controllers, respectively. Eventually, the active power controller receives its set-point from the active power optimization functions (such as the maximum power point tracker in photovoltaic systems).
On the contrary, the set-point of the reactive power controller can be generated according to various schemes classified in Section 1.
The inverter presents very fast dynamics in terms of current responses, which can be neglected because the related pulse frequencies are far beyond the bandwidth of the closed-loop control of active/reactive powers [23,25]. Moreover, the inverter allows an independent control of the two current components, i d i and i q i , and, consequently, of active and reactive powers [26]. On the basis of these premises, the DG model can take into account only the closed-loop response of the reactive power control, which is given by [27]: where Q DG i is the reactive power injected into the distribution grid by the ith DG system, Q ci is the set-point that is determined by the voltage droop controller, and G i (s) is the transfer function of the reactive power control loop, being G i (0) = 1. Concerning the distribution system, its model must represent the effect of the reactive powers injected by DG systems (Q DG i ) on the voltage amplitudes of the nodes to which the DG systems are connected (V m i ), see Figure 1. To this aim, the steady-state operation of the distribution network and of the substation transformer, which is connected to the slack bus, must be modeled. The classical load-flow equations can be used, in particular the DistFlow equations in the case of radial topology. However, the result is a large set of nonlinear equations that are difficult to handle in view of the control design. Several linear approximations of the DistFlow equations have been proposed in the literature. In the following, reference is made to [28], which provided the following closed-form relation: with: (2) linearly relates the variations of the squared amplitudes of the nodal voltages at the PCCs of the DG systems with the reactive power injections. The voltage variations are expressed with respect to the squared voltage values without reactive power injections; Γ is a non-singular matrix of sensitivity coefficients belonging to IR × , whose elements γ ij depend on the operating conditions of the distribution system. The linear expression of the voltage amplitudes can easily be derived from (2): with: where the voltage amplitudes V 0 are evaluated assuming the distribution system operating without the reactive power injections, that is Q DG = 0; the elements of the sensitivity matrix Θ are evaluated as: In practice, the matrix Γ is evaluated by closed-form expressions according to [28], for the given operating conditions of the distribution system, and in turn, matrix Θ is obtained and its determination indicated by Θ 0 . When the operating conditions of the distribution system change, the elements of the matrix Θ are subject to variations with respect to the corresponding values assumed in Θ 0 . Typical changes of the operating conditions of the distribution system that affect the elements in Θ are: • variations of the loads; • changes in the power injections by the DG; • variations of the voltage amplitude imposed by the slack bus, due to variations of the operating conditions of the main system.
Defining: (1) for the case of DG systems is rewritten in compact form as: Substituting (4) in (3) yields the overall system model:

Droop Constants' Design
It is well known that the SISO design of each decentralized voltage controller focuses on a single control loop, but neglecting the interaction with the other voltage control loops, it may lead to instability. To avoid this drawback, a MIMO approach should be adopted to take into account the interaction among the control loops. However, in the case of a distribution system, the number of droop controllers is very large, and a MIMO approach may lead to a more complex design. For this reason, it would be convenient to account in the SISO design of each droop constant for an additional constraint that guarantees stability when all the controllers are connected to the distribution system. To this aim, let the decentralized control structure shown in Figure 2 be considered, in which V des is the vector of reference signals and K c = diag{k c1 k c2 . . . k c } is the diagonal matrix of the voltage droop controller's gains. Without loss of generality, in the remainder, it will assumed V des,i = 1.
According to the SISO approach, each gain k ci is designed using only the model P ii (s) = θ ii G i (s) of the ith DG, and each single loop is stable. Since the off diagonal elements P ij (s) = θ ij G i (s) of the matrix plant P(s) have been neglected, it is important to derive a condition that has to be satisfied by each droop constant in order to preserve overall control system stability.
The concept of the interaction measure for MIMO systems proposed in [24] is a suitable approach to meet the target.
whereΘ = diag{θ 11 θ 22 . . . θ }. Model (5) shows how the interaction is modeled as an additive uncertainty on the model plant; see Figure 3. The ith droop constant k ci is designed with a SISO technique that guarantees stability for the single closed-loop with transfer function: Consequently, the diagonal matrix W(s) is also stable, being: The interaction measure expresses the constraint imposed on the values of the controller's gains k ci in order to guarantee that the full closed-loop system is stable with transfer function: In the interaction analysis, the key role is played by matrix E, defined as: E can be viewed as a relative error with respect to the diagonal plantP(s), that is the error induced when the full plant P(s) is approximated by the diagonal plantP(s) employed in the SISO design stage.
It is important to note that E is a constant matrix depending only on matrices (Θ −Θ) andΘ. Theorem 2 in paper [24] states that the closed-loop system W(s) is stable if: where ρ denotes the spectral radius. Based on the generalized Nyquist criterion, Condition (7) is satisfied if the characteristic loci of the eigenvalues of EW( ω) do not encircle the critical point (−1, 0). Now, since: withσ the largest singular value, it results that (7) is satisfied if: where ||W( ω)|| ∞ max ωσ (W( ω)).
Usually, each transfer functionW i ( ω) presents a frequency response that is monotonically decreasing; then, ||W i ( ω)|| ∞ = |W i ( 0)|. Consequently, Condition (8) becomes: with: Substituting (10) in (9) and recalling that G i ( 0) = 1, it has: Inequality (11) then represents the constraint imposed by the interaction measure on k ci to guarantee closed-loop stability with full plant P(s). It is worth pointing out that (11) is a sufficient condition that defines a conservative bound on k ci , whereas Condition (7) does not give any information about it.
The application of the proposed criterion requires the evaluation of the upper bound in the condition (11). To this aim, several different operating conditions of the distribution system have to be considered. In each operating condition, matrix Θ in (3) is firstly calculated using the linear method in [28], and then, matrix E in (6) and its infinity norm are evaluated. At the end of this procedure, a set of matrices Θ and of corresponding values of ||E|| ∞ are obtained. Eventually, in this set, the element that corresponds to the smallest value of the right-hand side term in (11) is chosen and indicated as Θ , whereas ||E|| ∞ is the infinity norm of the corresponding E. Each droop constant must comply with the following constraint derived from (11): with i = 1, . . . . Fulfilling the additional constraint (12) in the SISO design guarantees that the chosen k ci does not cause instability due to the interaction with the other − 1 loops.
Concerning the possibility that one or more DG is switched off, it is important to notice that the proposed design still guarantees closed-loop stability. In fact, as a general rule, the loss of a DG system causes a reduction of the dimension of matrix E and a consequent decrease of its infinity norm. In Figure 4, an example of the variation of ||E|| ∞ with respect to the variation of the number of DG systems is reported. It is apparent that the relationship is approximately a direct proportionality. Then, the choice of k ci according to (12) guarantees that the stability condition (11) stands also in the case of the loss of a DG system.
Conversely, the case of the installation of a new DG represents a rare and known in advance event, which is planned by the DSO. In this circumstance, it is necessary to update the values of the droop constants k ci . In fact, the dimension increases by unity, and the DSO must recalculate a new matrix Θ and the corresponding value of ||E|| ∞ in (12).
As a final consideration, it is worth noting that an increase of the number of DG systems yields an increase of ||E|| ∞ (see Figure 4) and, consequently, a reduction of the values of the droop constants k ci in (12).

Case Study
The 24-node MV/LV distribution system shown in Figure 5 was considered. A MV/LV substation with a 20/0.4 kV-0.25 MVA transformer supplied two LV main feeders with subfeeders. The electric parameters of the lines and the rated loads are reported in Table 1 on a 25 kVA power basis. Concerning the DG, six 20 kW PV systems were connected at different nodes along the feeders and subfeeders.
The distribution system was simulated by PSCAD/EMTDC [29]. Each PV system was simulated in detail as in [26], including the PV generator non-linear equations, the DC bus, the VSI, and the connecting transformer, both with rated power equal to 25 kVA. Each VSI was equipped with current control loops, which imposed the current components i d and i q in a frame of dq coordinates. This frame was synchronized with the AC system voltages so that the voltage lied along the d axis; synchronization was achieved by means of a Phase-Locked Loop (PLL), which forced v q = 0. The reference values of the current components along the d and q axes were determined by, respectively, the active and reactive power control loops. The set-point active power control loop was determined by the DC voltage controller, which regulated the voltage at the PV generator busbar to the optimal value determined by a Maximum Power Point Tracking (MPPT) algorithm. The set-point Q ci of the ith reactive power control loop was determined by the local voltage droop controller. More details about the VSI control architecture were reported in [23].  In the following, firstly, the design of the droop constants k ci according to the proposed procedure is presented, showing the main results concerning the stability of the voltage control systems. Then, some results giving evidence of the performance of the designed voltage droop controllers are presented in various operating conditions of the distribution system and of the PV systems.

Design of the Voltage Droop Controllers
According to the proposed procedure, in the first step, the overall system model P(s) = Θ G(s) (see (5)) was determined. To this aim, the matrix Θ in (3) was evaluated by using the linear method in [28]. It was determined for various operating conditions of the distribution system and of the PV systems, obtained by combining the following variations: • the voltage of the slack node was varied in the range 0.95-1.05 p.u., • the loads were varied in the range 30-90% of their rated values, • the active powers injected by the six PV systems were varied in the range 0-20 kW.
Then, a specific matrix Θ 0 was chosen for the SISO design, in particular the one characterized by the unitary voltage for the slack bus, loads at 70% of their rated values, and 15 kW active powers injected by the six PV systems. Concerning G(s) in (4), the same transfer function was assumed for all six PV systems, in particular: with T = 0.020 s and τ = 0.035 s, being identified by analyzing the response of the reactive power Q DG i to a step variation of the set-point Q c i .
In the second step of the procedure, the SISO design of the voltage droop constants k ci was performed. Considering each PV system at a time, its droop constant was determined by imposing a phase margin of 40 • at the crossover frequency of the open loop transfer function k ci θ 0ii G i (s). The resulting values for the droop constants are reported in the first row of Table 2. The SISO design was also tested by numerical simulations. To this aim, the distribution system was simulated in various operating conditions with only one PV system connected at a time, verifying its stability. As an example, let the test of the PV system #1 be considered. Concerning the operating conditions, the loads and the voltage of the slack node were assumed equal to their minimum values and the active power injected by the PV system to its maximum value. Figure 6a reports the time evolution of the voltage at node #6, at which the single PV system #1 was connected, when the load at node #9 was disconnected at time instant 10 s.   Table 2, (b) with all six PV systems connected to the network and the chosen values for k ci ; see Table 2, (c) with all six PV systems connected to the network and doubled values for k ci with respect to (b). The third step of the procedure consisted of evaluating the maximum values of k ci imposed by (12). The obtained upper limits are reported in the second row of Table 2. Comparing the first and the second row of Table 2, it is apparent that the values obtained from the SISO design were approximately five times larger than the limits imposed by (12). Since such limits represent sufficient yet not necessary conditions, it can be stated that adopting the SISO design may not guarantee a stable behavior of the voltage droop controllers when all six PV systems were connected to the distribution system.
The final step of the proposed procedure was to choose, for each PV system, the value of k ci , which was obtained by the SISO design subject to the criterion (12). The chosen values are reported in the third row of Table 2 and guaranteed a stable interaction among all six voltage droop controllers. For the sake of comparison, Figure 6b) reports the time evolution of the voltage at node #6 assuming the same operating conditions of the distribution system as the ones of Figure 6a), but connecting all six PV systems with the chosen constants for the voltage droop controllers.
An ultimate issue concerned how much margin from instability was guaranteed by the proposed procedure. In this specific case study, underdamped oscillations arose when the values of k ci were increased up to about two times the chosen values (see Table 2), denoting that the stability limit was being reached. It is evident from Figure 6c, which reports the same transient as the one in Figure 6b, but doubling the values of k ci for all six connected PV systems. It is important to underline that in this case, underdamped oscillations appeared for values of the droop constants that were much smaller (about 2/5) than the corresponding values obtained by the SISO design. Obviously, a further increase of the values of k ci brought instability to the system.

Performance of the Designed Voltage Droop Controllers
To assess the performance of the chosen voltage droop constants, various numerical simulations were performed. In the following, the results are reported referring to the following two simulations: 1.
the voltage of the slack bus was set to 1.0 p.u. for the whole simulation; the solar irradiation was equal to 800 W/m 2 at the beginning and changed to 1000 W/m 2 at time instant 6 s; an additional 15 kW/7.5 kVAr load was connected to Node #8 along Feeder 1 at time instant 1 s and disconnected at time instant 11 s; 2.
the voltage of the slack bus was set to 1.0 p.u. at the beginning of the simulation and was subject to a step reduction to 0.99 p.u. at time instant 1 s; the solar irradiation was equal to 800 W/m 2 for all the simulations; the PV System #2 was disconnected from Node #11 at time instant 11 s.
In both simulations, the loads were equal to 70% of their rated values, and the reference values for the voltage droop controllers were fixed equal to 1.0 p.u. Figures 7 and 8 report, for the first simulation, the time evolution of the voltage amplitudes at the network node and of the corresponding PV reactive powers (with the same color), respectively for Feeder 1 and Feeder 2. The graphs give evidence of the stable operation of the voltage droop controllers and of their support to the system voltage regulation. In particular, the transients due to the load connection (disconnection) at time instant 1 s (11 s) caused large voltage variations along Feeder 1, where the load varied; in response, the voltage droop controllers of the PV systems along Feeder 1 were large and counteracted the voltage variations. On the contrary, along Feeder 2, the variations of the voltage amplitudes and of the reactive powers were smaller and opposite in sign with respect to the ones along Feeder 1; the reason was that the two feeders were coupled by the voltage at the LV busbar in the supplying substation, which was affected more by the reactive power variations of the PV systems along Feeder 1 than by the load variations. At time instant 6 s, the increase of the solar radiation caused an increase of the active powers injected by all PV systems. Consequently, all the voltage amplitudes at the network node increased, and the voltage droop controls reacted with a reactive power reduction.  Similarly, Figures 9 and 10 report, for the second simulation, the time evolution of the same quantities as the ones in Figures 7 and 8. In this simulation as well, the graphs give evidence of the stable operation and of the support to voltage regulation by PV systems. In particular, at time instant 1 s, the step decrease of the voltage at the slack bus caused a sudden reduction of the voltage amplitudes at all the PCCs; the voltage droop controllers reacted by increasing the reactive power of all the PV systems. At time instant 11 s, the loss of PV System #2 caused a sudden reduction of the voltages along Feeder 1 and the consequent reaction by the voltage droop controller of PV Systems #1 and #3, which significantly increased their reactive powers. The combined effects of the perturbation and of the controllers' actions caused small perturbations on the voltage amplitudes along Feeder 2.

Conclusions
The voltage-reactive power droop controllers of DG systems are commonly used because of their simplicity and are typically available in commercial DG systems as the requirement of international standards. A simple constraint was proposed, which guaranteed that the droop constants, designed according to a SISO technique, preserved the stability of the active distribution system. The constraint was derived using the concept of the interaction measure under decentralized control. The proposed procedure was applied to a case study with six DG systems, and the results of numerical simulations gave evidence of its effectiveness. In future work, the proposed approach will be: (i) applied to voltage-active power droop controllers used in the presence of energy storage systems; (ii) extended to large active distribution networks with a huge number of DG systems requiring approaches based on voltage control zones.
Author Contributions: Conceptualization, G.F. and M.R.; methodology, G.F. and M.R.; formal analysis, G.F. and M.R.; validation, G.F. and M.R.; writing, G.F. and M.R. All authors have read and agreed to the published version of the manuscript.
Funding: This research was funded by the Italian Ministry of University and Research by the special grant "Dipartimenti di eccellenza".

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