Influence of parametric uncertainties and their interactions on small-signal stability: A case example of parallel-connected active loads in a DC microgrid

Classical stability analysis techniques based on nominal models do not consider the uncertainty of system parameters, their interactions, and nonlinearity, which are important characteristics of practical highly coupled microgrids. In this work, variance-based sensitivity analysis is used to identify parameter combinations that have a significant impact on the smallsignal stability of a microgrid featuring two parallel active loads. The analysis indicates that the effectiveness of source-side damping is reduced when resonant frequencies of load input filters become matched. Further results using derivative-based sensitivity analysis reveal that source-side resistance can exhibit drastically different effects on the stability if load input filter resonant frequencies are matched with respect to the case when they are well separated. These behaviours are verified using time-domain switching models.


Introduction
In small-signal modelling of conventional electrical power systems (EPSs), network dynamics are generally neglected since their time constants are much smaller than those of synchronous generators and their control loops. In conventional EPSs, the effects of individual loads are insignificant (due to their small capacity), and thus loads can be aggregated and collectively modelled, whereas in microgrid (MG), individual loads of EPSs play a critical role in system stability and have to be included in stability analysis [1].
Classical techniques to identify dynamic interactions through the use of participation factors (PFs) in EPSs [1,2] are based on the analysis of nominal models. The work reported in [1] provides important insights into the complexity of coupling phenomena between the sources, loads, and network. Although these techniques are prevailing, they are insufficient for describing how the interactions change under possible variations of parameters. For example, the recent work [3] identifies a coupled resonance phenomenon between multiple parallel load-side filters in an MG. The observed filter coupling effects depend strongly on parameter values. This observation highlights the necessity of modelling system uncertainties, such as the uncertainties in the length of lines between different nodes in an MG. Similarly, Wu and Lu [4] identified mutual coupling phenomena between multiple parallel load-side filters that must be considered at the design stage of the proposed stabilising controller. Also, Wan et al. [5] illustrated reduction in stability margins due to coupling effects from between transmission line impedances and the source and multiple loads, including impedance in-between loads in an AC system. Recent publication [6] quantifies effects of multiple simultaneous parametric uncertainties on the stability of a power electronics system by applying the structured singular value-based stability analysis method.
To identify and analyse this type of dynamic phenomena, Luo and Ajjarapu [7] proposed tracing eigenvalues over parametric changes to observe the movement of oscillatory modes. This work shows that modal resonance effects in power systems, such as subsynchronous resonance, can be caused by the interaction between two distinct modes, i.e. between torsional modes and subelectrical modes, and can be viewed as a precursor to system instability. For problems with numerous multiple uncertainties, exhaustive testing of all possible combinations quickly becomes computationally infeasible. As a result, there is an increasing interest in applying global sensitivity analysis (SA) techniques which aim to quantify the behaviour of an output over the entire range of uncertainties. Probabilistic methods can be used to examine the scenarios from the whole range and are easily implemented using Monte Carlo (MC) techniques. Global SA techniques are reviewed in [8] and some considered for the importance ranking of uncertain parameters in power networks in [9,10]. Of interest are those that can help identify different interactions between system parameters, such as variance-based sensitivity analysis (VBSA) [11].
The main idea of this paper is to: (i) demonstrate the application of VBSA to MG EPSs under uncertainty to identify the critical parameters that exhibit nonlinear effects through their interactions on the small-signal stability; (ii) using the identified parameters from (1), present local SA results for selected uncertain system scenarios. Interactions are shown by critical eigenvalue sensitivity to distribution line resistance for different load conditions. The paper is organised as follows. Section 2 describes the model of the system featuring two parallel active loads. Two unique modelling cases considered are as follows: load input filter resonant frequencies are matched, and load input filter resonant frequencies are separated. In Section 3, the VBSA technique is defined and applied to the analysed system to identify interaction effects among multiple uncertain system parameters. Section 4 performs derivative-based SA to identify effects of line resistance on the stability. Section 5 demonstrates these behaviours via switching model simulations, and Section 6 summaries main findings of the paper. Fig. 1 shows the configuration of a distribution system with two active loads sharing a common voltage bus. The source is modelled as an ideal voltage source with a finite line resistance (R line ). Each active load is modelled as a tightly regulated synchronous buck converter with an LC low-pass input filter, as shown in Fig. 2 [12]. The stability analysis of each active load subsystem individually is inapplicable in the presence of other loads in parallel. R line can not only provide a damping effect for the system under study [4], but it can also contribute to the coupling between the loads, as shown in [3].

