Sensitivity Analysis for Voltage Stability Considering Voltage Dependent Characteristics of Loads and DGs

The sensitivity analysis becomes particularly critical for voltage stability analysis due to the fluctuation in power outputs of renewable energy resources. Besides, impacts of different load modeling and the operation mode of Distributed Generations (DGs) are not addressed in the well-known sensitivity analysis methods. Therefore, this work presents a new sensitivity analysis approach to find the relation between the Voltage Stability Margin (VSM) and the control variables of power systems, considering the voltage dependent characteristics of loads and DGs. The sensitivity analysis is performed on VSM, defined from equivalent nodal analysis, via its differential equation. To include the voltage dependent characteristics, loads are modeled as polynomial function (ZIP model) and DGs are considered to be operated with constant current and constant power modes. Based on this analysis, the sensitivity of VSM can be directly obtained by taking the derivatives of nodal voltages with respect to control variables. The validity of the developed approach is demonstrated on the IEEE 118 bus system.


I. INTRODUCTION
The continuous growth in power demands and renewable distributed generation have posed challenges on voltage stability of large-scale power networks. Voltage instability or collapse makes the system insecure and may cause blackouts [1]. To keep the system secure, a sufficient loading margin for network stability is required. Thus, Loading Margin (LM) has to be monitored closely to ensure voltage stability. It gives an indicator about the amount of additional power consumption that would cause a voltage collapse [2].
Preventive or corrective control is normally used to mitigate voltage instability. Voltage stability control methods mainly depend on the relationship between LM and control variables. Sensitivity analysis is usually utilized for this pur- The associate editor coordinating the review of this manuscript and approving it for publication was Tariq Masood . pose to manage control variables for voltage instability mitigation. Various sensitivity approaches have been proposed for voltage stability analysis.
A LM sensitivity method based on Continuation Power Flow (CPF) is developed in [3] for load shedding in power systems. In [4], LM sensitivity based on the optimal multiplier power flow is used for load shedding purposes. A method depending on multi contingency sensitivity using CPF is proposed for preventive stability control [5]. In [6], the sensitivity analysis method depends on singularity of Jacobian matrix (J) for voltage stability control. Linear and quadratic approximation approach is used in [2] to show the relationship between LM and network parameters (or controls). This approach requires to find the eigenvalues of the Jacobian at maximum loadability point. The voltage stability control developed in [7] uses LM sensitivity analysis based on a complex modeling of systems components and VOLUME 9, 2021 This work is licensed under a Creative Commons Attribution 4.0 License. For more information, see https://creativecommons.org/licenses/by/4.0/ depending on the equilibrium tracing method to find the maximum loadability point. In [8], sensitivity calculation is used with CPF for voltage stability analysis to take into account the low voltage conditions. In [9], tangent vector method is developed to trace PV and PQ curves. In [10], linear sensitivities are combined with eigenvalue analysis for voltage contingency ranking. A sensitivity method based on modal analytical approach for LM is presented in [11]. However, the aforementioned methods depend on CPF via iterative process [12] or on the singularity of J to find the critical point. They require remarkable computational efforts and thus they are not suitable for online applications of modern networks hosting renewable energy resources. Indeed, the fast response of FACTS devices and DG units can change the operation conditions of modern power systems [13].
In [14], a probabilistic LM sensitivity method based on bootstrap technique is developed. A trajectory sensitivity approach for sensitivity analysis considering load uncertainty is proposed in [15]. This method is based on the base case condition without capturing the system changes. To provide fast sensitivity analysis, a new LM sensitivity approach based on Look-Ahead method is proposed in [16], [17]. However, these methods suffer from inaccuracy.
Measurement-based LM sensitivity methods can avoid the complex calculation and reduce the computational time. In the literature, many methods are developed for voltage stability analysis as: VSI-index [18], Voltage Instability Predictor [19], Tellegen's theorem [20], decision tree [21], VCPI Indicator [22], Corsi Identification Algorithm [23], L-index [24], [25], coupled single-port circuit concept [26]. However, such works have no information to provide about the decisions of preventive control. The information of which control variables to use for preventive or corrective voltage stability control is equally as important as the voltage stability assessment itself. The authors in [27], [28] have recently performed sensitivity analysis on voltage stability indices derived from the Thevenin theorem. The assumption that the collapse occurs when load impedance equals the equivalent impedance of a load bus makes the analysis not enough rigorous. Besides, they modeled the loads as constant impedances through sensitivity analysis, which is not accurate for practical power systems. Indeed, the voltage dependent characteristics of loads have an important impact on load power consumption. The sensitivity coefficients are normally obtained by approximating the change in the voltages and power flows with respect to changes of preventive controls. This means that such sensitivities suffer from inaccuracy because they are not able to include the variation of load powers with voltage [29]. Furthermore, such characteristics of DGs can also affect LM sensitivity analysis according to their operation modes. Since voltage dependent characteristics of loads and DGs play an inevitable role in LM sensitivity analysis, they would definitely affect the preventive or corrective control.
In this regard, a new sensitivity analysis approach considering the voltage dependent characteristics of loads and DGs is presented. The sensitivity analysis is performed on VSM, defined from equivalent nodal concept, via its differential equation. The polynomial model of loads and the constant current and constant power modes of DGs are used in the analysis to reflect the voltage dependent characteristics. The analytical model for sensitivity analysis is developed as a function of the derivatives of nodal voltages with respect to control variables.
Other advantages are also associated with the proposed method. By this method, the sensitivities to control variables can be updated to take into account the continuous and fast variations caused by DGs and fast-response devices. The proposed method is fast and requires less computation time and thus it can meet the requirements of smart grid applications, especially in the context of optimization techniques. Compared to the well-known techniques, the proposed approach does not require iterative process and there is no need to calculate the critical point via CPF or the singularity of the Jacobian. The proposed method can also accurately guide the network operators to rank the control variables and provide the most effective ones.
The key contributions of this work are: (a) to include the voltage dependent characteristics of loads and DGs; (b) To the best of our knowledge, sensitivity analysis is performed for the first time on VSM, defined from nodal equivalent concept, without any approximation. The sensitivity of VSM is obtained as a function of the derivatives of nodal voltages with respect to control variables.
The rest of this work is organized as follows. Section II presents equivalent nodal-based VSM. Section III presents the sensitivity analysis method. Extension to include the voltage dependent characteristics of loads and DGs is presented in section IV. Section V shows simulation results and Section VI states the assumptions and limitations.

