Analysis of Control Strategies Based on Virtual Inertia for the Improvement of Frequency Stability in an Islanded Grid with Wind Generators and Battery Energy Storage Systems

: Rising levels of non-synchronous generation in power systems are leading to increasing difﬁculties in primary frequency control. In response, there has been much research effort aimed at providing individual electronic interfaced generators with different frequency response capabilities. There is now a growing research interest in analyzing the interactions among different power system elements that include these features. This paper explores how the implementation of control strategies based on the concept of virtual inertia can help to improve frequency stability. More speciﬁcally, the work is focused on islanded systems with high share of wind generation interacting with battery energy storage systems. The paper presents a methodology for modeling a power system with virtual primary frequency control, as an aid to power system planning and operation. The methodology and its implementation are illustrated with a real case study.


Introduction
In the last two decades, electric power systems have experienced very significant changes in their operational requirements, mainly due to a growing share of non-conventional renewable energy in the electricity generation mix, along with the operation of grids with a higher degree of complexity [1][2][3]. From the beginning of the massive electrical generation industry at the beginning of the 20th century, conventional synchronous generators have been responsible for guaranteeing frequency stability [4], among other tasks within the operation of a power system. However, mainly in weak or islanded networks, synchronous machines do not always have the capacity to respond fast enough to keep the system stable in the face of power disturbances caused by imbalances between generated and demanded power [5][6][7].
If an increasing share of electricity generation based on non-conventional renewable energy resources is added to this type of power systems, one of the consequences is a significant reduction of the equivalent inertia of the system [8]. This is due to the fact that such renewable generators have an interface for connection to the grid based on the use of electronic power converters, in which a frequency variation in the grid side does not entail a modification of the active power supplied, as the latter keeps a value close to the reference set by its corresponding control system. From the point of view of the power system, this is seen as a lack of inertial response. To address this issue, previous works have outlined the design and implementation of additional control schemes that provide certain complementary regulation capabilities to certain renewable generators (e.g., photovoltaic [9,10], wind [11,12], hybrid of both [3,13]) and battery energy storage systems [14,15], with the aim of modifying the active power injected into the network in response to variations in frequency, similar to how conventional synchronous generators do; hence, the use of the term "virtual inertia". Further, very interesting review papers addressing inertia and frequency control strategies can be found in [8,16,17].
One of the most used resources to improve frequency stability in island-type microgrids is a battery energy storage system (BESS), with an increasing degree of utilization in electrical systems worldwide, for example, in South Australia or China, with 100 MW/ 129 MWh and 200 MW/800 MWh, respectively [18]. However, despite the effectivity of BESS in technical terms, it is not possible to ignore its significant economic cost [19,20]. This has motivated the exploration of less cost-intensive alternatives or complementary technologies.
In this sense, scientific contributions presented in [21,22] propose the use of additional control schemes in variable speed wind turbines that make use of the inertial resources existing in their rotating masses, giving support to the frequency recovery of the system in the event of power disturbances. This ancillary service is also provided by wind turbines by controlling the pitch angle of their blades to limit the incoming wind power and to contribute to the power balance of the grid [23,24]. The application of these solutions allows a grid-friendly integration of wind turbines in terms of the frequency recovery capability of the whole system. In these works, detailed models are employed to accurately represent the time-domain behavior of the power system elements in a simulation environment. However, their application in multi-machine systems may be impractical due to the huge computational effort required. In inertial response and frequency control studies, the dynamics of different agents of the system can be represented by a load frequency control scheme [10,25] which allows to achieve a noticeable modeling simplification without penalizing the accuracy of numerical results. This paper proposes a methodological guide for modeling the elements of a wind-diesel-energy storage islanded grid, based on the load frequency control scheme along with its implementation in a computational environment. The resulting model can be used to represent other power systems with similar characteristics, and it facilitates the integration of additional energy agents: generation, storage, and consumption.
This work presents a comparative study of the potential benefits obtained by the implementation of two control strategies designed to improve the frequency control actions of the system: one, based on the use of BESS, and the other, on the use of the inertial resources existing in a variable speed wind turbine. The hypothesis of this research is that the combined action of these two solutions will lead to more effective results in the mitigation of frequency fluctuations and the reduction of the frequency nadir, compared to applying such strategies separately.
It is intended that, with the results of this study, the operators of weak and isolated power systems with a significant penetration of wind generation can have an analysis tool that allows them to define new operating procedures and expansion plans of the system in order to improve the quality of the frequency response to disturbances caused by power imbalances. This paper is organized as follows: Section 2 details the methodology used to carry out the analysis of the grid frequency control, describes the electrical system considered as a case study, and presents the computer simulation modeling of the system under study; Section 3 shows the results obtained when analyzing the behavior of the system in four operating scenarios along with a critical discussion; and Section 4 draws the main conclusions of the study.