System structure
A single linearised model based on the nominal operating condition is not sufficient to describe system dynamics under uncertainty; therefore, it is necessary to have an analytical nonlinear system model through parameterisation. This is modelled using a set of differential equations f and output equations : where x are state variables, u are system inputs, and p are all other parameters describing the system. For conciseness in this paper, a k-number of parallel active loads can be represented by deriving a set of general equations based on a single active load in Fig. 2. Firstly, let the parameters of each subsystem of an active load be denoted by subscripts, i.e. the active load i is formed by combining the filter F i , averaged buck converter B i , controller C i and load R i . By applying the circuit averaging technique and Kirchhoff laws to the dynamic models in Figs. 1 and 2, each active load has the differential equations in (2)-(8) and bus voltage in (9). To note, (3) shows the sharing of the common bus through the i L,Fi state variable, i.e. R line contributes to a voltage drop of V bus based on the sum of the currents flowing into each active load: (see (5)) For the two active load system in Fig. 1, this results in a state space model of 14th order, such that Two case studies are considered for the system in Fig. 1 Fig. 2, are identical for both active loads, with parameters shown in (12). For buck converter switching frequency of 50 kHz and under nominal parameters, the regulator voltage loop gain is designed to a 1 kHz crossover frequency: low-side switch (Q2) on state resistances 5 m

Variance-based sensitivity analysis (VBSA)
The main idea of VBSA is the quantification of the influence of uncertain input variables X = [X 1 , X 2 ,..., X n ] on the model output Y = f (X), without the assumption of model linearity. This model can be approximated as a combinatorial expansion of 2 n additive components, known as the high-dimensional model representation [13]: where the component f i (X i ) represents the contribution to Y by the independent variation of the ith variable X i the second-order term f ij (X i , X j ) represents the contributions to Y due to the interactions between inputs X i and X j , not attributable to either f i (X i ) or f j (X j ), and so on. Each component is orthogonal to each other. This can further be decomposed in terms of the variances of the output quantity and normalised [14], such that where each S term is an associated sensitivity measure of input parameter(s) denoted by the subscript(s). The first-order sensitivity index, S i , or main effect is an indicator of the significance of X i on the variance of Y, and can be calculated using [11] where the total variance of Y under all variations is denoted by and the variance of the conditional expectation is: Homma and Saltelli [11] introduced the total effect indices, S Ti , that measure the effect of single variable (X i ) including its higher order effects with other variables: where X ∼i is a set that includes all variables except X i . Both main effects S i and total effects S Ti are computed through the estimation of the multi-dimensional integrals using Monte Carlo integration [14,11]. Estimation of the variance of conditional expectations is achieved by keeping certain input variables common, for example, S i is calculated by keeping X i constant while sampling all other variables independently as indicated by (17). The difference between the total effect indices and the first-order indices (S Ti − S i ) accounts for the interactions of parameter X i with other variables.

Modelling for VBSA
The stability of the system can be assessed from the real part of the eigenvalues of the linearised (small-signal) model. Due to the parametric uncertainty, the steady-state operating point (x 0 ) and system eigenvalues are also uncertain. Eigenvalues must be numerically computed and are dependent on the uncertain input variables. This process can be considered as a model input-output mapping denoted by: f: X Y, where X is the vector of uncertain input parameters and Y is the real-part of the eigenvalue(s). Random sampling from the probability distributions of the uncertain variables generates a unique scenario. where Matrix A is an n × n matrix which has a set of n distinct eigenvalues ( 1 , 2 , …, n ), representing the modal response of the system, and is used to find the stability margin by assessing the real part of the critical eigenvalues which correspond to the oscillatory modes with the least damping ratio. Under uncertainty, a difficulty lies in the situation where the critical eigenvalues are related to different subsystems. In this system, the most unstable response can either be from the poorly damped filter in the active load 1 (AL1) or active load 2 (AL2). PF analysis is used then to associate each eigenvalue with their most dominant state variables. PFs show the influence of a state variable (i) on any given eigenvalue (j) and are defined in (21) [1]. This allows us to reliably identify and group the eigenvalues to unique subsystems: The two most critical eigenvalue pairs (in further text referred to as 1,2 and 3,4 ) in the system under study are related to the LC filter state variables of AL1 (V C,F1 , I L,F1 ) and of AL2 (V C,F2 , I L,F2 ), respectively. Table 4 shows the values of 1,2 and 3,4 for Cases A and B under nominal conditions. Cases A and B differ only in nominal values of input filter parameters of AL2, as defined in Table 3. However, it is observed that under different scenarios within the uncertainty ranges defined in Table 2, the real parts of both the critical eigenvalues, where each critical eigenvalue is related to either AL1 or AL2, do change. This indicates that nonlinear coupling effects between subsystems AL1 and AL2 are the result of the interactions among a set of parameters of AL1 and AL2. VBSA is used here to estimate the sensitivity (output variance) of the real part of critical eigenvalues (Re( 1,2 ), Re( 3,4 )) to uncertain system parameters (R line , C F1 , L F1 , R L,F1 , R 1 , C F2 , L F2 , R L,F2 , R 2 ) in the uncertainty range defined in Tables 3. The procedure is outlined below: (i). Sample the probability distribution of all uncertain parameters using the pseudo-random sampling methodology for VBSA as in [14] to generate an input vector of n unique scenarios. (ii). For each scenario, use the procedure described in Section 2. (a) Calculate a steady-state operating point for each scenario; (b) linearise the model to obtain the A matrix and calculate eigenvalues; (c) identify the two critical eigenvalue pairs, 1,2 and 3,4 , and associate them with load AL1 or AL2 using PF analysis. (iii). Calculate S i and S Ti for both Re( 1,2 ) and Re( 3,4 ) using the estimator presented in [14].
The total calculation time was ∼246 s to perform 110,000 unique model evaluations using an Intel Core i7-4790 @ 3.60 GHz with 16.0 GB of memory.