II. VOLTAGE STABILITY MARGIN (VSM)
For multi-bus power network, the node voltages can be written in terms of line impedances and node currents as: where Z is the system impedance matrix. V and I are vectors of bus voltages and currents, respectively. The system nodes are usually divided into generator buses G and load buses L. Assuming that G = {1, . . . , M } and R = {M + 1, . . . , N }, the relation in (1) can be replaced by (2), as shown at the bottom of the next page.
From (2), the load voltage at bus i R can be found as: Equation (4) is equivalent to: where The terms V eq,i , I eq,i , and Z eq,i are the equivalent voltage, current and impedance referred to node ''i'' respectively.
where * denotes conjugate. v i is the voltage magnitude at bus ''i''. By defining V i I * eq,i = S eq,i as equivalent power of node i, we obtain: The formula illustrated in (7) is similar to the formula developed in [24] for defining L-index. According to [24], [25], a voltage stability index, , for any load bus i R can be written as: where Y ii = Z −1 ii is the self-shunt admittance of bus ''i''. The value of can vary within the range [0,1]. The instability point occurs at = 1.0. The formula expressed in (8) can be rewritten as: It is clear from (9) that the index of a particular bus depends on its own voltage, load currents, and system impedance matrix. The value of can be considered to represent the voltage stability margin. The greatest index value among network bus indices can be used to represent the system LM [26].

