Optimal Energy Management in Autonomous Power Systems with Probabilistic Security Constraints and Adaptive Frequency Control

The decarbonization of many heavy power-consuming industries is dependent on the integration of renewable energy sources and energy storage systems in isolated autonomous power systems. The optimal energy management in such schemes becomes harder due to the increased complexity and stability requirements, the rapidly varying operating conditions and uncertainty of renewable sources, the conflicting objectives across different timescales, the limited amount of reliable power sources and energy storage. The state of charge management when energy storage is used for multiple services, such as optimal scheduling and frequency support, is one of the most notorious problems in this context. To address this issue, an optimal energy management system is proposed in this paper. It co-optimizes the primary frequency control layer and the dispatch schedule of conventional generators and energy storage by taking advantage of an algorithm that provides adaptive active power demand uncertainty quantification, theoretical guarantees for frequency stability, and bounds for the reserves for frequency support assigned to the energy storage system. A convex reformulation is derived enabling the efficient solution of the involved optimization problem, being a test case of an isolated offshore oil and gas platform presented for validation.


I. INTRODUCTION
T AKING optimal decisions to reduce i) fuel consumption, ii) greenhouse gases emissions, iii) equipment degradation and iv) system insecurity in power systems with a high share of renewable energy sources (RESs) is an intricate task.This endeavor requires the solution of problems such as the unit commitment (UC) [1], economic dispatch (ED) [2], and allocation of reserves for frequency control [3], resulting in complex non-convex optimization formulations.Those have typically conflicting objectives, continuous and binary decision variables, and high uncertainty from particular nodes, such as loads and RESs.Reliability can, for instance, be increased where more dispatchable units, such as gas turbines (GTs), are kept online for longer periods, as they respond well to fast variations of loads and RESs.This action may, however, increase fuel consumption, emissions, and equipment degradation, affecting negatively operational expenditure (OPEX).
When compared to traditional bulk power systems (BPSs), this type of optimization can have different characteristics in autonomous power systems (APSs), such as isolated industrial plants, oil and gas (O&G) platforms, ships, islands, and community microgrids.Due to the size and complexity of the problem, several challenges exist in BPSs to implement real-time, advanced algorithms for optimal dispatch and realtime allocation of frequency reserves [4].Such problems can however be addressed in many APSs due to the limited number of dispatchable power sources, which typically provide simultaneously several ancillary services.The methods for allocation of frequency reserves and tuning of the system damping may also differ considerably between BPS and APS.The system damping is, in the former, adapted using binary decision variables which switch on and off units having a fixed active power-frequency droop, as usually the individual contribution of each unit to the total damping is small.In the latter, the droop of individual units may represent a large portion of the system damping, being therefore necessary to readjust them in real-time using integer variables to obtain optimal results.
Optimization objectives and decisions in APSs may also be coupled and affected by constraints in different time scales, being energy storage systems (ESSs) remarkable examples of equipment that enhances such dependencies.ESSs can, for instance, provide increased flexibility towards the optimal scheduling of fossil-based energy sources and avoid prolonged operation in partial loads, where emissions are much higher.They can also be assigned as spinning reserves for frequency control, which would require less GTs on for the same system security requirements and increase environmental gains even further.Where these two grid services are provided simultaneously, the scheduled trajectory of the ESS' state of energy (SoE) may be disturbed, affecting the optimality or even the feasibility of the original schedule.To decide in BPSs the effect that the provision of frequency reserves has on the ESS' optimal schedule, scenarios such as the worstcase active power disturbance (APD) or the N-1 criterion are many times used [3], [5].These criteria can however arXiv:2208.08953v1[eess.SY] 18 Aug 2022 be over-conservative in APSs [6] and the capital expenditure necessary for a fully-fledged ESS may not justify its benefits, being some probability of load shedding and generation curtailment acceptable many times [7].When it comes to energy management in APSs with high penetration of RES, there is therefore a need and potential for better assessment of frequency stability requirements as well as coordination and scheduling of reserves.