System Description
As a case study, the electrical system of San Cristobal Island (Galapagos, Ecuador) has been considered, whose electricity demand is supplied by a hybrid wind-diesel generation system [26]. Given its geographical conditions, an interconnection with a continental power Energies 2021, 14, 698 3 of 18 system is technically and economically unfeasible, so this constitutes a typical case of a microgrid operating in island mode, with diesel generation plants and high penetration of wind generation. The simplified single-line diagram of this system is illustrated in Figure 1.
This operating situation constitutes a latent risk for the continuity of the power supply, since such frequency excursions could lead to the tripping of certain over/under-frequency protection relays installed in different pieces of equipment of the system, a situation that has actually occurred on occasions in the system under study. An immediate solution to this problem would be to increase the number of diesel units in the microgrid, in order to strengthen the equivalent inertial response of the system. However, this measure is difficult to apply given the island's enormous ecological interest [27,28].
The analysis of this system has been carried out with the Matlab/Simulink ® (Math-Works, Natick, Massachusetts, USA) computer simulation tool, in which the elements of interest have been modeled and parameterized taking as a reference the work presented in [26]. An imbalance between generated and demanded power of ±0.1 p.u. has been considered as a disturbing event to quantitatively evaluate the impact of incorporating each of the proposed strategies on frequency stability in the system under study. WT   Regarding the operating conditions of the case study, historical records of operation and dispatch of the real system have been analyzed to select a particularly problematic case in terms of frequency control, with a wind penetration of 66% in a peak demand situation (2424 kW) [26]. Under these high wind conditions, the system tends to undergo pronounced frequency fluctuations, given the inability of the existing thermal generation to fully absorb the fluctuations of injected wind power.
This operating situation constitutes a latent risk for the continuity of the power supply, since such frequency excursions could lead to the tripping of certain over/under-frequency protection relays installed in different pieces of equipment of the system, a situation that has actually occurred on occasions in the system under study. An immediate solution to this problem would be to increase the number of diesel units in the microgrid, in order to strengthen the equivalent inertial response of the system. However, this measure is difficult to apply given the island's enormous ecological interest [27,28].
The analysis of this system has been carried out with the Matlab/Simulink ® (Math-Works, Natick, MA, USA) computer simulation tool, in which the elements of interest have been modeled and parameterized taking as a reference the work presented in [26].
An imbalance between generated and demanded power of ±0.1 p.u. has been considered as a disturbing event to quantitatively evaluate the impact of incorporating each of the proposed strategies on frequency stability in the system under study.
In addition, four operational scenarios have been analyzed: in the first one, the grid frequency control tasks are performed only by diesel generators; in the second scenario, a virtual inertia algorithm is incorporated in the variable speed wind turbines to contribute to the frequency support of the system together with the operational synchronous generation; in the third one, the corrective action provided by a battery storage system is evaluated by implementing an emulated inertial response in the power controller of its electronic converter; and, in the fourth scenario, the impact of both improvement strategies acting simultaneously, and in interaction with conventional thermal generation, is analyzed.
For assessing the magnitude of frequency deviations experienced by the system in the different scenarios, an actual national grid code is taken as a reference framework. In this sense, the Ecuadorian Electricity Regulation and Control Agency establishes, in Dispatch and Operation Procedure V2, Chapter 5.3 [29], that the target frequency of the National Interconnected System is 60.00 Hz, with a range of variation in normal operation conditions between 59.85 and 60.15 Hz, except for emergency states, faults, and restoration periods.
It is worth noting that, in the case of island systems, grid codes tend to be more permissive in their minimum and maximum limits in respect to what is established for continental systems, given their constructive and operational characteristics. Table 1 summarizes the frequency ranges within which the calibration of the frequency relays protecting the Ecuadorian power system generators is recommended.

System Modeling
In traditional electrical systems, the inertial power response is directly associated with the kinetic energy stored in rotating machines. In a synchronous machine, this energy can be expressed as [4] E kinetic = 0.5J [kg·m 2 ] ω 2 where J is the moment of inertia of the turbine-generator set and ω is the angular speed of the rotor of the machine. On the other hand, the rotor speed variation depends on the balance between electrical and mechanical torque, as follows: where T m and T e represent the electrical and mechanical torque and P m and P e are the mechanical and electrical power, respectively. The rate of change of frequency (RoCoF) can be expressed as: where H is the inertia constant of the set, defined by the relationship between the kinetic energy stored in the rotating masses and the nominal power of the machine S: Since the frequency stability of the system will be evaluated, it is convenient to introduce the simplified models of each one of the generating units in the representation given by the oscillation equation. As a function of incremental magnitudes, in values per unit, according to [4,30], this equation is: where H T is the equivalent inertia constant of the system and D is the coefficient of variation of demand with frequency. From the previous equation, it is possible to obtain the model of the system as shown schematically in Figure 2 [4,30], where: • ∆P L is the deviation in the value of the load that the system experiences in per unit. It will be negative if the charge decreases and positive when the opposite happens. • ∆P dn represents the power compensation contribution generated by each diesel generator to contribute to the recovery of the system frequency. • v n represents the incident wind speed at each wind turbine. • w i is the participation factor of each generation unit, both in the frequency support of the system and in the coverage of demand. In per unit values, this factor is defined as the relationship between the installed capacity of each generation unit and the base power considered for the analysis. according to [4] and [30], this equation is: where HT is the equivalent inertia constant of the system and D is the coefficient of variation of demand with frequency. From the previous equation, it is possible to obtain the model of the system as shown schematically in Figure 2 [4,30], where:  ΔPL is the deviation in the value of the load that the system experiences in per unit. It will be negative if the charge decreases and positive when the opposite happens.  ΔPdn represents the power compensation contribution generated by each diesel generator to contribute to the recovery of the system frequency.  vn represents the incident wind speed at each wind turbine.  wi is the participation factor of each generation unit, both in the frequency support of the system and in the coverage of demand. In per unit values, this factor is defined as the relationship between the installed capacity of each generation unit and the base power considered for the analysis. The modeling of the diesel and wind generation units is described below, as a summary of what is described in greater detail in [26], as well as the energy storage unit.   The modeling of the diesel and wind generation units is described below, as a summary of what is described in greater detail in [26], as well as the energy storage unit.

