Allocation of Spinning Reserves in Autonomous Grids Considering Frequency Stability Constraints and Short-Term Solar Power Variations

Low-inertia, isolated power systems face the problem of resiliency to active power variations. The integration of variable renewable energy sources, such as wind and solar photovoltaic, pushes the boundaries of this issue further. Higher shares of renewables requires better evaluations of electrical system stability, to avoid severe safety and economic consequences. Accounting for frequency stability requirements and allocating proper spinning reserves, therefore becomes a topic of pivotal importance in the planning and operational management of power systems. In this paper, dynamic frequency constraints are proposed to ensure resiliency during short-term power variations due to, for example, wind gusts or cloud passage. The use of the proposed constraints is exemplified in a case study, the constraints being integrated into a mixed-integer linear programming algorithm for sizing the optimal capacities of solar photovoltaic and battery energy storage resources in an isolated industrial plant. Outcomes of this case study show that reductions in the levelized cost of energy and carbon emissions can be overestimated by 8.0% and 10.8% respectively, where frequency constraints are neglected. The proposed optimal sizing is validated using time-domain simulations of the case study. The results indicate that this optimal system is frequency stable under the worst-case contingency.


I. INTRODUCTION
The depletion of mature oil and gas (O&G) fields reduces the energy return on investment of the field and increases the emission of nitrogen oxides and greenhouse gases (GHG) [1]. The use of renewable energy sources (RESs) in O&G operations has therefore become an active research area in recent years. Groups have analyzed this challenge from a broad range of perspectives [2], [3], [4], [5], [6], [7], [8], [9]. It is, however, widely accepted that the planning or operation of hybrid energy systems in the O&G industry requires the solution of an optimization problem. The objective of such optimization can differ considerably. The power balance of the energy system is, however, a physical constraint that must always be satisfied.
Most O&G plants, platforms, and vessels are supplied ac power by low-inertia, isolated systems [10]. Large active power variations and considerable deviations from the rated frequency can, however, occur where the penetration of RESs in such systems is high. This is a problem that is also found in other autonomous power systems (APS), such as islands [11], [12] and community microgrids [13]. Frequency stability constraints and the proper sizing of spinning reserves, such as frequency containment reserves (FCR) and frequency restoration reserves (FRR), are therefore important dimensional aspects and should be given special attention in optimization algorithms [14].
Short-term variations caused by cloud passage across solar photovoltaic (PV) parks or wind gusts at wind farms can furthermore saturate the ramping capabilities of other generators. The relationship between frequency stability, RESs penetration and variability, size of energy storage systems (ESSs), and scheduling decisions, therefore, is a topic that should be investigated in more detail in such APS. Large interconnected and meshed electricity networks also are progressively behaving as low-inertia systems, as the energy transition progresses worldwide and penetration of RESs increases [15], [16]. Interest in the topic of this paper therefore is wide.

A. LITERATURE REVIEW
Reliability-constrained unit commitment (UC) is described, for example, in [17], [18], and [19], and addresses the UC problem during unexpected events such as the severe loss of load, generation or transmission capacity. The strategies VOLUME 11, 2023 applied to mitigate this include n+1 redundancy, and spinning reserve for frequency containment and restoration. Reference [20] for example proposed a two-stage stochastic mixed integer linear programming (MILP) formulation for the calculation of UC and reserve scheduling, at frequency deviations of 0 and ±10 mHz. Reference [21] also discusses formulations for addressing wind power variations and the correct allocation from ESSs of spinning reserves. Generator contingency was, in this, considered to be the largest wind power fluctuation from one optimization time-step to another. Power balances were handled before and after contingencies in these works, and without generator transient period dynamics or ramping capacities being taken into account.
Reference [22] has proposed frequency-constrained UC problem implementation alternatives to circumvent these limitations. References [23], [24] thoroughly discuss the challenges of a frequency-constrained model. Options such as the use of maximum ramping capacity or minimum value of inertia were evaluated, but can lead to a nonlinear problem. This could, however, be handled by a Benders decomposition. Reference [25] introduced a nonlinear model that includes frequency security constraints and the dynamics of the transient period. This, however, requires the use of a genetic algorithm to solve the optimization problem. Reference [6] also use a genetic algorithm to optimally size the energy system of an offshore platform interconnected with a wind farm. The maximum ramping capacity of generators and frequency stability limits were furthermore evaluated using time-domain simulations of worst-case scenarios for each candidate solution.
Reference [26] integrated the power drops observed in solar PV, within 15-minute time windows, into the constraints of an optimal planning problem, this accounting for the stochastic variation of solar irradiance. A recent work [27], however, highlighted the varying duration and magnitude of the short-term variability of PV plants, and also that maximum perturbation may rise to 60 % of the plant's rated power in less than 30 s. These short-term variations must be taken into account if realistic decisions on storage investment or unit scheduling in optimization algorithms are to be achieved. Reference [14] furthermore provided an updated review of the frequency-constrained UC problem, a classification of the different approaches documented in the literature, and a discussion of their many shortcomings. Also proposed a linear model that includes frequency constraints directly in a MILP algorithm, and that therefore does not require the use of external time-domain simulations in the assessment of solution security. They did not, however, consider the effect of fast short-term variations in RESs, combined with generator ramping restrictions and use of ESSs on frequency stability.
Research gaps therefore still exist in the efficient sizing, planning, and operation of autonomous, low-inertia power systems with high RESs penetration, especially wind and solar PV.

