An Enhanced Two-Stage Grid-Connected Linear Parameter Varying Photovoltaic System Model for Frequency Support Strategy Evaluation

Securing a future sustainable decarbonised economy involves moving towards a system with rising penetration levels of distributed photovoltaics (PVs) within the low voltage distribution network (LVDN). This power system evolution is displacing conventional generators and has resulted in a decline in inertia that is essential for frequency stability. Emerging network codes require PV generators to maintain a scheduled curtailed active power (CAP) reserve for under-frequency contingencies. In this paper, the development, verification and application of an enhanced, two-stage grid-connected, state-space, linear parameter varying (LPV) PV system model is presented. The LPV model provides accurate and efficient modelling for PV systems over the wide range of operating points associated with curtailed active power and is suitable for power systems with large numbers of distributed PV systems capable of frequency support in the LVDN, to be simulated within reasonable simulation times. In addition, the LPV model can be used to investigate voltage rise due to reverse active power. The model performance is evaluated using recorded experimental data with step changes in irradiation and active power curtailment. The measured data is generated from a power hardware in the loop (PHIL) testbed. The model’s performance is investigated on an adapted radial European LVDN benchmark with several distributed PV systems to present some of the challenges, opportunities and benefits. Step changes in solar irradiation are used to evaluate the dynamic behaviour of the LPV model compared to a discrete-time electromagnetic transient (EMT) model. A frequency droop control characteristic for frequency support is demonstrated. The results show a computational burden reduction of 132:1 compared to the EMT model and demonstrate the voltage rise due to reverse active power from providing frequency support during under-frequency contingencies.


Introduction
Electricity generation using low carbon technologies such as solar PV and wind turbine systems is beginning to play an essential part in achieving a sustainable decarbonised economy. Over the past decade, PV module prices have dropped by 79% leading to an exponential growth from a globally installed capacity of 6 GW in 2006 to 638 GW in 2019 [1][2][3]. In this transition, the existing electricity network is evolving towards a decentralised infrastructure with the decline of large fossil-fuelled conventional generators and the proliferation of millions of small geographically distributed PV systems, presenting multiple challenges to system operators [4]. One of the main challenges relates to This paper presents the development of an LPV PV system model to accurately capture the steady-state and dynamic response of a two-stage PV system capable of operating away from maximum power point (MPP) for under-frequency contingencies. To the author's knowledge, this is the first presentation of an LPV model of a two-stage PV system. The development of an LPV model is significant for simulation-based investigations of large networks with multiple distributed PV systems that otherwise would be time-consuming and computationally expensive using detailed switching EMT models. The LPV model is an addition to previous work modelling a comprehensive linear state-space two-stage (LSS) PV system [27]. In that work, the PV system is linearised around a single equilibrium point at MPP, offering fast and accurate modelling. However, the LSS model is only valid for a very limited operating range with curtailed active power (CAP) near the equilibrium point of linearisation. In order to extend this operating range to anticipate the range required for frequency support services, multiple linear models are combined in the LPV model. The model is derived from an EMT model and both are validated using recorded experimental data for dynamic and steady-state performance at several solar irradiation intensities with an active power reserve. The recorded experimental data was created using a scaled PV system with a power hardware in the loop: real-time simulator (RTS).
The LPV PV system model is implemented on an altered radial European LVDN benchmark with low synchronous inertia to investigate the impact of declining system inertia and implement possible frequency support solutions [17,19]. The network is comprised of a steam turbine generator, eight network loads and seven 32 kWp distributed PV systems at each load centre to represent 100% penetration. Instead of using load perturbations [15,17,19,21], the effect of variable solar irradiation is used to assess the effectiveness of a droop controller response to offset adverse impacts at high penetration levels.
The development of the LPV PV system model, its validation with experimental data and its application in a distribution network are presented as follows. A description of the two-stage PV system under investigation is presented in Section 2. The development of the LPV model is described in Section 3. A comparison between the experimental measurements and the created LPV model are discussed in Section 4. In Section 5 the LV distribution network under investigation is described and the operation of the created LPV model is evaluated for both frequency response and voltage transients for various step changes in irradiation. The LPV model's dynamic and steady-state response is compared with the full EMT model from which it is derived. The conclusions are outlined in Section 6.