Diesel Power Plant
The thermal power plant of the system under study is made up of 3 units, each of which has a synchronous generator with salient poles and an exciter with rotating diodes. Its general configuration is schematically represented in Figure 3 [26] and its main characteristics are listed in Appendix A. The thermal power plant of the system under study is made up of 3 units, each of which has a synchronous generator with salient poles and an exciter with rotating diodes. Its general configuration is schematically represented in Figure 3 [26] and its main characteristics are listed in Appendix A. As a diesel engine model, the one proposed by [31], particularized for this case study, has been used. Figure 4 shows the schematic representation of the proposed model and the values assigned to its different parameters are indicated in Appendix B, where  ωr is the angular velocity of the generator.  Tmax and Tmin are the mechanical torque limits, Tm, that are input to the actuator.  ΔPm represents the deviation of the mechanical power compared to a pre-established initial value, Pm0.
It is important to note that the physical processes involved in the control actions are represented by first-order systems with time constants T2 and T3 and by a delay function of value T1. Additionally, KP and KI are the PI controller gains.

Wind Power Plant
The San Cristóbal wind farm is made up of three variable speed wind turbines with full power electronic converter (type IV, from the manufacturer MADE, model A-59). Each one has a brushless synchronous generator driven by a horizontal axis turbine. It is schematically represented in Figure 5 [26].  As a diesel engine model, the one proposed by [31], particularized for this case study, has been used. Figure 4 shows the schematic representation of the proposed model and the values assigned to its different parameters are indicated in Appendix B, where As a diesel engine model, the one proposed by [31], particularized for this case study, has been used. Figure 4 shows the schematic representation of the proposed model and the values assigned to its different parameters are indicated in Appendix B, where  ωr is the angular velocity of the generator.  Tmax and Tmin are the mechanical torque limits, Tm, that are input to the actuator.  ΔPm represents the deviation of the mechanical power compared to a pre-established initial value, Pm0.
It is important to note that the physical processes involved in the control actions are represented by first-order systems with time constants T2 and T3 and by a delay function of value T1. Additionally, KP and KI are the PI controller gains.

Wind Power Plant
The San Cristóbal wind farm is made up of three variable speed wind turbines with full power electronic converter (type IV, from the manufacturer MADE, model A-59). Each one has a brushless synchronous generator driven by a horizontal axis turbine. It is schematically represented in Figure 5 [26]. P  • ω r is the angular velocity of the generator. • T max and T min are the mechanical torque limits, T m , that are input to the actuator. • ∆P m represents the deviation of the mechanical power compared to a pre-established initial value, P m0 .
It is important to note that the physical processes involved in the control actions are represented by first-order systems with time constants T 2 and T 3 and by a delay function of value T 1 . Additionally, K P and K I are the PI controller gains.