B. PAPER CONTRIBUTIONS
This paper addresses some of the modeling challenges and shortcomings described above, and introduces a set of algebraic frequency stability constraints that can be directly applied to linear optimization formulations of the UC problem in low-inertia APS with high RESs penetration. Practical issues such as limited generator ramping capacities and variability of RESs are taken into consideration. Short-term frequency variations are also used to reduce the conservatism of a UC under uncertainty, which permits complementary FCR and FRR action and a reduction in the power required from ESSs.
These ideas are exemplified in a case study of the integration of a solar PV plant and a battery ESS into an existing isolated industrial installation fed by gas turbines (GTs). A MILP algorithm is used to solve a rolling frequencyconstrained UC problem, and to optimally size the installed capacity of solar PV and battery ESS. The RES short-term variability is evaluated by using a method recently published in [27], that identifies equivalent solar ramping scenarios during cloud passage. The novel proposal, described in detail in Section II-C, are linear constraints for frequency stability that model tolerated frequency dynamics and the complementary FCR and FRR action taking place during shortterm power variations, such as cloud passage events in solar PV installations or wind gusts in wind farms. This allows a significant reduction in the amount of FCR and consequently the rated power of an ESS designed to support the grid in such events, as demonstrated by the results of the case study in Section III-B. The frequency-stability results for the constrained and non-constrained implementation were compared, and major differences were highlighted, being the solution data obtained from this also used in a timedomain simulation to validate the results under the worst-case scenario.

II. ADDRESSING FREQUENCY STABILITY WITH LINEAR CONSTRAINTS
This section reviews the equations that describe the frequency dynamics of a power system, and the types of spinning reserves required to ensure frequency stability. The rationale for determining spinning reserves by obtaining a set of algebraic expressions from the dynamic equations, and applying these as constraints to the linear optimization problem, is also introduced. The first proposed formulation treats the sudden disconnection of loads or generators. The next formulation addresses the problem of the short-term variation of loads or RESs, this requiring gradual compensation by dispatchable generators with limited ramping capacity.