III. VSM-SENSITIVITY ANALYSIS A. SENSITIVITY OF VSM TO CONTROL VARIABLES
The dependence between VSM and control variables can be explained by considering the small power network presented in Fig. 1. The system has 4 load buses and one generator bus. It is assumed that a control variable (marked in red color) is available at bus 3. Fig. 1a represents the current system state before performing control actions. Fig. 1b shows the system state after performing control actions. Any action (positive or negative) made by the control variable will result in changing the nodal voltages and currents by V, I respectively. By referring to (9), we can see that the index of any bus will then change by . The control action at any specific node will not only affect VSM of its own node, but also the VSM of other nodes.
If it is assumed that i = |f |, the sensitivity of the index of any bus to the change in the control variable ''u x '' can be then obtained as (where f is a function in Cartesian coordinates): The can be extracted from the derivative of the term f with respect to the control variable u x as: df Taking the derivative of (9) with respect to the control variable ''u x '', we obtain: It is clear form (11) and (12) are functions of sensitivity of load voltages and currents with respect to control variable u x .

B. SENSITIVITY OF LOAD CURRENTS TO CONTROL VARIABLES
The load current I * j can be expressed in terms of complex power and node voltage as: Taking the derivative of (13) with respect to u x , one can obtain: At constant load impedance, the terms S j and dS j du x can be calculated as: Substituting (15) and (16) into (14), we obtain: It is clear from (17) that the current sensitivity can be expressed in terms of voltage sensitivity. By referring to (10), (11), (12) and (17), we can then conclude that the index sensitivity d i du x can be expressed in terms of voltage sensitivity.

C. SENSITIVITY OF LOAD VOLTAGES TO CONTROL VARIABLES
The dependence between network load voltages and control variables can be found using conventional power flow equations. Let us consider the system power flow equations at fixed operating condition represented as [30]: where y is the network states (voltage magnitudes and angles) and u denotes for control variables. Equation (18) assumes that the only disturbance in the system is the control variables u. If the conventional power flow equation is expanded into Taylor series, we obtain (assuming that high-order terms at the operating condition (y o , u o ) are ignored): where y and u are the change in the state and the preventive control, respectively. H y (y o , u o ) and H u (y o , u o ) are the partial derivatives of power injections with respect to node voltages and preventive controls, respectively. Thus, the sensitivity of the operating states to control variables can be calculated as: The term H y (y o , u o ) represents the Jacobian matrix while H u (y o , u o ) represents the sensitivity of power changes to control variables.
It is worth mentioning that the control variables can be the generators' terminal voltage, capacitor switching, tap changer transformers and reactive power injections by DGs.

IV. SENSITIVITY ANALYSIS INCLUDING VOLTAGE DEPENDENT CHARACTERISTICS A. EXTENSION TO LOAD MODELS
Many models can be used to represent the loads. However, VSM and the sensitivity analysis can be extended to the polynomial load models (ZIP model) illustrated in (21). The model shows that the active P and reactive Q powers vary with the magnitude voltage v.
where P 0 and Q 0 are specified active and reactive power of the loads. v 0 is a specified voltage of the system nodes. P P , I P , Z P and Q q , I q , Z q are constant parameters of the model and satisfy the following equations: It is clear from (21) and (22) that the polynomial model (ZIP model) represents a combination of loads (constant power PQ, constant current I, and constant impedance Z). They are modeled by assigning a percentage of the total load to each of the three aforementioned load models.
To include the ZIP model in sensitivity analysis, (13) can be replaced by: where P j and Q j are the real and reactive power of load ''j''. Taking the derivative of (23) with respect to u x , one can obtain: The derivative of (21) with respect to u x is: Substituting the derivatives dP du x and dQ du x referred to load ''j'' into (24), we obtain: It is clear from (26) that the current sensitivity can be expressed in terms of the voltage sensitivity to control variables. For constant power mode, the active and reactive power output by DGs are P CP p and Q CP q , respectively. For constant current mode, the active and reactive power output by DGs are I CC p v and I CC q v, respectively. The terms P CP G , Q CP G , I CC p and I CC q are constants. Accordingly, the active and reactive power generation (P DG,j and Q DG,j ) at bus ''j'' can be expressed as: If the DG unit installed at bus ''j'' is operated at constant current mode, P CP p,j and Q CP q,j equal zero. If the DG operated at constant current mode, I CC q,j and I CC q,j equal zero. If no DG is installed at bus ''j'', the generation powers P DG,j and Q DG,j equal zero.
If it is assumed that the load bus ''j'' in (2) has a DG unit, the current I * j illustrated in (9) can be then replaced by: where the subletters L and DG denote the load and DG unit, respectively. Taking the derivative of (28) with respect to u x we obtain: The derivative of (27) with respect to u x is: Substituting (16) and (30) into (29), we obtain: It is also clear from (30) and (31), that the overall load current sensitivity can be expressed in terms of the voltage sensitivity.
To include both effects (i.e. effect of load models and DG modes) in the analysis, (28) and (29) can be written as: Substituting (25) and (30) into (33), we obtain: (34), as shown at the bottom of the next page.
A Flowchart of the proposed sensitivity method is presented in Fig. 2.