Wind Power Plant
The San Cristóbal wind farm is made up of three variable speed wind turbines with full power electronic converter (type IV, from the manufacturer MADE, model A-59). Each one has a brushless synchronous generator driven by a horizontal axis turbine. It is schematically represented in Figure 5 [26].
Its representation derives from some of the author's proposal [21], consisting of a simplified electromechanical model of a variable speed wind turbine conceived particularly for the study that is the object of this work. The relevant technical information for the implementation of the wind turbine model is summarized in Appendix C. Given that the objective of this study is based on the application of virtual inertia, a brief presentation of the concept and its application in a variable speed wind turbine will be given next.
This concept consists of introducing changes in the active power command signal dispatched by the wind turbine in response to the variation of some variable of interest; in this case, the network frequency. This characteristic, supported by the kinetic energy stored in the wind turbine, is implemented in the active power controller of the machine. The virtual inertia control scheme used in this study has been taken from [21] and is outlined later in Figure 8. Before describing the mechanics of the method, some aspects concerning the control criteria of the active power dispatched by the wind turbine are addressed. Its representation derives from some of the author's proposal [21], consisting of a simplified electromechanical model of a variable speed wind turbine conceived particularly for the study that is the object of this work. The relevant technical information for the implementation of the wind turbine model is summarized in Appendix C. Given that the objective of this study is based on the application of virtual inertia, a brief presentation of the concept and its application in a variable speed wind turbine will be given next.
This concept consists of introducing changes in the active power command signal dispatched by the wind turbine in response to the variation of some variable of interest; in this case, the network frequency. This characteristic, supported by the kinetic energy stored in the wind turbine, is implemented in the active power controller of the machine. The virtual inertia control scheme used in this study has been taken from [21] and is outlined later in Figure 8. Before describing the mechanics of the method, some aspects concerning the control criteria of the active power dispatched by the wind turbine are addressed.
In Figure 6, the characteristic P-ω curve of a three-bladed wind turbine is shown, particularized for the commercial model adopted in this study. The maximum power point tracking (MPPT) characteristic has been superimposed on the curves of the mechanical power developed by the turbine. In order to model the dynamics of the MPPT curve, the optimal power in each operating zone is characterized, as shown in the following set of equations: In Figure 6, the characteristic P-ω curve of a three-bladed wind turbine is shown, particularized for the commercial model adopted in this study. The maximum power point tracking (MPPT) characteristic has been superimposed on the curves of the mechanical power developed by the turbine. Its representation derives from some of the author's proposal [21], consisting of a simplified electromechanical model of a variable speed wind turbine conceived particularly for the study that is the object of this work. The relevant technical information for the implementation of the wind turbine model is summarized in Appendix C. Given that the objective of this study is based on the application of virtual inertia, a brief presentation of the concept and its application in a variable speed wind turbine will be given next.
This concept consists of introducing changes in the active power command signal dispatched by the wind turbine in response to the variation of some variable of interest; in this case, the network frequency. This characteristic, supported by the kinetic energy stored in the wind turbine, is implemented in the active power controller of the machine. The virtual inertia control scheme used in this study has been taken from [21] and is outlined later in Figure 8. Before describing the mechanics of the method, some aspects concerning the control criteria of the active power dispatched by the wind turbine are addressed.
In Figure 6, the characteristic P-ω curve of a three-bladed wind turbine is shown, particularized for the commercial model adopted in this study. The maximum power point tracking (MPPT) characteristic has been superimposed on the curves of the mechanical power developed by the turbine. In order to model the dynamics of the MPPT curve, the optimal power in each operating zone is characterized, as shown in the following set of equations: In order to model the dynamics of the MPPT curve, the optimal power in each operating zone is characterized, as shown in the following set of equations: Zone BC : Zone CD : From D : P opt = P max (9) where • P opt = Optimum power delivered by the wind turbine.
According to the criterion proposed in [21], the changes in the value of the setpoint active power injected by the wind turbine into the grid are carried out by shifting the MPPT curve above or below its original position, as illustrated in Figure 7. This is achieved by manipulating the optimization constant, K opt , in Equations (6)-(9), by multiplying by a factor, k VIC , which presents sensitivity to changes in the frequency of the system, ∆f, and its derivative, df/dt.
According to the criterion proposed in [21], the changes in the value of the setpoint active power injected by the wind turbine into the grid are carried out by shifting the MPPT curve above or below its original position, as illustrated in Figure 7. This is achieved by manipulating the optimization constant, Kopt, in Equations (6)-(9), by multiplying by a factor, kVIC, which presents sensitivity to changes in the frequency of the system, Δf, and its derivative, df/dt.
Depending on the deviation of the frequency from its rated value, the multiplication factor, kVIC, will generate a series of maximum power follow-up curves, which will be called VIC curves (Virtual Inertia Control curves). Next, the mechanics of the method are explained, using Figure 7. Assuming that the wind turbine is subjected to a constant wind speed (v = 8 m/s), according to the MPPT algorithm, its operating point on the P-ω characteristic of the machine will be located at point A. If, at a certain moment, the network experiences an increase in the demanded Depending on the deviation of the frequency from its rated value, the multiplication factor, k VIC , will generate a series of maximum power follow-up curves, which will be called VIC curves (Virtual Inertia Control curves).
Next, the mechanics of the method are explained, using Figure 7. Assuming that the wind turbine is subjected to a constant wind speed (v = 8 m/s), according to the MPPT algorithm, its operating point on the P-ω characteristic of the machine will be located at point A. If, at a certain moment, the network experiences an increase in the demanded power (or a loss of generated power), the frequency will begin to decline. To make the wind turbine support the grid frequency, it must respond with an increase in the active power it injects.
This method of virtual inertia forces the MPPT curve to shift upwards to drive the wind turbine operation to point O immediately. In this new operating condition, the rotor decelerates, and the kinetic energy stored in its rotating masses is released to support the restoration of the grid frequency. During the deceleration process, the operating point moves from O to B along the modified MPPT characteristic.
Once the frequency support is no longer needed, the wind turbine must return to its optimal operating point, A, therefore, the method will return the MPPT curve to its original position. The control action provided by the virtual inertia method described here, whose deductive process is detailed in [21], is defined in the following equation: • W VIR = weighting factor, defined as the ratio of the kinetic energy developed by the wind turbine to the maximum kinetic energy that it can store, as follows: Figure 8 shows a general schematic representation of the virtual inertia algorithm that must be implemented in the wind turbine speed controller to enable its participation in frequency control tasks. The details of the implementation of the scheme and the values of its parameters are described in detail in [21].
wind turbine support the grid frequency, it must respond with an increase in the active power it injects.
This method of virtual inertia forces the MPPT curve to shift upwards to drive the wind turbine operation to point O immediately. In this new operating condition, the rotor decelerates, and the kinetic energy stored in its rotating masses is released to support the restoration of the grid frequency. During the deceleration process, the operating point moves from O to B along the modified MPPT characteristic.
Once the frequency support is no longer needed, the wind turbine must return to its optimal operating point, A, therefore, the method will return the MPPT curve to its original position. The control action provided by the virtual inertia method described here, whose deductive process is detailed in [21], is defined in the following equation: dt (16)  WVIR = weighting factor, defined as the ratio of the kinetic energy developed by the wind turbine to the maximum kinetic energy that it can store, as follows: Figure 8 shows a general schematic representation of the virtual inertia algorithm that must be implemented in the wind turbine speed controller to enable its participation in frequency control tasks. The details of the implementation of the scheme and the values of its parameters are described in detail in [21].