A. FREQUENCY STABILITY AND SPINNING RESERVES
The simplified model of a flywheel spinning at angular speed ω s can be used to express the dynamics of frequency in a power system. The rotating masses of all synchronous generators and motors are represented by an equivalent moment of inertia J . Generators deliver, at one end of the shaft, energy to the flywheel via torque T G . Load removes energy through T L at the other shaft end. The natural and controlled damping of the system is represented by coefficient B. Applying Newton's second law of motion to this simplified model, and multiplying both sides by the angular speed ω, gives Eq. (1) [28], where P G and P L are respectively the active power delivered by generators and consumed by loads.
Eq. (1), which is also known as the Swing Equation, represents the power balance that is required to maintain this system at its rated angular speed ω s . The frequency deviatioñ f = ω−ω s ω s can, after this equation has been normalized by ω s and by the total apparent power of the generators S n , be estimated by Eq. (2) [28], where M , D, p G , and p L are respectively the normalized inertia, damping, and active powers delivered by generators and consumed by loads.
The system is considered to be frequency stable [29] where any power imbalancef remains within intervals ±r tr during transient conditions and ±r ss during steady-state conditions. These intervals are often defined in grid codes or industry standards. In Europe, r ss = ±1% and r tr = ±3% are specified for generators connected to transmission [30] or distribution systems [31], [32]. Wider limits are normally specified for low-inertia, autonomous systems, for example r ss between 2 and 5%, and r tr between 5 and 10% for fixed and mobile offshore units [33].
The minimum controlled damping D min that generators must provide when assuming a steady-state in Eq. (2) (ḟ = 0), andf < 1% and null natural damping from loads, can be approximated by Eq. (3,) where P b = (p G − p L ) represents the maximum continuous imbalance that the power system can withstand and still maintain frequency stability. P b also defines the minimum level of spinning reserves that are required for frequency variations to be kept within the permitted range of ±r ss . This is therefore referred to as frequency containment reserves (FCR) [34]. The assumptions used to derive Eq. (3) may underestimate D min whenf values of more than 1% are possible [35]. A more robust approximation may therefore be given by Eq. (4).
Note that the FCR strategy (proportional control) cannot restore frequency to its rated value. Achieving the zero steady-state error (f = 0) therefore requires p G and p L to be matched by generator re-dispatching, through integral control [28], [36]. This requires additional spinning reserves, which are referred to as frequency restoration reserves (FRR) [34]. The action of FRR, FCR constraining frequency during disturbances to defined bands, therefore can be slow. This is desirable, not just to avoid new distur-bances, but also to accommodate the ramping limitations of generators.
A single generator can provide both FCR and FRR. This may not, however, be the optimal solution from a technical or an economic perspective. GTs can, for example, usually provide both types of reserves simultaneously. Excessive activation of FCR on GTs can, however, increase actuator wear and tear and increase maintenance costs. More suitable FCR could be batteries designed to supply large amounts of power for short periods of time. Their participation in FRR will, however, usually require large energy storage capacity, which will increase investment costs and space requirements.
The optimal allocation of spinning reserves can, even for small power systems, quickly become a complex process. An optimization algorithm is therefore often used in planning and operational decisions. The following sections present, for such optimization problems, the proposed algebraic constraints which can determine the optimal allocation of FCR and FRR based on the frequency dynamics of a power system.

B. FORMULATION 1: CONSTRAINTS FOR SUDDEN LOAD OR GENERATION LOSSES
Eqs. (5) and (6) present the constraints that must be met by generators participating in FCR when considering the assumptions and conditions imposed by Eq. (3), where P FCR m,h , P min m , P max m , P m,h are respectively the assigned FCR contribution, the minimum, the maximum, and the current commitment of generator m in time step h. The sum of FCR provided by all generators is to also, at each time step h, be greater than the worst-case power disturbance P b,h , as shown in Eq. (7).
P b,h varies with time, because it depends on the worst expected sudden load or generation loss for the current state of the power system. This value is normally, in a security assessment, obtained from the evaluation of possible contingencies. This complex topic is, however, beyond the scope of this work, a simplified approach however being presented in Section III-A1 and used in this paper. A deeper discussion about security assessments can be found in [37]. Note that Eq. (5) can be rewritten based on the assumptions and conditions imposed by Eq. (4), wheref values of more than 1% are permitted.