V. SIMULATION RESULTS
To examine the validity of the developed sensitivity method, the IEEE 118 bus system shown in Fig. 3 is used in this study. The network data and parameters can be found in [31]. Different scenarios (normal and contingency scenarios) are used in this work. The pilot bus, which has the smallest VSM is used as a representative bus for the analysis. To cover all possible controls, the analysis discusses the sensitivity of VSM with respect to active and reactive power injections at network buses, not only to preventive controls. Indeed, most of preventive controls act as power injections (i.e. generators' active power, power generation of DGs, shunt capacitors).

A. COMPARISON WITH CLASSICAL METHODS
To examine the accuracy of the proposed sensitivity method, the numerical values of the obtained sensitivities of pilot bus with respect to active and reactive power injections are compared with the results obtained using linear sensitivity analysis [2]. To make the system more stressed, the loads are multiplied by 1.7 of their nominal values. The comparison is investigated at two scenarios: at no contingency and at line 98-100 outage. Fig. 4 shows the relative errors in the sensitivities obtained using the proposed method and the results of linear sensitivity analysis.
The results in Fig. 4 shows that the relative errors in the sensitivities (at all conditions) are very small. Besides, the sensitivity ranking with respect to power injections is kept unchanged, compared to the sensitivity ranking obtained using the linear sensitivity analysis. This provides a demonstration on the accuracy of the proposed method for VSMsensitivity analysis.

B. IMPACT OF LOAD MODELS
To study impact of voltage dependent characteristics of loads on system VSM, five load models are considered in this work. Model 0 is the base model and assumes that the loads are represented as constant impedances. The other models are assumed as shown in Table 1.
Remark: For comparison purposes and without loss of generality, this study considers five different cases of the polynomial coefficients of load modeling. However in practical systems, the polynomial coefficients can be obtained for a very good fit to the measured data (i.e. using a curve-fitting procedure).
To validate the proposed method under different scenarios, network loads are multiplied by 1, 1.7, and 2 of their nominal values. Besides, one contingency (line 98-100 outage) is also considered for the analysis.
The sensitivity of VSM of pilot bus with respect to active and reactive power injections by considering the five cases of load modeling is calculated. Sensitivity results to active power injections during different loading conditions at no contingency and during line 98-100 outage are presented in Fig. 5 and Fig. 6, respectively. Sensitivity results to reactive power injections at no contingency and during line 98-100 outage are presented in presented in Fig. 7 and Fig. 8, respectively.
It is worth mentioning that the nodes with very high sensitivities are not shown in the Figures 5-8 to improve the visibility of the comparison. From these figures, it is clear     that by considering the voltage dependent characteristics of loads (i.e. the polynomial modeling), the sensitivity of VSM VOLUME 9, 2021    to active and reactive power injections can vary from the one obtained using constant load impedance (base case). Some variations can be small while others can be larger. It is generally clear that as load level increases, the difference in the sensitivity also increases. Besides, the difference in the sensitivity with respect to active power injections is greater than those with respect to reactive power injections.
To classify the differences in the sensitivities, Histograms of relative errors between the sensitivities obtained at base model (constant load impedance) and the ones obtained using the other load models (models 1-4) during different loading conditions are carried out. The relative error can be found as: where δ is the sensitivity obtained by assuming constant load impedance while δ is the sensitivity obtained using other load models. Histograms of relative errors in the sensitivity to active power injections at no contingency and during line 98-100 outage are shown in Fig. 9 and Fig. 10, respectively. Histograms of relative errors in the sensitivity to reactive power injections at no contingency and during line 98-100 outage are shown in Fig. 11 and Fig. 12, respectively. From the figures 9-12, it is clear that the errors in the sensitivities are less than 0.2 in case of normal loading conditions. As load demands increase, the sensitivity errors increase. The relative errors of some sensitivities can reach for more than 1 (i.e. 100%). For some load models, around 80% of the sensitivities to active power injections and around 60% of the sensitivities to reactive power injections have generally a relative error more than 1. This high normalized density occurs in case of using constant current load models (model 1) which assumes that no constant impedance loads are available. This validates the accuracy of the proposed method and demonstrates the need to consider the voltage dependent characteristics on loads during sensitivity analysis of VSM.
By computing the average relative error for each load model, we can generally notice that as load level increases, the average relative error is also increase. This notice is clear in the average relative errors of the sensitivity to both active and reactive power injections and for all scenarios (no  contingency and line outage scenarios). Refer to Tables 2 and 3, which summarize the average relative errors.
From Tables 2 and 3, it is clear that the average relative errors of the sensitivity to active power injections are greater than the sensitivity to reactive power injections during different loading conditions. For active power injections, we can also notice that the average relative errors in case of line outage are greater than the average relative errors in case of no contingency. In contrast, the average relative errors at both VOLUME 9, 2021  cases (line outage and no contingency) are almost close for reactive power injections.
From both tables, we can also conclude that the maximum average errors occur at load model 1 (i.e. constant current model). This means that when I P orI q = 1, the loads have a higher effect on the sensitivity.
From this comparison, it is clear that the voltage dependent characteristics of loads can significantly affect the sensitivity analysis of VSM. Therefore, considering such characteristics in the sensitivity analysis is necessary. The variation of the average relative error in the sensitivity of VSM of pilot bus as a function of Z P and Z q can be shown in Fig. 13 and Fig. 14, respectively. It is clear that the maximum average relative error occurs when Z P or Z q = 0 (i.e. I P or I q = 1). As Z P or Z q starts to increase, the average relative error will decrease to reach zero at Z P or Z q = 1, which represent the base model. This characteristic is clear for all scenarios, which also proves the validity of the proposed method. A negative slope of the curves can be obtained by taking the average relative errors as a function of I P or I q .
The variation of the average sensitivity of VSM of pilot bus as a function of Z P or Z q can be shown in Fig. 15 and Fig. 16, respectively. A linear function is fitted to each average sensitivity curve along all the scenarios. It is clear that the value of Z P or Z q affects the average sensitivity. The significance of the effect of Z P or Z q on the sensitivity depends on the system operating condition. Moreover, the curves show that the average sensitivity increases as the value of Z P or Z q increases. It can also be seen that the scenario with 2 of nominal loading conditions is more sensitive to the value of Z P or Z q compared to the other loading conditions. Besides, the value of Z P or Z q has a higher impact on the average sensitivities with respect to active power compared to the ones with respect to reactive power injections.

C. IMPACT OF DG OPERATION MODES
To study the impact of the voltage dependent characteristics of DGs on sensitivity analysis of system VSM, it is assumed that 12 DGs are installed at the buses: 3,11,29,35,45,53,60,78,88,98,106, and 118 (refer to Fig. 3). The rating for each DG is 70% of the load power of its own bus. To make the system more stressed, the load demands are multiplied by 2. Besides, it is assumed that the loads are modeled as constant  impedances. The validation for this section is achieved at no contingency and under one line outage.
The sensitivity of VSM of pilot bus with respect to active and reactive power injections is calculated at both constant power mode and constant current mode. The coefficients I CC p , I CC p = 1 and P CP p,j , Q CP q,j = 0 for constant current mode. In case of constant power mode, the coefficients I CC p , I CC p = 0 and P CP p,j , Q CP q,j = 1. From sensitivity calculations, it is noticed  that by considering the voltage dependent characteristics of DG units, the sensitivity of VSM to active and reactive power injections can vary. The relative errors between the sensitivities obtained by considering the two modes are also carried out. Histograms of relative errors of the sensitivity with respect to active power injections at no contingency and during line outage are shown in Fig. 17. Histograms of relative errors of the sensitivity with respect to reactive power injections at no contingency and during line outage are shown in Fig. 18. From these figures, it is clear that the relative error of the sensitivity to both active and reactive power injections and at all scenarios can reach more than 1. The maximum normalized density of the errors generally occurs when the relative errors are in the order of (0.6 -0.8) or (0.8 -1). This also demonstrates the necessity to include the voltage dependent characteristics of DGs in sensitivity analysis. By computing the average relative errors, we can notice that the average relative errors of the sensitivity to active power injections are also greater than the sensitivity to reactive power injections. It is also noticed that the average relative errors in case of line outage is greater than the average relative errors in case of no contingency.
From this comparison, we can conclude that the voltage dependent characteristics of DGs has a significant impact on the sensitivity analysis of system VSM. Therefore, considering such characteristics in the sensitivity analysis is also necessary.

D. IMPACT ON VOLTAGE STABILITY MARGIN
Impact of the proposed method considering the different load models on VSM is investigated in this section. The sensitivities obtained in section V-A for each load model are used to obtain the margin of pilot bus. The margin can be obtained by multiplying the sensitivity of to control variable (i.e. d du x ) with the amount of control change (i.e. u x ) as: The control change u x is assumed to be 1.0 p.u of active and reactive power injections at pilot bus. The sensitivity d du x are changed based on the load models. The relative errors between the margin obtained using the base mode and the ones obtained using other models are then calculated. The percent relative errors in the margin of pilot bus due to active and reactive power injections are illustrated in Table 4 and Table 5, respectively. Although some errors are small and can be acceptable for power system operators, other errors are high and cannot be ignored. Some percent relative errors exceed 6.7% for the study system, which may be higher for practical larger systems. As a result, ignoring the load models can significantly overestimate/underestimate the voltage stability margin of power system.

VI. ASSUMPTIONS AND LIMITATIONS
To include the voltage dependent characteristics of loads and DGs for voltage stability analysis, ZIP model is used in this work. However, there are some assumptions and limitations for the proposed method.

Assumptions:
In real networks, the voltage dependent characteristics of loads are not known. Such characteristics can be obtained by using a polynomial model or the exponential model of loads. These models are accurate enough to model the load characteristics and can help engineers to fit a curve to the measured data to represent the voltage dependence of power demand.

Limitations:
The polynomial coefficients can be obtained for a very good fit to the measured data. This can be achieved using a curve-fitting procedure for load power to minimize the error between the fitted approximation and the measurement values of voltage points. Such approach of determining the polynomial coefficients requires high effort and thus it represents a limitation for the proposed method. This limitation can be mitigated by considering a finite number of voltage points within a range of operating voltages (for example, 15%).
In the context of smart grids, the advanced measuring infrastructure (AMI) can be used to collect extensive measurements of the power demand.
Moreover, the voltage stability margin developed in this work is obtained through a centralized approach. It is wellknown that centralized methods may increase the burden of computation and communication. A decentralized-based voltage stability margin can be used to mitigate such problem. However, the same procedure discussed in this work can be performed to consider the voltage dependent characteristics of loads and DGs.

VII. CONCLUSION
This work considers the voltage dependent characteristics of loads and DGs in the sensitivity analysis of VSM of power systems. The method does not require iterative process and directly depends on the derivatives of nodal voltages with respect to control variables.
Simulation is conducted on the IEEE 118 bus system under different operating conditions. The voltage dependent characteristics can be specified based on the polynomial coefficients. Results show that the voltage dependent characteristics of load and DGs can significantly affect the sensitivity values. It is clear that the relative errors between the sensitivities can reach more than 100%. The Histograms of relative errors also show that up to 80% of sensitivities have a relative error more than 100%. This demonstrates the necessary to pay attention to load and DG modeling in the sensitivity analysis. The results also show that as the load level increases, the average relative error also increases. Besides, the average relative errors of the sensitivity to active power injections are greater than the sensitivity to reactive power injections. Furthermore, as the value of Z P or Z q of load modeling increases, the load has a higher effect on the sensitivity and a lower impact of the relative error.
Impact of voltage dependent characteristic on VSM is also studied. The results show that the percent relative errors due to load modeling can exceed 6.7% for the study system. Ignoring the load and DGs models can significantly overestimate/ underestimate the VSM of power system.
Our future work is to develop a decentralized sensitivity analysis approach for voltage stability assessment. Considering different types of load models for voltage stability will also be studied.