Battery Energy Storage System
Among the potential benefits of incorporating battery energy storage systems (BESS) into electrical power systems are: improvement of the transmission of electrical energy,

Dead-band
Wash-out filter [pu] f   Figure 8. Scheme of the virtual inertia algorithm implemented in the variable speed wind turbine.

Battery Energy Storage System
Among the potential benefits of incorporating battery energy storage systems (BESS) into electrical power systems are: improvement of the transmission of electrical energy, damping of power oscillations, dynamic voltage stability, short-term power reserve, leveling of the demand profile, reduction of power supply interruptions due to low-frequency load shedding processes, sub synchronous resonance damping, and improvement of power quality [32,33].
A complementary application of interest is frequency regulation, the object of this work, in which the BESS is loaded or unloaded in response to the active power compensation needs in the system in the event of an increase or decrease of frequency, respectively [34].
In general terms, a BESS consists of three main components: the battery module, which stores electrical energy through an electrochemical process; the electronic power converter, normally referred to as an inverter, which converts direct current into alternating current to be injected or absorbed from the system as appropriate; and the power control systems, whose function is to establish active and reactive power setpoints based on the role adopted by the BESS within the electrical power system. A typical general diagram of a BESS is illustrated in Figure 9.
eling of the demand profile, reduction of power supply interruptions due to low-frequency load shedding processes, sub synchronous resonance damping, and improvement of power quality [32,33].
A complementary application of interest is frequency regulation, the object of this work, in which the BESS is loaded or unloaded in response to the active power compensation needs in the system in the event of an increase or decrease of frequency, respectively [34].
In general terms, a BESS consists of three main components: the battery module, which stores electrical energy through an electrochemical process; the electronic power converter, normally referred to as an inverter, which converts direct current into alternating current to be injected or absorbed from the system as appropriate; and the power control systems, whose function is to establish active and reactive power setpoints based on the role adopted by the BESS within the electrical power system. A typical general diagram of a BESS is illustrated in Figure 9. In order to incorporate the model that represents the dynamics of the energy storage unit within the scheme considered for the representation of the electrical power system, in the time horizon covered by the frequency control actions, the model proposed by [35] is adopted, as illustrated in Figure 10. In this, the virtual inertia control is implemented by adding two control loops: one dependent on the frequency variation, ∆f, like a powerfrequency droop control, and the other, dependent on the frequency variation rate, in the time domain, df/dt, or RoCoF. Thus, when the system frequency varies, this scheme will respond with a change in the active power setpoint ΔPBESS in order to emulate an inertial response.
A common practice in this type of study is to represent the dynamics of the battery by means of a first-order transfer function, with a TBESS time constant depending on its technology. The values of the parameters of the scheme of Figure 10 are summarized in Appendix D. In order to incorporate the model that represents the dynamics of the energy storage unit within the scheme considered for the representation of the electrical power system, in the time horizon covered by the frequency control actions, the model proposed by [35] is adopted, as illustrated in Figure 10. In this, the virtual inertia control is implemented by adding two control loops: one dependent on the frequency variation, ∆f, like a powerfrequency droop control, and the other, dependent on the frequency variation rate, in the time domain, df/dt, or RoCoF. Thus, when the system frequency varies, this scheme will respond with a change in the active power setpoint ∆P BESS in order to emulate an inertial response.
work, in which the BESS is loaded or unloaded in response to the active power compensation needs in the system in the event of an increase or decrease of frequency, respectively [34].
In general terms, a BESS consists of three main components: the battery module, which stores electrical energy through an electrochemical process; the electronic power converter, normally referred to as an inverter, which converts direct current into alternating current to be injected or absorbed from the system as appropriate; and the power control systems, whose function is to establish active and reactive power setpoints based on the role adopted by the BESS within the electrical power system. A typical general diagram of a BESS is illustrated in Figure 9. In order to incorporate the model that represents the dynamics of the energy storage unit within the scheme considered for the representation of the electrical power system, in the time horizon covered by the frequency control actions, the model proposed by [35] is adopted, as illustrated in Figure 10. In this, the virtual inertia control is implemented by adding two control loops: one dependent on the frequency variation, ∆f, like a powerfrequency droop control, and the other, dependent on the frequency variation rate, in the time domain, df/dt, or RoCoF. Thus, when the system frequency varies, this scheme will respond with a change in the active power setpoint ΔPBESS in order to emulate an inertial response.
A common practice in this type of study is to represent the dynamics of the battery by means of a first-order transfer function, with a TBESS time constant depending on its technology. The values of the parameters of the scheme of Figure 10 are summarized in Appendix D. It should be noted that this linear model, common in frequency stability studies, has certain limitations since it obviates non-linearities. For example, the response capacity of the storage system to frequency variations has a certain dependence on the state of charge of the battery. With very low state of charge, the power injection capacity can be limited to values below the nominal power of the electronic converter in order to protect the battery.
Similarly, with very high state of charge, the power absorption capacity of the battery may be diminished. However, it must be taken into account that the inertial emulation, on the one hand, takes place during very short periods of time, which implies really reduced energy exchanges between the battery and the electrical energy system, and, on the other hand, as the sign of the frequency deviation is constantly changing, the energy exchange is also changing, so the average energy balance is practically zero.
Finally, to particularize the model for the case study, given that there is currently no storage unit installed on San Cristóbal Island, the one existing on a neighboring island (Isabela), in which a 660 kW system is already installed, will be taken as a reference. Thus, the participation factor of this technology within the primary frequency control scheme proposed in this study (Figure 2) is quantified as:

Results and Discussion
In this section, the performance of the control strategies designed to improve the effectiveness of the system frequency recovery actions in contingency situations are evaluated. The analysis has been divided into the following case studies: Next, the results obtained, with time domain simulations, for some variables that are of interest in this study, are illustrated: ∆f and RoCoF are illustrated in Figures 11 and 12 and are summarized in Tables 2 and 3. Energies 2021, 14, x FOR PEER REVIEW 12 of 18 Figure 11. Grid-frequency response to a sub-frequency event. Figure 12 shows a significant reduction of the RoCoF peak with respect to the value reached in the base case: 23%, 27%, and 48% for cases 1, 2, and 3, respectively. From these results it can be inferred that the participation, whether of wind generation or BESS, in the frequency control is beneficial for the system, the advantages being more evident when these two agents carry out their work in combination.  Figure 13 illustrates the active power response of the group of wind turbines in each of the scenarios. Note that the wind turbines respond with an increase in active power in the first moments of the contingency (in cases 1 and 3), which contributes to reduce the maximum frequency deviation and the maximum value of RoCoF (Figures 11 and 12). In this figure, the curves that represent the base case and case 2 are overlapping, since in Figure 11. Grid-frequency response to a sub-frequency event.
Energies 2021, 14, x FOR PEER REVIEW 12 of 18 Figure 11. Grid-frequency response to a sub-frequency event. Figure 12 shows a significant reduction of the RoCoF peak with respect to the value reached in the base case: 23%, 27%, and 48% for cases 1, 2, and 3, respectively. From these results it can be inferred that the participation, whether of wind generation or BESS, in the frequency control is beneficial for the system, the advantages being more evident when these two agents carry out their work in combination.  Figure 13 illustrates the active power response of the group of wind turbines in each of the scenarios. Note that the wind turbines respond with an increase in active power in the first moments of the contingency (in cases 1 and 3), which contributes to reduce the maximum frequency deviation and the maximum value of RoCoF (Figures 11 and 12). In this figure, the curves that represent the base case and case 2 are overlapping, since in these cases the participation of the wind turbines in the frequency control tasks is not considered. Figure 14 shows the active power response of the BESS for the scenarios studied. Given that in the base case and in case 1 the action of this agent in the frequency support is not expected, the resulting curves are overlapped at a zero level. However, in cases 2 and 3, the active power contribution of the BESS contributes to the frequency recovery efforts. The power injection value is lower when this unit acts in coordination with the    Figure 11 shows that, by promoting the participation of wind turbines and BESS in frequency control tasks by adding virtual inertia control schemes, a significant reduction in the peak of frequency variation following the disturbance is obtained, with values around 45% in case 1, 41% in case 2, and 57% in case 3. Figure 12 shows a significant reduction of the RoCoF peak with respect to the value reached in the base case: 23%, 27%, and 48% for cases 1, 2, and 3, respectively. From these results it can be inferred that the participation, whether of wind generation or BESS, in the frequency control is beneficial for the system, the advantages being more evident when these two agents carry out their work in combination. Figure 13 illustrates the active power response of the group of wind turbines in each of the scenarios. Note that the wind turbines respond with an increase in active power in the first moments of the contingency (in cases 1 and 3), which contributes to reduce the maximum frequency deviation and the maximum value of RoCoF (Figures 11 and 12). In this figure, the curves that represent the base case and case 2 are overlapping, since in these cases the participation of the wind turbines in the frequency control tasks is not considered.
results it can be inferred that the participation, whether of wind generation or BESS, in the frequency control is beneficial for the system, the advantages being more evident when these two agents carry out their work in combination.  Figure 13 illustrates the active power response of the group of wind turbines in each of the scenarios. Note that the wind turbines respond with an increase in active power in the first moments of the contingency (in cases 1 and 3), which contributes to reduce the maximum frequency deviation and the maximum value of RoCoF (Figures 11 and 12). In this figure, the curves that represent the base case and case 2 are overlapping, since in these cases the participation of the wind turbines in the frequency control tasks is not considered. Figure 14 shows the active power response of the BESS for the scenarios studied. Given that in the base case and in case 1 the action of this agent in the frequency support is not expected, the resulting curves are overlapped at a zero level. However, in cases 2 and 3, the active power contribution of the BESS contributes to the frequency recovery efforts. The power injection value is lower when this unit acts in coordination with the wind generation (case 3) than when it does it alone (case 2).  Figure 14 shows the active power response of the BESS for the scenarios studied. Given that in the base case and in case 1 the action of this agent in the frequency support is not expected, the resulting curves are overlapped at a zero level. However, in cases 2 and 3, the active power contribution of the BESS contributes to the frequency recovery efforts. The power injection value is lower when this unit acts in coordination with the wind generation (case 3) than when it does it alone (case 2).  In Figure 15, the temporal response of the synchronous generation (diesel groups) is presented, in terms of power, for the four simulated cases. In the base case, the frequency control actions are governed by synchronous generation through the actuation of their speed controllers in response to frequency deviations, Δf. With the participation of wind generation, BESS or both in frequency control tasks, the frequency stabilization and re- In Figure 15, the temporal response of the synchronous generation (diesel groups) is presented, in terms of power, for the four simulated cases. In the base case, the frequency control actions are governed by synchronous generation through the actuation of their speed controllers in response to frequency deviations, ∆f. With the participation of wind generation, BESS or both in frequency control tasks, the frequency stabilization and recovery process improves compared to the base case due to the extra contribution of active power provided by the other agents ( Figure 11). It means slower performance of the speed controllers of the synchronous generation, since, by decreasing ∆f, which is the control variable, their part of active power increase is now a shared task with the other participating agents. This situation is reflected, of course, in a process of restoring the frequency in a slightly longer time compared to the base case. In Figure 15, the temporal response of the synchronous generation (diesel groups) is presented, in terms of power, for the four simulated cases. In the base case, the frequency control actions are governed by synchronous generation through the actuation of their speed controllers in response to frequency deviations, Δf. With the participation of wind generation, BESS or both in frequency control tasks, the frequency stabilization and recovery process improves compared to the base case due to the extra contribution of active power provided by the other agents ( Figure 11). It means slower performance of the speed controllers of the synchronous generation, since, by decreasing Δf, which is the control variable, their part of active power increase is now a shared task with the other participating agents. This situation is reflected, of course, in a process of restoring the frequency in a slightly longer time compared to the base case. Next, the benefits offered by the frequency control strategies implemented in the different generation agents that make up the electrical test system (base case, cases 1, 2, and 3) when it is subjected to more realistic operating conditions, characterized by a high variability of the wind resource, are evaluated. To do this, a wind profile like the one shown in Figure 16 has been applied to the three wind turbines. This time series of wind speed has been obtained from real records taken in the San Cristóbal island wind farm in a windy day. As in the previous study, the contingency takes place 10 s after the start of the simulation. Next, the benefits offered by the frequency control strategies implemented in the different generation agents that make up the electrical test system (base case, cases 1, 2, and 3) when it is subjected to more realistic operating conditions, characterized by a high variability of the wind resource, are evaluated. To do this, a wind profile like the one shown in Figure 16 has been applied to the three wind turbines. This time series of wind speed has been obtained from real records taken in the San Cristóbal island wind farm in a windy day. As in the previous study, the contingency takes place 10 s after the start of the simulation. The results presented in Figures 17 and 18 show that, under a variable wind scenario, the injected wind power introduces continuous power imbalances in the system, which translates into high frequency fluctuations under normal operating conditions (base case). A similar thing happens with the dynamics of the RoCoF. The results presented in Figures 17 and 18 show that, under a variable wind scenario, the injected wind power introduces continuous power imbalances in the system, which translates into high frequency fluctuations under normal operating conditions (base case). A similar thing happens with the dynamics of the RoCoF. The results presented in Figures 17 and 18 show that, under a variable wind scenario, the injected wind power introduces continuous power imbalances in the system, which translates into high frequency fluctuations under normal operating conditions (base case). A similar thing happens with the dynamics of the RoCoF.  In these figures, the effect of the participation of the different generation agents in the frequency control (cases 1, 2, and 3) can be appreciated, which produces a notable reduction in the variability of the system frequency and the RoCoF following the contingency. This is because the control strategies implemented in the wind turbines, the BESS, and the speed controller of the synchronous generation itself are permanently acting due to the existence of a constant deviation of the system frequency, not only in contingency conditions. This situation can be verified in Figures 19 and 20, which show the active power injected into the system by the wind farm and by the BESS, respectively. The power curves shown in Figure 19 for the base case and case 2 are overlapped because the participation  The results presented in Figures 17 and 18 show that, under a variable wind scenario, the injected wind power introduces continuous power imbalances in the system, which translates into high frequency fluctuations under normal operating conditions (base case). A similar thing happens with the dynamics of the RoCoF.  In these figures, the effect of the participation of the different generation agents in the frequency control (cases 1, 2, and 3) can be appreciated, which produces a notable reduction in the variability of the system frequency and the RoCoF following the contingency. This is because the control strategies implemented in the wind turbines, the BESS, and the speed controller of the synchronous generation itself are permanently acting due to the existence of a constant deviation of the system frequency, not only in contingency conditions. This situation can be verified in Figures 19 and 20, which show the active power injected into the system by the wind farm and by the BESS, respectively. The power curves shown in Figure 19 for the base case and case 2 are overlapped because the participation In these figures, the effect of the participation of the different generation agents in the frequency control (cases 1, 2, and 3) can be appreciated, which produces a notable reduction in the variability of the system frequency and the RoCoF following the contingency. This is because the control strategies implemented in the wind turbines, the BESS, and the speed controller of the synchronous generation itself are permanently acting due to the existence of a constant deviation of the system frequency, not only in contingency conditions. This situation can be verified in Figures 19 and 20, which show the active power injected into the system by the wind farm and by the BESS, respectively. The power curves shown in Figure 19 for the base case and case 2 are overlapped because the participation of the wind farm in the frequency control is not contemplated. A similar fact occurs in Figure 20 regarding the power injected by the BESS for the base case and case 1.  Figure 20 regarding the power injected by the BESS for the base case and case 1. Note in Figure 17 that when the contingency occurs, the drop in frequency is successfully mitigated in all the cases studied. However, reducing both the frequency deviation and the maximum value of RoCoF is more effective in case study 3, which considers the combined participation of wind generation and BESS in the provision of this ancillary service. Similar conclusion is reached in the study under constant wind conditions (Figures 11  and 12).   Figure 17 that when the contingency occurs, the drop in frequency is successfully mitigated in all the cases studied. However, reducing both the frequency deviation and the maximum value of RoCoF is more effective in case study 3, which considers the combined participation of wind generation and BESS in the provision of this ancillary service. Similar conclusion is reached in the study under constant wind conditions (Figures 11 and 12).  Finally, Figure 21 shows the dynamics of the power provided by synchronous generators, where it is possible to see how the demand for their participation in frequency control tasks is reduced when the other generation agents (wind, BESS or both) do too.