C. FORMULATION 2: CONSTRAINTS FOR SHORT-TERM POWER VARIATIONS
In this paper, short-term power variations are considered to be power deficits or surpluses that can be compensated for by generators that provide FRR without infringing their permitted ramping rates, and by generators that provide FCR without infringing the permitted frequency variations in VOLUME 11, 2023 steady-state. Examples of these variations in plants fed by RES include cloud passage across solar PV parks or wind gusts at wind farms.
This formulation requires the following assumptions: 1) A power system in balance and in steady-state before RES variations at t = t 0 occur. In other words, in Eq.
(2),f (t 0 ) = 0,ḟ (t 0 ) = 0, and p G (t 0 ) = p L (t 0 ). 2) Worst-case short-term variations can be approximated by a maximum amplitude p RES , a minimum duration T RES , and a constant rate of change rr RES = p RES / T RES . This assumption has been corroborated by recently published analysis of shortterm wind [38] and solar [27] variability, the determination of these parameters from high-resolution wind speed and solar irradiance measurements also being discussed in detail in these papers. 3) p RES is compensated for without infringing the maximum ramping rate rr FRR m of the generators participating in FRR. 4) Short-term RES power variations are sufficiently smooth. It is therefore reasonable to assume minimal influence of system inertia (Mḟ ≈ 0) during interval T RES .
Assumption 1 implies that any power variation in the system can be arbitrarily attributed to either p L or p G in Eq. (2). Eq. (8) is obtained based on the assumption that RES variations affect p L and on Assumption 2. Eq. (9) also accounts for the contributions of all generators that provide FRR in term p G , being Assumption 3 here applied. Substituting Eqs. (8) and (9) into Eq. (2) and using Assumption 4 gives Eq. (10).
An algebraic constraint that allocates FCR and FRR for compensation of short-term RES power variations is finally obtained in Eq. (11) (11) Note that the constraint in Eq. (11) implies that FCR and FRR are activated concurrently when rr RES h > m rr FRR m . It is therefore implicitly assumed that sudden load variations and simultaneous worst-case RES power variations were taken into consideration when determining P b,h in FCR sizing. This leads to the additional constraints expressed in Eqs. (12) and (13), i.e., that generators m participating in FRR have a total available up and down capacity of P b,h at each time step h. ∀h, P b,h may eventually become asymmetrical, the positive worst-case power variation being different from the negative worst-case power variation. The up and down-ramping capacities of generators can also be distinct, FCR and FRR being asymmetrical in this situation. Separate parameter sets r ss , r tr , should then be considered for up and down reserves. The formulations proposed in this and the previous section can be developed further for this general case but will not, in the interest of brevity, be presented in this paper.
Regarding Assumption 4, note that short-term power variations of RES typically occur in ramp, as described later in Section III-A1, and detailed in Appendix A, [27] for solar PV, and [38] for wind power. A cloud passage would likely take many seconds to complete shade the total area of a MW-scale plant, meaning that the power output would drop approximately in a ramp. Where such ramps are smooth enough, the contribution of the GTs inertia for the power balance represented in Eq. (2) would indeed be minimal, as turbine governors would have enough time to react to the grid frequency changes and compensate the power imbalance. This is exactly what is being assumed in this formulation.
This concept is better understood in Fig. 1, which sketches what happens when a cloud passage occurs in an APS and where P 1 , P 2 , P max represent respectively the power delivered by two GTs and their maximum limit, P PV and P bat are the power delivered by a solar PV farm and a battery ESS, and the superscripts DC and PC indicates respectively the power variations during and post cloud passage.
The dynamics of short-term power variations are, in summary, different from sudden load or generation loss, the latter being addressed by the constraints proposed in Section II-B. The GTs inertia will, in the latter case, be very important to guarantee the power balance represented in Eq. (2) and avoid extreme frequency variations, a topic well explored in [39], for instance.

III. CASE STUDY OF AN INDUSTRIAL INSTALLATION
The case study of an isolated O&G installation with a peak electric power demand of 160 MW is presented in this section, to exemplify the use of the proposed formulations. This plant is currently equipped with four GTs each of 45 MW. 30 MW of the total load is considered to be non-essential, and can therefore be shed during an extreme contingency. This design therefore allows the critical power demand of 130 MW to be supplied, and a n+1 redundancy for GTs. The operator would like to integrate a solar PV farm and a battery ESS into this plant, to reduce emissions of nitrogen oxides and GHG. FCR and FRR should be supplied by the existing GTs, the battery ESS being sized to provide only FCR. Fig. 2 presents the proposed architecture and Table 1 the main technical parameters.