PV System
The two-stage PV system under investigation is presented in Figure 2a. The system comprises of a PV array, step-up DC/DC converter, DC/AC converter, an LCL filter configuration and controllers. The red outlined section represents the PV system components and controllers derived as an equivalent state-space model in [27]. The PV array is represented by the equivalent single diode circuit model equations [31]. The instantaneous operating voltage of the array is used to solve the array current using iterative methods such as the Newton-Raphson method [32,33]. Grid support services for frequency support are implemented using the control structure, as shown in Figure 2b. At any time, the value of curtailed active power (CAP) P cap for power reserve can be a constant predefined by the network operator or dynamically updated as required. Based on the measured frequency, the frequency support controller (FSC) determines P f s which for under frequency is the portion of P cap delivered to the grid and for over frequency is the portion of P cap curtailed from maximum power point grid delivery. This results in an instantaneous curtailed active power for the PV system of P ic = P cap − P f s . Based on this instantaneous curtailed active power, the active power reference for the PV system is P re f = P mpp − P ic . The active power controller as shown Figure 2b can operate at MPP (P mpp ) or operate below MPP with an instantaneous curtailed active power. As the PV array power is a non-linear function of the operational voltage v pv , it is difficult to predict the maximum power available at any instance, and when operating with a curtailed capacity, this can be achieved at voltages below and above the maximum power point voltage on the P-V curve. Also, this relationship continuously changes with variations in irradiance and panel temperature. Therefore, a PI controller regulating power is not effective as power does not increase linearly with duty cycle [14]. For this reason, the active power controller regulates the operating voltage v pv rather than power. Traditional maximum power point tracking methods such as perturb and observe (P&O) and incremental conductance (IC) do not require a measurement input of irradiation and temperature. To implement grid support services, as proposed in this work, it is essential to measure these meteorological parameters so alternative MPP techniques are considered. The approach adopted is to generate 2D look-up tables (LUTs) that use these measured inputs to determine the PV array maximum power P mpp . A 3D LUT provides a voltage reference v pv re f for the PV array [16]. The LUT method permits regulation of the PV array power to a given reference power P re f . The maximum power available P mpp is calculated using a 2D LUT dependent on irradiation (G) and temperature (T). A 3D LUT provides a voltage reference v pv re f for the PV array [16]. The reference voltage for the PV array v pv re f is a function of the power reference P re f , irradiation (G) and temperature (T). For operations other than the calculated points, the data is interpolated using trilinear interpolation. The LUT data is generated offline using the PV panel parameters from NREL System Advisor Model (SAM) [34] and solving the transcendental equation for a single PV panel over the expected operating space. The addition of panels in series or parallel changes the operating voltage v pv re f (v pv re f = v pv re f · N s ) and the power reference P re f P re f = P re f N s ·N p . In a physical implementation the 3D LUT and DC/DC controller will set PV operating voltage v pv re f . Continuous comparison of the modelled and measured PV array output power will allow the PV plant operator to identify faults and partial shading effects and notify the network operator. The model presented here does not implement this function and would also require further development to include partial shading effects.

PV System Modelling
In modelling the PV system there are several approaches that can be taken as shown in Figure 3. Each method has its own advantages and disadvantages related to simulation speed and accuracy. At the top level, the EMT model represents a component level system description that includes discrete switching devices that require the generation of PWM signals. In this work, these signals are generated at 10 kHz and require a 1 µs sampling time to provide a 1% duty cycle resolution. The system controllers use a 30 µs sampling time and data is recorded every 60 µs. The PHIL model uses a sampling time of 30 µs for both the controllers and the recording of data. In this case, the PWM signals are generated external to the simulation environment using a gate array. While the EMT model provides the most detailed information for comparison to the PHIL model, it imposes a large computational burden. In order to increase simulation speed and maintain a high level of accuracy, the switching device behaviour can be averaged over a switching period, creating an average model. This eliminates the need for a small sampling step time to generate the PWM signals and allows the use of variable step solvers that allow for increased sampling steps at steady-state operation. The computational burden can be further reduced by linearising the nonlinear average model. The linear model is time-invariant and is suitable for controller design, stability analysis and parametric sensitivity analysis and yields accurate results over a limited operating range near the equilibrium point [27]. This operating range is extended using the LPV model framework containing multiple linear models. In this work, the linear and LPV models are compared to the EMT model; the sampling time and data recording is reduced to 300 µs with a linear model sampling time of 3 ms which provides 27 samples in a step response of 80 ms for the PV system. The solvers are applied with both a maximum step size and relative tolerance of 0.003. The effect of different data logging sampling times used in the EMT model has an insignificant impact on simulation time.