Conclusions
This work presents a study of the effect of two control strategies based on the concept of virtual inertia (through batteries and the incorporation of a virtual inertia algorithm in wind generation units) on the frequency regulation in an electric power system, within Finally, Figure 21 shows the dynamics of the power provided by synchronous generators, where it is possible to see how the demand for their participation in frequency control tasks is reduced when the other generation agents (wind, BESS or both) do too.  Finally, Figure 21 shows the dynamics of the power provided by synchronous generators, where it is possible to see how the demand for their participation in frequency control tasks is reduced when the other generation agents (wind, BESS or both) do too.

Conclusions
This work presents a study of the effect of two control strategies based on the concept of virtual inertia (through batteries and the incorporation of a virtual inertia algorithm in wind generation units) on the frequency regulation in an electric power system, within

Conclusions
This work presents a study of the effect of two control strategies based on the concept of virtual inertia (through batteries and the incorporation of a virtual inertia algorithm in wind generation units) on the frequency regulation in an electric power system, within the time frame of primary frequency control under a scenario of high penetration of wind generation. By introducing each of the control strategies, significant improvements in the system's frequency response can be evidenced in a disturbing event, mainly on the variables: frequency and RoCoF. The simulation results also reveal that such improvement is much more noticeable when the two strategies are applied in combination. This fact allows to corroborate the fulfillment of the research hypothesis.
The proposed model is outlined as an analysis and diagnostic tool that allows studying the dynamics of certain variables of interest of an electrical power system of these characteristics, offering the network operator a tool to study different scenarios that pose a risk to the network frequency stability and being able to take the necessary measures that lead to its prevention or correction.
It is important to note that this is only a first step of the analysis that aims at seeing the big picture of the interaction among different elements of a power system. The advantage of this load frequency control approach is that it allows the analysis of sensitivity issues and the drawing of general conclusions, allowing to define the general behavior of the system and to tune the controllers accordingly. However, this approach has some limitations related to the necessary simplifications, e.g., a unique frequency is considered, parts of the models are described as small-signal ones, and simplified controllers are used. This can lead to inaccuracies for some variables that can have an impact on the detailed description of the system behavior when subject to some particular disturbances. A second or complementary step of the analysis, out of the scope of this work because of its specificity, is based on the approach of detailed computational simulations, even including hardware-in-the-loop or real-time simulations. Although it is difficult to obtain general conclusions with this second approach, it allows a more accurate modeling and is adequate for studying specific situations of interest in a given case, thus complementing the load frequency control approach described in the paper.

Conflicts of Interest:
The authors declare no conflict of interest. Table A1. Technical data of the diesel generator [26].

Characteristic
Value and Unit

Parameter Symbol Value and Units
Rated power (base power) P AG,base 800 kW Maximum/minimum generator power P g,max /P