A. MILP ALGORITHM FORMULATION
A MILP algorithm was developed to support the operator's investment decision. The objective of optimization was to find the best performing pair of solar PV installed capacity P inst PV and battery ESS installed capacity P inst bat for the plant. Eq. (14) shows the objective function, c PV , c bat being the installation costs for solar PV and battery ESS in $/kW, e denoting the discount rate of the project and c f aggregated fuel costs including CO 2 and NO x penalties. Operational costs across the plant's lifetime Y inst in years take into consideration the fuel consumption FC m,h of each generator m at each time step h, as defined in Eq. (15).
∀m, ∀h, Eqs. (16) to (25) express the system operational constraints. Eq (16) ensures that system power is in balance at each time step h. Note that the ESS does not take part in permanent load balancing as it is used only for FCR, that is ensuring frequency stability in case of short-term PV drop or large contingency. The energy flows related to charge and discharge therefore are considered negligible compared to the industrial load. The injected solar PV power P inj PV ,h is calculated in Eq. (17) using the available PV area A PV ,h , the average irradiance I h in time step h, and derating factor d PV . Eq. (18) integrates P inst PV into the objective function of Eq. (14).
Eqs. (19) to (25) Note that Eqs. (20) to (21) allocate the FCR for each generator m defined in Eqs. (5) and (6) according to a predefined damping D m , calculated based on Table 1 data. D m = P max m /(β m P base ). β m is therefore the GT droop, P base = 45 MW is an arbitrary value adopted for normalization. The total damping of the GTs m D m may not, however, be sufficient to comply with Eq. (7), battery ESS power P inst bat in this case being sized to provide the remaining FCR. The constraints necessary for this are discussed in the next section, which are based on the formulations proposed in Sections II-B and II-C.

1) BATTERY SIZING
The first step in sizing spinning reserves is to apply sudden load or generation loss constraints, a simplified security assessment in this study being used. The loss of a GT in operation was therefore considered to be the worst contingency of this type. The system was therefore designed for n+1 redundancy of the GTs, the sudden power disturbance VOLUME 11, 2023  term P sud b,h therefore being defined by Eqs. (26) and (27).
Cloud passage events at the solar PV farm should also be taken into consideration. Integration of these into the MILP formulation is however, in this study, based on the assumption that the optimization algorithm has access to a compact set H max of worst-case solar irradiance ramp events. Appendix A summarizes a procedure for building H max from high-resolution irradiance time series. Solar irradiance data provided by NREL [40] was used. The wavelet variability model [41] was also employed to take into account the geographical smoothing effect, using Sandia's PV_LIB Toolbox for Matlab [42]. The datasets used to build H max in this study are available in [43], being H max for h = 11 presented in Table 2. The duration and irradiance drop associated with ramp event r in time step h are denoted T r,h and I PV r,h . The worst-case solar PV disturbance term P PV b,h at each time step h can, where this framework is assumed, be specified using Eqs. (28) and (29).
Spinning reserves can be specified where both P sud b,h and P PV b,h have been defined and the formulations proposed in Sections II-B and II-C are applied. FRR is sized using Eqs. (30) and (31), which are obtained when the constraints in Eqs. (12) and (13) ∀h, FCR is then calculated by combining the minimum damping requirements defined in Eq. (7) and the constraint in Eq. (11), which assumes that P PV r,h is partially compensated for by FRR during a cloud passage. The battery ESS is therefore sized, at each time step h, to simultaneously compensate for any power deficit between a) P sud b,h , FCR provided by GTs; and b) P PV r,h , FRR provided by GTs for all ramp events in H max . Eqs. (32) and (33) describe these requirements. Finally, P inst bat is integrated into the objective function of Eq. (14) via Eq. (34). Note that the proposed sizing is only addressing the maximum power of the ESS, not its energy capacity. To calculate the energy flows demanded by FCR, it is necessary to  simulate the grid using more detailed models at sub-second time steps, which cannot be performed using the methodology here proposed. This limitation is later exemplified and further discussed in Section IV.

B. OPTIMAL SIZING WITH STATIC AND DYNAMIC FREQUENCY CONSTRAINTS
The MILP algorithm described in Section III-A was implemented in Gurobi 9.1, using an optimality gap tolerance of 1%. For reproducibility, the final results and values of intermediate steps of the optimization model described in this section are available in [43]. Table 3 lists the technoeconomic parameters used in the optimization. Four scenarios were considered, these being selected to highlight the importance of frequency stability constraints: 1) Baseline is the current operation of the case study plant, power generation being based on 4 GTs. 2) No FC includes the integration of the solar PV farm, but without any frequency stability conditions. The constraints in Section III-A1 are therefore omitted.