A. Literature review
The integration of frequency stability constraints in the scheduling phase of BPSs with high penetration of RES has been an active area of research in recent years and is well described in the literature [8].The use of linear constraints had been a typical approach to increase the fidelity of the frequency response in the UC and ED problems, which are obtained by linearizations and/or analytical solutions of the swing equation model [9]- [16].These procedures provide good approximations in BPSs, where the largest APD is a small fraction of the installed capacity and transient and steady-state frequency deviations are usually required to be below 2% of its rated value [5].Frequency deviations may, on the other hand, be much higher in APS due to low inertia and limited amount of frequency reserves when compared to the worst-case APDs, causing the effects of non-linear dynamics to be sizeable during large disturbances [17].Where these nonlinear effects are considered, bilinear terms are introduced in the frequency stability constraints, those requiring specialized reformulation-linearization techniques [18]- [20].This increase of the optimization model complexity can play an important role in APS, while it may bring only negligible frequency stability improvements in BPS.
The challenge of leveraging optimal system operation and security in APSs was also addressed recently, being alternatives presented to the deterministic evaluation of the worstcase APD or the N-1 criterion.References [13], [21]- [24], for instance, considered the effect of anticipated net load variations and applied dynamic constraints for sizing frequency containment reserves and required inertia.These works did not evaluate however the impact that re-adjusting the droop of online units, instead of switching on and off units with fixed droops, would have on frequency stability and OPEX.References [9], [15], [20], [25] tried moreover to tackle the problem by applying simplified and static uncertainty models, such as non adaptive and arbitrary uncertainty intervals, distributions, and scenario selection, which have the limitation of not providing probabilistic guarantees.

B. Paper contributions
The main contribution of this paper is an algorithm for the energy management system (EMS) of a generalized APS equipped with an ESS that is capable of simultaneously achieving optimal scheduling and securing system operation under dynamic uncertainty considerations and bounded impact on its optimal scheduling.The proposed algorithm introduces realistic frequency constraints with general applicability, those not being limited to small scale power variations, combined with a probabilistic security analysis based on a novel adaptive power variation uncertainty quantification scheme.This new method allows an EMS to allocate time-varying optimal frequency control reserves with bounded divergence from its optimal SoE.A simple strong mixed integer linear programming (MILP) reformulation is derived to efficiently implement the proposed algorithm, whose effectiveness is verified by simulations using the case study of a wind-powered offshore O&G platform.
The rest of the paper is organized as follows.The proposed EMS algorithm and its numerical implementation are presented in section II, while simulations in Matlab/Simulink are used to validate and discuss the various features of the proposed algorithm using the case study of a wind-powered offshore O&G platform in section III.The main conclusions are finally presented in section IV.

II. METHOD
The concept and the methodology proposed in this paper is presented in detail in this section.It assumes the existence of a centralized EMS in an APS that is capable not only to dispatch a set of g generators and a single ESS b, but also to decide for the primary frequency controllers the proportional gain (droop) D g of each generator, the proportional gain (droop) D b and derivative gain (virtual inertia) of the ESS, as summarized in figs. 1 and 2. To take these decisions, an adaptive and probabilistic step-like algorithm quantify net APD at discrete time intervals T is applied, as shown in fig. 2