PV System EMT Model
The developed EMT model of the three-phase PV system is shown in Figure 4. The model is used as the basis for derivation of the linear state-space PV system model's differential equations and provides a reference for comparison. The parasitic losses related to non-ideal passive components are incorporated as resistor elements, derived from measurements and measured experimental data. The switching devices (IGBTs and diodes) are shown in red to display losses and are modelled as ideal switches with characterised losses depicted as a comparable voltage source [27].

LPV Model Development
In this section, an LPV model for the two-stage PV system is developed. The concept of LPV modelling is based on the traditional gain-scheduling framework [35]. In gain-scheduling, the non-linear system model is linearised at several operating points of interest resulting in a collection of local linear time-invariant (LTI) models. For each linear model, an LTI controller is designed. These controllers are then interpolated to create a global control solution for the non-linear system model [30]. The interpolation function is called the scheduling function and depends on the current operating point of the system. The changes in the operating point are described by an introduced signal called the scheduling signal (ρ). The parameters of the controller now depend on the varying scheduling signal. The controllers are now parameter varying.
Using the gain scheduling approach, an LPV model is created using the framework shown in Figure 5 by using a finite gridded domain of linear state-space models from several equilibrium points. Initially the LPV model is developed using the equivalent averaged time-invariant PV system differential Equations (1) and (2) in terms of z states x(t) ∈ R z , k inputs u(t) ∈ R k , q output measurements y(t) ∈ R q related by n first-order state derivativesẋ for continuous-time systems, t ∈ R where f pv is the vector of PV system functions containing n differential equations, y is the system output measurements and h is the vector of output functions.
The linear time-invariant (LTI) state-space representation of (1) and (2) in continuous time is represented by (3) and (4) for t ≥ t 0 and the initial conditions, x(t 0 ), where A is the state matrix A ∈ R z×z , B is the input matrix B ∈ R z×k , C is the output matrix C ∈ R q×z and D is the feed-forward matrix D ∈ R q×k .
The main advantage of using a linear parameter varying system S(ρ) instead of an LTI system model is the ability to account for dynamic changes using a time-varying scheduling signal ρ(t), this allows for system analysis and controller design subject to time-varying trajectories. The scheduling signal is a function ρ = f (y, r) of the independent variables y and external signals r, this time variation is not known in advance but it is a measurable quantity at each instance in time. These independent variables are referred to as scheduling parameters. The LPV system is described as an array of parameter-dependent linear systems, whose state space representation depends on the scheduling parameters. The scheduling parameters are described as freely varying parameter values in the region P, ρ ∈ P. The scheduling parameters chosen for the PV system are irradiation and instantaneous curtailed active power (G, P ic ), where P ic = P cap − P f s . The LPV equivalent of a nonlinear system is presented in (5) and (6).

Linearisation
The linear equations are obtained from the non-linear differential equations by using the first order terms of Taylor's expansion series (7) [29] from an equilibrium point (x e , u e ) by solving each derivative of the vector f pv (x e , u e ) = 0 where x e is the vector of equilibrium points for each state variable, u e is the vector of equilibrium points for each input, (δu = u − u e , δx = x − x e ) denotes the variation from the equilibrium points and ∇ f pv is the gradient of f pv evaluated at (x e , u e ).
Linearising the PV system about a single equilibrium point yields the linear state-space model (8) and (9), these linear models are used for the LPV framework as shown in Figure 6b: where The LPV model is created by using a finite gridded domain of linear state-space models from several equilibrium points (x i e , u i e ). The collection of LTI models describes the behaviour of the nonlinear system throughout the scheduling space P. Each model is a local approximation that describes the nonlinear system for small deviations. The entire set of LTI models (x i e , u i e ) ∈ P yields a linear parameter-dependent linearisation (ρ e ). A regular grid is chosen to schedule between local LTI models, as shown in Figure 6b. The scheduling parameters' irradiation and instantaneous curtailed active power (G, P ic ) define the local model.