3) Static FC is the full implementation given in
Section III-A. It adopts P FCR m,h = 0 in Eq. (32), which means that the frequency dynamics and the FCR contribution during cloud passage are ignored. 4) Dynamic FC includes the FCR contribution in Eq. (32). The results are given in Table 4. The highest total costs are in the Baseline case, minimum costs being in the No FC case at a discount of 8.9 % on Baseline. This can be explained  (32). This assumption is also used in Eq. (20) and (21) in No FC. The system may, in No FC, face a power deficit of 50.2 MW. Critical PV ramps are r 2 and r 3 , which would undoubtedly lead to a grid blackout. P FCR bat,r is obtained, in the Static FC case, by adding 22.5 MW of one GT loss to 10.8 MW for the imbalance between FRR and PV ramp rate. The generator loss is fully compensated for in the Dynamic FC case, by the damping of the remaining GTs. This shows the importance of taking into consideration the damping capabilities of generators participating in FCR and the tolerated frequency deviations.

C. VALIDATION WITH TIME-DOMAIN SIMULATIONS
Time-domain simulations implemented in OpenModelica v1.19.2 were used to corroborate the battery sizing by the proposed MILP algorithm, [43] containing the OpenModelica model used for the validation described in this section and the Python script used to plot Fig. 5. Fig. 4 presents a schematic diagram of the simulation model. The dynamics of the electricity grid are reduced into active power flows, Eq.
(2) being used to calculate grid frequency deviation. Battery and solar PV models consist of ideal active power sources, battery ESS adopting a droop control strategy including saturation and only providing FCR. GTs and their governors are represented by a droop model that contain a first-order filter to represent the delay of the actuator, being a time constant of 0.5 s employed. Ramp-rate saturation is added for the FRR component. Other relevant parameter values are given in Tables 1, and 4.
The validation consists of simulating the worst-case cloud passage event and the simultaneous loss of one GT, as presented in Table 2, a check of whether frequency deviation is lower than or equal to r ss then being carried out. Remark that all other cases will produce smaller frequency drops or lower FRR setpoint rate of change. The power perturbation is generated by applying a load step that corresponds to P sud b,h = 22.5 MW at t = 10 s, and by simultaneously applying a negative power ramp corresponding to r 2 in Table 2, i.e. 27.7 MW. The system inertia M is equivalent to three GTs. Fig. 5 shows  the results of this simulation, including the battery ESS, the GTs active power and the grid frequency. Battery support means the frequency drops to 49.5 Hz after the worst-case event, meeting the target of up to 0.5 Hz deviation. The grid frequency has also an oscillatory behavior after the sudden loss of one GT due to the time delay of the turbine governors, which is compensated by the fast action of the battery ESS. A further discussion of this topic is beyond the scope of this paper, but the interested reader can obtain more information in [44].