A. Frequency Control and Reserves Allocation
The non-linear dynamics of the frequency deviation for the center of inertia (COI) in an ac power system are given in Fig. 2. Visualization of the proposed adaptive and probabilistic step-like net active power disturbance quantification at a time instant in the continuous range between the discrete points in time where decisions are taken.
It is shown in [17] that if the total system damping D is selected as in eq. ( 2), then X will be bounded post-disturbance by r ss for an instantaneous net active power imbalance P nl b , where r tr is a pre-defined constant bound for the nadir/zenith of X during the transient period.Note that the step-response of frequency for active power disturbances is, in most ac power systems, underdamped due to delays of actuators and their non-linearities [26].Those effects are not modelled in eq. ( 1) because they influence only the frequency nadir or zenith after the disturbance P nl b occurs, being the post-disturbance steadystate value of X unaffected [17].The variables r ss and r tr in eq. ( 2) represent, in other words, the allowed steady-state and transient frequency deviations in per unit of X , being usually defined in grid codes or by the system operator.
It is typical to assume r tr > r ss to provide a safety margin which avoids triggering protection schemes such as under frequency load shedding or over frequency generation curtailment during large disturbances.When r ss = r tr , there exists no safety margin between steady-state and transient values.It shall also be noted that eq. ( 2) is valid for any value of X , and not only for small deviations around the operating point.
The minimum required inertia M is retrieved from eq. ( 1), assuming that a) the disturbance P nl b occurs at t = t 0 , b) the system was previously in balance (X | t=t 0 = 1), and c) the maximum rate of change of frequency (RoCoF) is Ẋ | t=t 0 = γ, which results in eq. ( 3).
Ensuring that enough power reserves are available in an APS is a necessary condition to satisfy the constraints for frequency control imposed by eqs.( 2) and ( 3).The values of D and M depend on the number of online units and their individual parameters, being expressed by eqs.( 4) and ( 5), where D g and M g denote the droop setting and the inertia of each online generator, and D b and M b are the virtual damping and inertia emulated by the ESS, respectively.The EMS must guarantee that each individual equipment has enough available power capacity to provide the frequency reserves assigned to it, what leads to the constraints in eqs.( 6) and (7).
B. Optimal and Bounded Energy Management under Uncertainty 1) Allocation of frequency reserves with bounded energy storage: The EMS must also ensure that enough energy capacity is available in the ESS.This is required not to violate the storage limits when providing frequency control in the current scheduling step t 0 while charge or discharge events had already been scheduled for the next step t 0 + T , where T is the time step of the EMS.To deal with this situation, the constraint in eq. ( 8) is proposed, where ∆E b (t) is a deviation around the optimal scheduled SoE E b (t), and λ is a userdefined (absolute) percentage value.The proposed constraint allows, in other words, an allocation of frequency reserves to the ESS that is proportional to its current energy level.
For energy calculations of frequency containment reserves in EMS algorithms, one can assume that X ≤ r ss , ∀ t ∈ [t 0 , t 0 + T ] when the system is frequency stable [17], this simplification resulting in eq. ( 9).Remark that, on one hand, r ss ≤ X ≤ r tr during the arrest and rebound periods of frequency containment.On the other hand, those periods combined last only a couple of seconds, while T may vary from 5 minutes to 15 minutes in typical EMSs.The integration error would therefore be minimal when adopting this simplification.The latter however allows ∆E b ≤ ∆E b (t) to be enforced in the optimization problem, as described later in section II-C, and to satisfy eq. ( 8), which avoids excessive deviations of the SoE from its optimal dispatch.
2) Adaptive Uncertainty Quantification: Load uncertainty is, in the proposed method, evaluated not only in rapidly varying loads but also in the intermittency of non-dispatchable power sources, being the net load variation defined as the combination of two random variables P p , where p = { , w} corresponds to load and non-dispatchable sources respectively).The time-varying cumulative distribution function (CDF) for each random variable and predication lead time presented in eq. ( 10) can be inferred when using ensembles of quantile values following the method in [27], where x r (t) is the regressor composed only by lagged values of the corresponding variable.The net load as a function of time is represented as ξ t+k (t) = P t+k (t) − P w t+k (t), where the two random variables are combined, resulting in the adaptive probabilistic forecasting capability which will be exploited in the construction of the security constrained energy management algorithm described in section II-C.