Interpolation of Grid Points
For scheduling parameters of irradiation G and instantaneous curtailed active power P ic that are located on-grid locations, as shown in Figure 6a, the LTI model circled in blue best describes the operating dynamics at that instance. Using these scheduling parameters, a linear parameter-dependent system Sρ e is capable of locally describing the nonlinear model (10) and (11).
When the scheduling parameter values ρ e are located off their grid locations as shown in Figure 6b, the scheduling function uses inverse distance weighted (IDW) interpolation of the 2n ρ nearest neighbouring LTI models in the scheduling space, where n ρ is the number of scheduling parameters. The nearest LTI models are obtained using the k-nearest neighbours (k-NN) algorithm by solving the Euclidean distance from the scheduled model ρ e as shown in Figure 6b. The IDW interpolation is based on the constant nodal function (12) where matrices A i , B i , C i , D i are interpolated as weighted sums, n m is the number of nearest models (n m = 2n ρ ) and w i is the weight assigned to each n m LTI model. The inverse distance weighted method is influenced most by the nearest model and less by the more distant models. The total sum of the weights is equal to one (13). Weights are calculated using (14) where h i is the distance from the scheduled model to the weighting LTI model and h j is the distance from the scheduled model to the nearest n m LTI models.
The LPV model is a family of LTI models linearised at various inputs of irradiation and reserve power (G, P ic ) so that each has their unique equilibrium points for states x, inputs u and outputs y that depend on these linearisation inputs. Taking this into account, the LPV model is extended to include the offsets in the state, input and output (ẋ, x, u and y) variables. Using the equilibrium offsets, the LPV model in (10) and (11) is represented by (15) and (16), which is an affine form of the LPV system: where δẋ(ρ e ), δx(ρ e ), δu(ρ e ) and δy(ρ e ) are offsets from the equilibrium points of irradiation and the instantaneous curtailed active power (G, P ic ) forẋ, x, u and y at the scheduling function ρ e = f (G, P ic ).

LPV PV System Model
The LPV PV system model is a linear system array approximation of the time average system dynamics. The PV system consists of a 1.25 kWp array made up of five series Conergy PH 250 Wp multi-crystalline panels. All PV system parameters are given in [27]. The PV system consists of 25 states and is considered to be of low complexity and can be efficiently solved without the need for truncation (the number of scheduling variables is at most 2 or 3 for this size system) [36].
The percentage error in output active power for a single LTI model linearised at 550 W/m 2 operating over the expected operating range with active power curtailment is shown in Figure 7 with a zoomed section over the operating range of 400-1000 W/m 2 , this emphases how a single models accuracy varies with curtailed active power. The percentage error is calculated using (17) where P EMT is the output active power of a full EMT model, and P LSS is output active power of a single linear state space model. The single model remains accurate when operating off MPP within -3.2% to +0.5% in the range of 800-1000 W/m 2 with a curtailed capacity of 500 W (40%). While operating at 600 W/m 2 for the same curtailed capacity, an error of 16% exists and further increases as irradiation lowers. The operating range of the model can be extended by using multiple linear models linearised at various operating conditions. As an example, Figure 8 shows the (P-V) characteristic curve evaluated at various irradiation levels at a constant temperature of 25 • C, MPP operation is highlighted on the knee of each curve by the dashed magenta line. The nonlinear curve for MPP voltage can be captured by using several linearised models. Three models linearised at [100, 550, 1000] W/m 2 using (7) are shown in Figure 9.
Increasing the number of linear models reduces the error highlighted in Figure 9 and provides a more accurate representation of output active power. Figure 10 shows how the percentage output error decreases for MPP operation using a total of up to 20 linear models, linearised in the range of 100 to 1000 W/m 2 with even spacing. A zoomed section over the operating range of 400-1000 W/m 2 shows that maximum error is 3% for a single model and the error rapidly decreases with just two models to 0.75%. In total, 10 linear models are chosen for changes in irradiance as the benefits above 10 models are minimal.   In order to capture the full operating envelope for this particular application, additional models are required for the inclusion of active power curtailment. Figure 11 shows how the PV array voltage varies for active power curtailment at 550 W/m 2 with a constant temperature of 25 • C. Depending on the irradiation level and active power curtailment, the operating range of the PV array varies widely, making it difficult for a single scheduling parameter to capture the expected operational range as shown in Figure 12.  Increasing the number of linear models for active power curtailment is shown in Figure 13 for a constant irradiance of 550 W/m 2 . This shows how the percentage output active power error decreases with an increasing number of linear models. The percentage error of the LPV model is calculated using (18) where P LPV is output active power of the LPV model. For active power curtailment, the maximum error in output active power decreases rapidly with an increasing number of models from 54.8% for a single model to 5.6% for two models and 0.2% for six models with an active power curtailment in the range of 0-500 W. The zoomed section in Figure 13 shows how rapid the error decreases when using three or more models with a curtailed capacity in the range of 0-500 W. From these results, six linear models are chosen to capture the PV systems' operating points when operating off MPP in the range of 0-500 W.
A combination of linear models for irradiation and active power curtailment is shown in Figure 14 for constant irradiation of 550 W/m 2 with a fixed number of six linear models for active power curtailment. For irradiation, various numbers of linear models are chosen (1, 2, 5, 10, 60), this gives a maximum combination of 120 linear models. The difference in simulation time from a single model to 120 models is small, and the most significant constraint in expanding the number of models applied is model generation time. From these results, considering model generation and accuracy, 60 linear state-space models were chosen for use with a regular grid, as shown in Figure 6. The scheduling parameter values irradiation G and instantaneous curtailed active power P ic consist of 10 and six values, respectively. The gridded array contains a linear state-space model at every combination of G and P ic .