IV. DISCUSSION
Optimal solar PV and battery ESS capacities were calculated for an industrial plant, by using a novel formulation of linear frequency constraints and their integration in a MILP sizing problem. Economic and environmental performance is reduced (at +10.8% of CO 2 emissions and +8.0% of total costs) in relation to optimization results without frequency constraints. It was, however, shown that the solution without frequency constraint is not resilient to worst-case events, such as a simultaneous cloud passage and loss of a GT. Taking frequency reserves into consideration in sizing optimization therefore avoids the over-estimation of the benefits of solar power integration, and provides a more robust architecture in terms of power supply security.
The battery ESS capacity requirement was reduced by 67.6% of the static calculation value that neglects the FCR contribution. This was achieved by taking into consideration the frequency deviation tolerance (r ss = 0.5 Hz) in the Dynamic FC scenario. Considering simultaneous PV drop and generator loss allows a robust approach, but also leads to conservative and costly ESS sizing. The problem of battery ESS capacity allocation might, in less sensitive applications, take into consideration the highest value of the worst-case PV ramp event and generator contingency, but not both simultaneously. The impact of this on overall costs and CO 2 savings would, however, still be relatively small, primarily due to the marginal effect that the proposed constraints would have on the fuel consumption of GTs, which is the main cost driver in this case study.
The battery energy storage capacity in MW h for FCR is not a dominant cost driver in this problem, and the cost function of Eq. (14), as consequence, do not consider it. The actuation time of FRR is usually less than 30 s in a typical energy management system (EMS) for APS, meaning that, if full FCR capacity (i.e. 10.8 MW) is delivered for 10 events every hour (i.e. 5 minutes in total), the required energy will be 0.9 MW h for the selected battery ESS in the case study. Typical losses of ESSs in the MW-range varies between 1 and 2%, the energy required by FCR therefore being in the same order of magnitude of the battery ESS losses. Remark also that the energy required by FCR represents only 0.56% of the load demand in the case study (160 MW h) and can be assumed negligible when taking into consideration the accuracy of the model employed in the UC problem and the optimality gap tolerance adopted. Including the battery energy storage capacity in the problem formulation would however be extremely relevant if the ESS was designed to provide FRR, or if the proposed constraints were employed in an operational UC problem designed to provide dispatch setpoints for the GTs and ESS in real-time.
The attentive reader may have noticed that the fuel and CO 2 costs presented in Table 3 are well above historical market prices. The reader may have also noticed that the project discount rate and installation costs for solar PV and battery ESS are lower than the prevailing values in the industry. This was to allow a high penetration of solar PV in the case study installation, to highlight the impact of frequency stability constraints in this scenario. The results of the MILP algorithm are, of course, sensitive to the economic parameters. For example, if the discount rate is doubled, then the installed solar PV capacity will drop to below 9 MW in the Static FC and Dynamic FC scenarios. The presentation of an economic sensitivity analysis in the case study is, however, beyond the scope of this work.
Just considering Eqs (20), (21), (30) and (31) is, in terms of spinning reserves allocation, equivalent to the formulations proposed in [14] and [26]. Introducing Eq. (32), however, contributes to security sizing problems, as it draws in both formulations. FCR modeling during short-term power variations therefore improves the optimal solution. The validation using a time-domain model that only considers active power flows, shows that this frequency stability proposal can produce valid results. The power management of larger systems is, however, a multi-objective optimization problem with coupled variables, and one in which not only active power flows, but also voltages and reactive power flows, are optimized  as discussed in [19]. The inclusion of frequency constraints in such problems, to deal with short-term power variations, is however a topic for future research.
Note as well that a droop with a fixed value of 10% for all GTs was considered in the case study. As briefly discussed in Section III-A, the droop value directly affects the total system damping. It would therefore also affect the size of the battery ESS. This observation suggests that this value should also be used as a decision variable in the optimization problem. This does, however, lead to a non-linear formulation. Further research on reserve allocation models and MILP formulation should therefore be carried out to circumvent this issue.
It is worth highlighting, at this point, that the constraints proposed in Sections II-B and II-C do not controlḟ , only assign the correct amount of FCR given bounded values of (p G − p L ) and M . The proposed constraints, in other words, guarantee that the frequency deviationf remains within intervals ±r tr during transient conditions and ±r ss during steady-state conditions (i.e. post-disturbance), which is the definition of frequency stability given in [29] and in Section II-A.
Note that large values off in the negative direction (i.e. lack of power generation) can cause magnetic over-flux in electrical machines and transformers. This situation not only increases losses and heating, but also induced voltage gradient between laminations that can break down the core insulation of these equipment and cause serious damages. This topic is discussed in depth in [39], [45], [46], and [47] and can be considered a more serious issue than large values ofḟ in low-inertia APS.
Despite not being explicitly mentioned in Section II-A, M influences the term r tr in Eq. (4), which represents the maximum value off during transient conditions, also known as frequency nadir or zenith. This term can therefore be used to include constraints on M . The model from Eq. (2), however, cannot capture the dynamics of the transient period, as it does VOLUME 11, 2023 not include the time delay of actuators. References [48], [49] have discussed this issue and introduced linearizations and/or analytical solutions that can be employed to include services such as fast frequency reserves (FFR) or inertia emulation in the UC problem. The constraints described in these references can be included in the formulation proposed in Section III-A where the dynamics of the transient period must be addressed.
In the case study case presented in Section III, however, a higher penetration of solar PV is mainly prevented by thef constraints post-disturbance and the installation costs of the battery ESS, as shown in Table 4 and discussed in Section III-B. The authors therefore consider reasonable to ignore the dynamics of the transient period in this specific case. For a more general problem formulation, the constraints introduced in [48] and [49] can be combined with those presented in Section II-C. This however will turn the problem into a nonlinear optimization, as the ratio D M determines thẽ f dynamics but each variable will become a decision variable in the optimization. A relevant topic for future work is, therefore, to address this issue and develop a convex reformulation of the problem that can be solved by MILP.