C. Probabilistically Constrained EMS
This section presents how the proposed EMS algorithm integrates into a single optimization problem including the various objectives and multiple time-scale requirements described in sections II-A and II-B.The system dynamics are first expressed in a hybrid state space system to achieve this integration as shown in eq.(11), where x + is the state at the next discrete time, u the corresponding control inputs and A, B the corresponding system matrices.Note that x SoC = E b / Ēb ∈ [x SoC , xSoC ] and x gt 1:Ng = {x gt g } 1:Ng ∈ {0, 1}.
The optimization problem presented in eq. ( 12) is then formulated and solved in each discrete time step t, where δ = {δ i k } 1:K is a random multi-sample, ξ = {ξ k } 0:K−1 is the deterministic mean net load forecast including the measured value ξ 0 at t = t 0 , and P b nl = {P b nl (δ k+1 |ξ k )} 0:K−1 is the net load perturbation for the whole prediction horizon as a function of the sample and the mean forecast.
In the cost function F, J (x + ) represents the cost of having generators online, J (z) penalizes the deviation from optimal operating conditions of online generators and the ESS, J (u) captures the start up cost of generators, and finally w T v accounts for costs of activating reserves for frequency control.The decision variables v = [D g , D b , M b ] T used for frequency reserves are weighted by w so that different contributions can be assigned to online generators and the ESS to ensure the system frequency stability.
The term C(•) encapsulates all the constraints related to the operation of the energy system, as described in [27].The last expression in eq. ( 12) is a chance constraint that, given the estimated distributions from eq. ( 10), relates potential instantaneous perturbations for a net load value ξ k and the next step sampled net load δ i k+1 .Equations ( 2) and ( 3) are, in other words, risk-constrained by , which act as a mechanism to leverage cases with rather pessimistic prediction intervals resulting from poor uncertainty range estimation.
The chance constraint complicates, however, the optimization problem.Note that P is a stochastic mixed integer nonlinear problem and therefore non-convex, so a standard scenario approach cannot be applied.Following [28], [29], a probabilistic set composed of a finite number of samples δ i ∈ ∆ N is computed as in eq. ( 13), being β and e user defined parameters that tighten the non-violation of the original chance constraint.Note that the number of random variables is 2|K| because two sources of uncertainty exist (load and renewable power injection) for the whole prediction horizon.The chance constraint can however be replaced by the deterministic set of linear inequalities presented in eq. ( 14), which transforms the stochastic problem P into a robust one, where the set ∆ N encapsulates the same risk-volume as the original chance constraint.When solving the corresponding robust program, solutions are not only feasible for the initial P but also satisfy the same probabilistic guarantees [30].
It is not feasible, on one hand, to implement the constraints expressed in eqs.( 2), ( 4) and ( 7) directly into a MILP formulation because of the bi-linear terms x gt g D g .It is possible, on the other hand, to take advantage of the following observation: with N g generators, there are exactly 2 Ng different possible system configurations depending on the values of x gt 1:Ng ∈ {0, 1}.The status variable x gt g,t takes a particular value (0 or 1) for each generator and configuration, being therefore feasible to represent it by binary strings of length N g .The feasible system configurations are then enumerated by taking all the possible permutations, which can be gathered in a data table A cf [2 Ng × N g ].The variables x gt g,k , ∀ k ∈ K can, in this way, be treated as constants given a selected system configuration, being indicator variables b j introduced to identify the configuration selected and associate it with the status of the generators by the set of constraints in eqs.(15) and (16).Equation ( 4) can, as consequence, be linearized using A cf as in eq. ( 17).
Equations ( 2) and (3) can also be reformulated using the sampled net APD terms from eq. ( 14).This is presented in eq. ( 18), where M B is a big-M value, b k , D k correspond to each lead time of the prediction horizon (k ∈ K) assuming that M g are constant values.The different possible configurations are restricted by a type 1 special ordered set constraint as in eq. ( 19).
Equation ( 7) is similarly linearized using A cf and M B as in eqs.( 20) and ( 21) where S b is the system's base power.
The ESS reserves from eq. ( 6) are reformulated in terms of discharging and charging power as in eq.(22), where s k is the discharging binary indicator variable.
The bound derived in eq. ( 9) can be easily integrated to P through the linear constraints in eq. ( 23), where ν = rssS b 3600 Ēb .
Equation ( 12) may therefore be expressed as a deterministic MILP robust program and solved efficiently by the EMS at each time iteration, as F(•) and C(•) can be expressed as linear combinations of the optimization variables when including the reformulations proposed in eqs.( 14) to (23).

III. SIMULATIONS
The various components of the proposed methodology are demonstrated in this section, namely the adaptive uncertainty quantification and the integration of frequency stability constraints with bounded use of the ESS stored energy in the optimal scheduling.To achieve this goal, the case study of a wind-powered offshore O&G platform is presented, where a reference operation time period is considered involving both regions of smooth, low magnitude and sudden, large net load variations.The effectiveness of the proposed methodology is then demonstrated with time domain simulations on the scheduling time scale.The optimization problem is solved with Gurobi 9.1.0in a 28 physical core multi-node cluster with Intel(R) Xeon(R) CPU E5-2690 v4 @ 2.60 Hz and 25 GB RAM.The solution time of eq. ( 12) using the proposed formulation is well below 15 minutes, which is assumed here as the minimum threshold for real-time power system scheduling.

A. Capabilities of adaptive uncertainty quantification
A representative example of how the adaptive uncertainty quantification framework works is illustrated in fig.3, where a case of a sudden APD was selected since these are the most interesting from a power imbalance perspective.The actual load values are presented with the solid black line and the blue cross indicates the time instant at which a probabilistic forecast is issued.From that moment and for a given prediction horizon (k ∈ K), the expected values of the prediction are plotted in red dashed line along with the prediction intervals around them in green at various quantile levels (α).Samples generated by the estimated CDFs are also plotted as purple dots.
The proposed adaptive uncertainty quantification algorithm generates samples that better describe the size of an APD.Observe that, as the blue cross moves forward in time, the prediction intervals adapt to capture the irregular event of the sudden load increase.Note also that, at the initial time (fig.3a), the prediction intervals are narrow and all the sampled values fall close and around the actual load values.From the next moment and later (figs.3b to 3g), however, the uncertainty increases and the prediction intervals expand to capture the possibility of an irregular sudden load increase.The deterministic forecast, in contrast, fails to capture the variation adequately, because it is dominated by the inertia of past values, a common drawback of auto-regressive models.