Results and discussion
VBSA results for the system under study are shown in Figs. 3 and 4. The first-order sensitivity indices S i indicate that in both Cases A and B, 1,2 responds to changes in AL1 parameters (R 1 , C F1 , L F1 , R L,F1 ), whereas Re( 3,4 ) responds mainly to changes in AL2 parameters (R 2 , C F2 , L F2 , R L,F2 ). This is expected as 1,2 and 3,4 have been classified by their dominant participation factors which associate them either with AL1 or with AL2. The results show that the most influential parameter for both Re( 1,2 ) and Re( 3,4 ) is R line , meaning that the reduction in the uncertainty of R line will reduce the output variance the most. Comparing Case A to Case B, it can be observed that they differ primarily in total effect indices S Ti and interaction effects (S Ti − S i ).
Under the matching filters case (Case A), the total variances of Re( 1,2 ) respond to changes in both AL1 and AL2 parameters, indicating interactions between the parameters of these two load subsystems, i.e. varying multiple parameters simultaneously produces a non-additive response on the movement of the eigenvalue. Of interest is the significant interaction effect of R line which suggests that the behaviour of R line on Re( 1,2 ) depends on the other interacting variables. This phenomenon is investigated further in the following section.
The results in Fig. 4 show that Case B exhibits minimal interaction effects and no interaction between AL1 and 3,4 , and AL2 and 1,2 ; meaning that parameters with zero interaction terms will have no effect on the sensitivity of any other parameter. For example, changes in the value of L F2 in the uncertainty range will not affect Re( 1,2 ) associated with AL1.
A limitation of using the total effect sensitivity measure is that it cannot identify which particular parameters are interacting. Modification of the adopted sampling procedure and estimator can be made to accommodate the 2nd order and higher order sensitivities [15], at the expense of increased total computation time. Additionally, it is important to note that VBSA indices only show the importance of particular variable(s) and do not quantify whether increasing or decreasing a variable will destabilise or stabilise the system.

Modal analysis of source-side resistance
To further investigate the effect of R line on the stability of Case A system, the local sensitivity ∂Re( 1,2 )/∂R line is calculated for different values of AL2 input filter resonant frequency. This is achieved by adjusting parameters L F2 and C F2 . Results presented in Figs. 5 and 6 indicate drastically different behaviours of the system under perturbations in R line . The results in Fig. 5 can be classified into three regions: Region 1) line resistance is stabilising due to equivalent series damping of the filters; Region 2) line resistance is destabilising due to increase in filter coupling effects; Region 3) critical eigenvalue insensitivity to changes in R line when the two filters are matched; For clarity, a cross-sectional view of Fig. 5 is shown in Fig. 6 when L F2 is nominal. From this, it can be seen that for R line to have a damping effect, the system has to operate in Region 1.

Validation via time-domain simulation
The detailed switching model of the system, as shown in Fig. 1, is simulated using Simulink/SimPowerSystems library. The effects of changes in R line under different operating conditions, as identified in Section 4, are verified using time-domain simulations. To investigate the effect on the stability, the time-domain simulation models are configured to operate in Region 1 (R line stabilising) and Region 2 (R line destabilising). Fig. 6 is used to select AL2 input filter parameters. In addition, the load resistances R 1 and R 2 are chosen so that the system operates near instability, as summarised in Table 5.
Referring to the simulation results in Figs. 7 and 8, at t = 0, the system operates in the steady state at the nominal line resistance (R line = 0.1 ). R line is then incrementally varied by ± 0.01 every 1 s to observe the effect on V bus . In Region 2, unstable oscillations occur when R line is increased to 0.14 at t = 4 (Fig. 7b), whereas in Region 1, R line can be safely increased to the upper bounded value of 0.2 without exhibiting instability (Fig. 8b). These results confirm the predicted system behaviour in Section 4. Decreasing R line results in the oscillatory behaviour in both Regions 1 and 2. In Region 2, the instability occurs when R line is decreased to 0.02 (Fig. 7a) compared to Region 1 instability at 0.06 (Fig. 8a). It can be concluded that R line has both damping and coupling effects on AL1 and AL2 and under certain operating condition, it can destabilise the system.

Conclusion
In this paper, a system with two active loads under parametric uncertainties is analysed. An improved analysis approach has been proposed using VBSA in order to identify system parameter(s) that have a significant impact on the small-signal stability. For the system under study, source-side resistance interacts with other uncertain parameters when input filter resonant frequencies of loads 1 and 2 are matched. Further analysis using local SA reveals that source-side resistance can exhibit different effects on the small signal stability: stabilising (damping), destabilising, and no impact (insensitivity). The predicted model behaviours are validated through time-domain simulation.