V. CONCLUSION
The problem of spinning reserves allocation for isolated grids with high penetration of RESs is addressed in this work. Linear frequency stability constraints were first formulated, to allow integration in a frequency constrained unit commitment (UC) problem, the proposed constraints ensuring the allocation of enough reserves to compensate for shortterm renewable variations and generator contingencies. This includes practical issues such as the limited ramping capacities of frequency restoration reserves (FRR). It also takes advantage of short-term frequency variations and the complementary action between frequency containment reserves (FCR) and frequency restoration reserves (FRR).
The reserve allocation strategy is exemplified by a case study in which a mixed integer linear programming (MILP) algorithm is formulated for the optimal sizing of a solar photovoltaic (PV) farm and an energy storage system (ESS). Linear frequency stability constraints are aggregated to ensure the resiliency of the grid, in the event of cloud passage and generator contingency occurring simultaneously. The results show that this method provides an optimal and secure architecture. It also shows that neglecting stability constraints leads to an inoperable solution, and overestimates system profitability and emission savings, by 8.0 and 10.8% respectively.

APPENDIX A GENERATION OF SOLAR PV POWER DROP SCENARIOS
Solar variability is typically measured by irradiance sensors, which generate high-resolution time series. The measurement sampling time must be lower than 5 s if cloud passage events are to be captured. Integrating these time series into a high-level optimization is however impractical, due to time increment discrepancies. Fig. 6 summarizes a procedure for extracting ramp scenarios from irradiance sensor time series, as was recently proposed in [27] and applied in Section III.
In step 1, high resolution data from sensors is decomposed into hourly time slices. In step 2, all ramp events within an hour slice are detected. A ramp event rr RES h,r is defined by its irradiance drop I r,h and ramp duration T r,h during cloud passage. In step 3, all events are grouped in a compact set. The hourly convex hull H max h is defined to gather the set of highest irradiance drops for each duration, to allow only the worst-case events to be extracted. This provides a subset of a limited number of ramps. Steps 2 and 3 are then repeated for each hour slice defined in step 1. Finally From 1992 to 2001, he was a Research Engineer with the Centre for Energy, Mines Paris-PSL, where he worked on modeling and control of energy systems and led and participated in several research projects at national and European levels. In 2001, he became a Professor with the Centre for Applied Mathematics, Mines Paris-PSL, where he is currently the Deputy Director and the Head of the post-graduate master program in energy systems optimization. He has published nine books and more than 100 journal and conference papers in the field of energy systems, energy economics, and applied mathematics. His research interests include energy systems optimization, modeling, control and forecast, energy economics, and applied mathematics. From 2009 to 2011, she was a Postdoctoral Fellow with the Norwegian University of Science and Technology (NTNU), Trondheim, Norway, where she worked on the grid integration of offshore renewable energies. She was a Researcher with Tecnalia, Bilbao, Spain, from 2011 to 2013, where she was the Principal Investigator in the FP7-Sea2grid Project, related to storage needs for the grid integration of wave energy converters. From 2013 to 2014, she was a Research Scientist with SINTEF Energy, Trondheim, and an Adjunct Associate Professor with NTNU. In 2014, she became a Full Professor of offshore grids with the Department of Electric Power Engineering, NTNU. Since 2020, she has been a Full Professor with the Department of Industrial Engineering, University of Trento, Italy. She has published two book chapters and more than 150 journal and conference papers in the field of marine energy and energy conversion systems. Her research interests include the design and control of energy conversion and transmission and distribution systems, with a focus on offshore energy and power-quality issues.
Dr. Tedeschi was a recipient of the Marie Curie Fellowship and a member of NTNU's Outstanding Academic Fellows Program.