B. Effect of optimal and bounded frequency support from ESS
To demonstrate the effect of optimally controlling the ESS to provide frequency support for an isolated APS, a step APD is considered at time t 0 = 2s representing a 0.4 pu load increase from a sudden motor startup, where a single GT is on in the platform of the case study.Simulations were run in Matlab/Simulink 2022a using the model given in eq. ( 1), which do not include delay of actuators.This simplification facilitate the interpretation of results, not implying however any loss of generality.
Fig. 3. Demonstration of adaptive uncertainty quantification by using autoregressive probabilistic forecasting for the load time-series.A case for a sudden step-like variation is presented for consecutive lead times (figs.3a to 3g).By updating the estimated prediction intervals it is possible to capture the sudden load variation and draw samples (purple dots) that span an appropriate range of values.
Results are illustrated in fig.4, where the frequency deviation (solid lines) and the RoCoF (dashed lines) for two cases are represented.In the first case (grey lines), the single online GT is the only source of primary frequency control, while, in the second case (black lines), the ESS supports in this task the GT, which has the same droop setting as in the previous case.Observe that the system's response in the first case violates not only the steady-state frequency bound r ss but also the maximum allowable RoCoF (γ), whereas in the second case both limits are respected.To respect the defined bound in the first case, the GT droop setting must be increased, leading to a larger deviation from the optimal GT operating In contrast with the case of non participation of the ESS to frequency regulation (solid and dashed grey lines), the optimally designed virtual inertia and damping results in a frequency response (solid black) and RoCoF (dashed black) that are bounded by their defined limit values.point, decreased efficiency, increased fuel consumption and emissions.The same droop setting for the GT, on the other hand, can be maintained and the r ss threshold observed where an optimal participation of the ESS in the primary frequency control (PFC) has been decided beforehand by the algorithm using an adaptive droop setting.Note also the compliance with the RoCoF limit where the ESS provides virtual inertia.To respect this limit without the ESS support would require an additional GT on, affecting significantly the overall efficiency of the system.The proposed algorithm, in other words, employs an adaptive droop and virtual inertia scheme, enabling not only the optimal scheduling of the power system by avoiding an additional GT on, but also guaranteeing that the services provided by the ESS to the grid will not cause it to deviate significantly from its optimal schedule, as the use of stored energy is bounded by eq. ( 23).As a matter of fact, the calculated energy for frequency support provision by the ESS to the described disturbance for T (15 minutes) is 0.9179 M W h while the bound for λ = 3% is 0.9369 M W h.