Model Validation
An experimental hardware testbed (EHT) as shown in Figure 15 was built to assess the EMT and LPV models. The measured experimental data is compared against the LPV model for dynamic and steady-state response. The testbed was built to emulate a scaled PV system with a power rating of 1.25 kWp. All system parameters remain the same for both the EMT and EHT models as described in [27].

Results
The results in this section are used to assess the accuracy of the developed LPV model against an experimental hardware testbed, a single linear state-space (LSS) model linearised at 500 W/m 2 with 0 W curtailed and a full EMT switching PV system model simulated in MATLAB. All models receive step changes in irradiation, as shown in Figure 16 with a constant temperature profile. An active power curtailment of 300 W with unity power factor (Q re f = 0 VAr) is applied at 3.5 s to demonstrate the capability of active power curtailment for frequency support. The power curtailment is achieved by making the PV array operating voltage exceed the MPP voltage, as shown in Figure 17 where v pv re f is the reference voltage to achieve active power curtailment and v pv re f mpp is the reference voltage for maximum power point operation. Operating at the right-hand side of the P-V curve in this way is preferred for improved dynamic performance, increased converter efficiency and controller robustness [14].  The active power extracted from the PV array is shown in Figure 18. All models generate similar steady-state operating values for PV array power as the PV array uses the equivalent single diode circuit model equations [31] and solves for array current based on the instantaneous operating voltage [32,33]. The active power output for all models is shown in Figure 19. Once operating off MPP with a reserve of 300 W the accuracy of the LSS model starts to decline as the irradiation decreases. The LSS model is unable to capture the full operation of the PV system with an active power reserve. The LPV model produces similar dynamics and steady-state operation for various irradiances with active power curtailment when compared to both the measured experimental data and EMT model. This validates the capability of the LPV model to capture the entire operation of a PV system for active power curtailment.  When operating at MPP, the duty cycle of the step-up DC/DC varies non-linearly by 2.15% as shown in Figure 20. As with Figure 9, the linearised model accuracy can be improved by including several linear models with three models linearised at (100, 550, 1000) W/m 2 as shown in Figure 20. Operating with a curtailed active power as required for GSS, the duty cycle has a broader range of nonlinear variation of 7.8%, as shown in Figure 21. The data for CAP 0 W is the same as Figure 20 while the effect of increasing levels of active power curtailment results in the high variation in duty cycle. The DC/AC converter for the same operation presents a linear relationship between the duty cycle and the output active power for control of both active and reactive power as shown in Figure 22 and Figure 23 respectively. These results show the non-linearity of the DC/DC step-up converter duty cycle variation over the full irradiation range and indicate that this may result in significant steady-state errors in a linearised model.  Step-up converter duty cycle with curtailed active power.

Network Simulations
The purpose of the developed LPV PV system model was to enable system-level simulations and its usefulness was tested on a radial LVDN. Two different scenarios are investigated: 1. Maximum power point tracker (MPPT) operation with no frequency support provided. 2. A P(f) droop characteristic for frequency support.
The results are validated against a full EMT model of the LVDN. The network disturbance considered focuses on the impact of variations in solar irradiation and all loads remain constant.
Step changes in solar irradiation are used for dynamic performance characterisation. Two irradiation profiles are used, as shown in Figures 24 and 25 for over and under frequency events. The temperature profile of the PV array is constant for all analyses (25 • C).