C. Comparative analysis and effect of bounds
To further demonstrate the capabilities of the proposed EMS, a reference operation period is simulated considering: I) the default EMS that does not integrate any bounds; II) the version that includes frequency variation bounds but not the energy bounds on the ESS; and III) the proposed version which includes all of them.
1) System inertia and damping evolution: The results for all EMS versions are aggregated in fig.5, being the system inertia and damping evolution depicted in fig.5a and fig.5b along with the net load signal P nl b (t) (solid black line), during a reference period of 8 hours.No frequency and RoCoF containment reserves are considered in case I, the droop settings not being optimization variables and the ESS not contributing with virtual inertia.The system inertia M(t) (red line in fig.5a) is therefore just the result of the online GTs based on the optimal techno-economic scheduling, not having any virtual component M b (t).Following the same color notation in both fig.5a and fig.5b, the evolution of the system inertia and damping is observed, being the contribution of the ESS in case II denoted by the superscript f (frequency bound) and in case III by f, e (frequency and energy bound).
Two main patterns are evident when observing fig.5a and fig.5b.Firstly, more inertia and damping are assigned for both cases II and III close to the time instants of the sudden net load variations (around 14:00 and 18:00 correspondingly).Secondly, the same inertia and damping are noticed during the relatively constant net load conditions (between 14:00 and 18:00).This effect demonstrates the adaptive capabilities of the proposed EMS algorithm to assign more or less inertia and damping in correspondence with the anticipated net load variations.The larger the variations expected, the more secure the system will be by properly deciding its power and energy reserves.Whether the ESS is not assigned to frequency control, securing the system for possible net load variations would require additional GTs to be on.This is confirmed in fig.5a during both sudden net load variations, where M(t) ≤ M f (t) (green line) and M(t) ≤ M f,e (t) (magenta line).Note also that around those instants, M f b (t) ≥ 0 (cyan line) and M f,e b (t) ≥ 0 (blue line), meaning that the additional GT can be avoided by properly assigning virtual inertia to the ESS.The comparison of the different methods (I, II, III) is also quantified through the cumulative results of specified key performance indicators in table I. Including bounds on the EMS results in slightly higher fuel consumption and operating costs, which is primarily attributed to higher cumulative operational hours of the GTs.This is in agreement with the results presented in fig.5, as securing the system against possible disturbances may be associated with increased inertia requirements.It is also noteworthy in table I that the incremental cost for including bounds on the use of ESS stored energy (III) compared to case II is negligible, meaning that the optimal GT scheduling is almost not affected by the inclusion of energy bounds in the ESS.The main difference is that the resulting optimal GT trajectory for II is associated with one less GT start up compared to III and slightly lower ESS cycling, reflected on the lower degradation, which is however almost equal in cases I and III.
2) State of charge evolution and energy bounds effect: An additional comparison demonstrating the additional benefits coming from method III over method II is illustrated in  Notice that both methods respect the upper and lower state of charge limits, as originally designed in the default method I.Both trajectories seem to follow similar patterns, i.e. initial discharge until k = 9, smooth re-bouncing and discharge until k = 25, and charging until the end.There are however small but important differences, those being clearly depicted (a) Without bounds on the energy deviation.IV.CONCLUSIONS Achieving optimal energy management in isolated power systems with energy storage and sudden load variations cannot be decoupled from ensuring their secure and robust operation, especially under the presence of intermittent renewable power sources.Even though decisions related to techno-economical operation are conventionally taken in the discrete time, those will inevitably affect the system's stability in its continuous operation and vice versa.To address this problem, this article proposed an energy management algorithm capable of integrating both higher time scale economic objectives and lower time scale stability constraints under adaptive uncertainty considerations.Additional constraints regarding the optimal use of the energy storage for providing flexibility and frequency support with bounded interaction between both services were proposed.A MILP formulation was derived and validated through time-domain simulations for an isolated offshore O&G platform integrating wind power.The results indicated that, under the proposed adaptive uncertainty framework, optimal decisions with dynamic frequency stability guarantees could be achieved and secured the system under an adaptive assessment of possible active power perturbations.This also reduced the conservatism from setting fixed damping and inertia requirements based on the expected worst-case and allowed better scheduling and operation of the GTs for longer periods.At the same time, the optimal sharing of primary frequency control contribution from conventional generators and the energy storage was found, while ensuring a tolerable impact on the storage optimal state of charge schedule and a negligible impact on the rest energy management objectives.

Fig. 1 .
Fig. 1.Hierarchical control system schematic where the upper layer optimal discrete time control is integrated with the lower time scale continuous adaptive primary frequency control.

Fig. 4 .
Fig.4.Effect of the optimally calculated participation of the ESS in providing frequency support for a step load change when only 1 GT is on.In contrast with the case of non participation of the ESS to frequency regulation (solid and dashed grey lines), the optimally designed virtual inertia and damping results in a frequency response (solid black) and RoCoF (dashed black) that are bounded by their defined limit values.
Similar observations are drawn from fig. 5b where D f b (t) ≥ 0 (cyan line) and D f,e b (t) ≥ 0 (blue line) around the sudden variation instants and D f b (t) = 0 and D f,e b (t) = 0 for the rest period.Remark that the effect of the energy bounds on the ESS are stronger in the damping terms than in the inertia.Figure 5a shows that the signals M f b (t) and M f,e b (t) are almost identical for the whole time period, whereas fig.5b depicts that sup D f,e b (t) ≤ sup D f b (t).This means that the peak contribution of the ESS for case III is smaller than the one for case II.It is indeed possible to observe that D f,e b (t) is almost always less than D f b (t), further demonstrating that version III is more cautious not to overuse the PFC of the ESS.

Fig. 5 .
Fig.5.Trajectories of the optimally designed system inertia (top) and damping (bottom) and optimal split between primary control provision by the GTs and the ESS for different bounds considerations.The simulated results are plotted against the net load signal (solid black) for a case where sudden step-like variations occur.
With bounds on the energy deviation.

Fig. 6 .
Fig.6.State of charge trajectories from the optimal scheduling and primary control design of the ESS, during the simulation period, for the cases of not including bounds on the resulting energy deviation from the participation in PFC (top) and the one when including the bounds (bottom). .

TABLE I PERFORMANCE
COMPARISON FOR THE WHOLE SIMULATION PERIOD.