Three-Phase Network Description
The impact of increasing PV system penetration is analysed on a modified radial, 20 nodes symmetrical three-phase commercial CIGRE European LVDN benchmark [37]. The network consists of seven distributed load centres totalling 217.8 kW. The LVDN is modified by including a steam turbine generator and by adding seven PV systems (7 × 32 kW) at each distributed load centre to represent 100% penetration. The modified one-line diagram of the European LVDN is shown in Figure 26. The detailed parameters of the distribution line are given in Table 1 [37]. As the lengths of the distribution lines are short, the capacitance is small and considered negligible.    The LVDN is modified by replacing the slack bus with a steam turbine generator, as shown in Figure 26. This provides a controllable system with a realistic performance for low inertia networks, allowing for the evaluation of high penetrations of PV with and without frequency support strategies on frequency dynamics. Since frequency dynamics are the primary objective, a single synchronous generator provides the full operating range of loads.
The electrical generator in this study is described by a sixth-order dq state-space model. The model takes into account the field circuit and a single equivalent damper winding [38]. A steam turbine is used to represent the prime mover without a reheat steam turbine as shown in Figure 27 and governor based on TGOV1 [39,40]. The model is characterised by the steam turbine, the main steam control valve and the control valve limits (V max and V min ). The steam turbine and the main steam control valve are described by first-order transfer functions with time constants T tur and T va respectively. The governor regulates the turbine's mechanical power (P mac ) into the synchronous generator based on the combined angular velocity (ω m ) and the generator power set-point (P gen setpoint ) with droop characteristic (R droop ). The droop characteristic provides primary frequency support speed regulation allowing for even generator load sharing based on each generators reference. A 4% droop characteristic is chosen as the default [41]. Restoration of the system frequency to the nominal value of 50 Hz after an event causing an increase or decrease in load requires a supplementary control action in the form of a secondary controller that adjusts the generator power set-point P gen setpoint from the scheduled power P sch . The automatic generator controller (AGC) is initiated 30 s after a frequency excursion outside the controller limits as shown in Figure 27. As the secondary frequency controller is not the focus of this work, an additional integral controller is added to speed up system restoration when the frequency is below 50.1 Hz or above 49.9 Hz after a frequency excursion. This addition controller reduces the overall time window for which the simulation must be run to demonstrate the various responses. Typically, primary frequency control responds within 30 s and secondary frequency control responds in a range of 5 s up to 10 minutes after the change in frequency [41].  The generator's terminal voltage is regulated using a static excitation system based on an IEEE Type ST1C voltage source [42,43] rectifier exciter as shown in Figure 28 with limited exciter voltage (V amax and V amin ). The voltage regulator gain and the excitation system time constant are expressed by K a and T a , respectively. The transient gain of the system is reduced by implementing a feedback gain expressed by the rate feedback gain (K f ) and the rate feedback time constant (T f ). Parameter values are given in p.u. in Table 2.

Response to Increasing Solar Irradiation
The LPV model's accuracy for over frequency events is evaluated in this section. The step changes in solar irradiation that increase generation by 10, 20 and 30% of rated capacity at constant temperature are shown in Figure 24. The droop characteristic for over frequency events is in the range of 50.2 to 51 Hz (1.6% droop). If the frequency is higher than 51 Hz the PV system fully curtails active power generation.
The system frequency response for the PV with an MPPT controller only and no frequency support is shown in Figure 29. When the PV systems are operating with an MPPT controller, the objective is to deliver as much active power as possible. Since each PV system has a fast dynamic response, the generator is unable to respond quickly enough to balance the system active power demand causing the frequency to rise. As a result, a maximum frequency excursion of 53.9 Hz and a RoCoF of 2.4 Hz/s occurs (exceeding rated limits) for a sudden increase of generation equivalent to 30% of capacity. The system frequency response for the droop controller is shown in Figure 30. The implementation of droop control alone from each PV system provides sufficient frequency support by maintaining frequency within limits, reducing the maximum frequency deviation by 79% (53.9 to 50.8 Hz). Using the droop controller does not reduce the initial RoCoF. The voltage at bus 17 is measured to assess model accuracy and voltage rise due to increasing solar irradiation. Bus 17 is chosen due to the distance from the generator and low local load to demonstrate worst-case voltage rise due to reverse active power flow. The voltage at bus 1 is regulated to 400 V. The voltage at bus 17 is shown in Figures 31 and 32 for MPPT and droop control respectively. As solar irradiation increases, PV active power generation becomes larger than local load demand resulting in reverse-active power flow, causing voltage rise. The sudden increase in voltage rise in Figures 31 and 32 is due to the step increase in solar irradiation and the fast response of the PV system. The transient response in Figure 32 is due to the droop controller working from a feedback loop using frequency. The event occurs before the controller responds, and active power cannot be curtailed until frequency increases. It can be seen in Figure 32 that the voltage settles briefly before active power is curtailed. Perhaps a Q(f) or P(Q) characteristic dependent on RoCoF could limit these voltage transients [44].

Response to Decreasing Solar Irradiation with Active Power Reserves
The LPV model's accuracy for under frequency events is evaluated in this section. The output is compared against an EMT model for step changes in irradiance. The step changes in solar irradiation decrease generation by 10, 20 and 30% for a constant temperature using the irradiation profile as shown in Figure 25. The droop characteristic of under frequency events is in the range of 49.8 Hz to 49 Hz (1.6% droop). If the frequency is less than or equal to 49 Hz, the PV system operates at maximum power. Once system frequency is above 49.8 Hz, all PV systems operate with the specified curtailed active power reserve of 16% for the droop controller.
The system frequency response for the PV with an MPPT controller only and no frequency support is shown in Figure 33. When each PV system operates with an MPPT, the maximum frequency deviation is 4.5 Hz, and a RoCoF of 2.77 Hz/s occurs for a sudden loss of PV generation equivalent to 30%. The system frequency response for the droop controller is shown in Figure 34. For a sudden loss of PV generation equivalent to an imbalance of 30% with each PV system operating with a droop controller, the maximum frequency deviation is 2.24 Hz which is a 50.22% reduction compared to MPPT operation. The droop controller, however, does not reduce the system RoCoF. The loss of PV generation due to decreasing irradiation causes bus voltage to fall. The voltage at bus 17 is shown in Figures 35 and 36. Without frequency support from each PV system, the voltage at bus 17 falls nearly instantly for each step decrease in solar irradiation. The implementation of the droop controller dampens the fall in voltage by delivering active power locally from an active power reserve. Active power plays an essential role in voltage regulation, where the network line's R/X ratios are high; this is especially true for LVDNs. The voltage profile at each bus can also be further improved by implementing a Q(V) droop controller offsetting the reactive component [45].  In summary, the results presented in this section show the validation of the LPV PV model for frequency support investigation and the potential of frequency support control with low synchronous inertia systems. The model was validated by comparing the dynamic performance with a full EMT switching model in MATLAB, as shown in Figures 29-36. From these time-domain simulations, the dynamic behaviour was assessed using system frequency and bus voltage. The LPV and EMT PV system models produce similar system responses and steady-state values for both increasing and decreasing step changes in irradiation with frequency support controllers. A summary of the simulation times and computational burden reduction between the EMT, average and LPV model is shown in Figure 37. The simulation time on average was reduced by 99.2%, from 37 h to 0.28 h; this is a computational burden reduction of 132:1 between the EMT and LPV model. These results demonstrate the benefits of the LPV PV model when investigating a high penetration of distributed PV systems in a large scale distribution network. Another advantage of this model is the ability to investigate bus voltage, which could violate voltage limits due to reverse active power when providing active power to support frequency.

Conclusions
In this paper, a linear parameter varying (LPV) model of a two-stage grid-connected PV system was introduced with flexible active power control for frequency support. The model builds on previously published work [27] to formulate a fast and accurate model, suitable for power system studies with increasing penetration levels of PV. The function of the LPV model is to accurately capture the non-linear behaviour of a two-stage PV system over the wide range of operating points associated with curtailed active power for GSS. The LPV model is validated with curtailed active power using recorded data from an experimental hardware testbed (EHT). The dynamic and steady-state performance at multiple emulated step changes in irradiation is assessed. The potential of the LPV model and the effectiveness frequency support using a P(f) droop control for PV systems were simulated on a 20 node isolated modified European LVDN with low synchronous inertia and compared against a full EMT switching model. These simulations demonstrate that the LPV model can emulate the behaviour of an EMT model with comparable dynamics and steady-state running values with a reduction in simulation time by 132:1.