A hierarchical frequency stability control strategy for distributed energy resources based on ADMM

With the large-scale distributed energy resources (DERs) interfaced by power electronic converters connected to the power grid, the traditional synchronous generation is gradually replaced, and the inertia and primary frequency regulation capability of the power system is weakened, which seriously threatens the frequency security of the power grid. In view of the above problems, a hierarchical control strategy by coordinating the distributed wind turbine (WT) and photovoltaic (PV) is presented to improve the transient frequency (TF) and steady-state frequency (SSF). For the sys-tem control level, the optimal frequency regulation coefficients of multiple WTs and PVs clustering systems are determined by the optimal control principle with minimum frequency regulation energy, while considering the constraints of TF and SSF. For the WT/PV control level, the frequency regulation coefficients of WTs and PVs are solved by the model predictive controller (MPC) and alternating direction method of multipliers (ADMM), with the minimum control cost of WT and PV. It achieves the optimal allocation of frequency regulation power and improves the computing efficiency. Finally, the validity of the proposed method is tested in the modified IEEE 10-machine 39-bus system.


Introduction
To construct a new power system with new energy as the main body in China, the traditional synchronous generators are gradually replaced by the power electronics devices, which contain photovoltaic power, wind power, energy storage system (ESS), and high voltage DC (HVDC) transmission, resulting in the insufficient inertia level and the poor primary frequency regulation capacity of the power system [1,2].Thus, the frequency stability of the power system is facing severe challenges.Especially, in recent years, some blackouts related to the frequency security have happened occasionally.For example, the tripping of wind farms in Great Britain [3] and South Australia [4] that happened in August 2019 and September 2016 caused the frequency drop and eventually led to a massive blackout.The frequency stability of new power systems has become a research hotspot.
Each country puts forward the related requirements for the TF and SSF.For example, the national grid of Great Britain mandated the limit that the maximum frequency deviation and steady-state frequency deviation should be maintained within 0.8 Hz and 0.5 Hz for events considered to be abnormal [5].In China, the technical regulations for under-frequency load shedding (UFLS) imposed the limit that the SSF is not less than 49.5 Hz after a large disturbance [6].In a new power system, the system frequency drops rapidly after an active power disturbance.When the frequency is less than the preset threshold, the UFLS is activated, and it will lead to shed an amount of load.Furthermore, if the frequency regulation capability of the system is deficient, UFLS alone will not be enough to recover the SSF to the stability range, and it may eventually lead to the separation of the system.
The distributed energy resources, which include the distributed WT/PV, energy storage system, and flexible load, are mostly located in the distribution network.Due to the flexible control and large quantity, DERs can be used as frequency regulation resources to participate in the system frequency control and ensure frequency safety.
Currently, there are many studies on the optimal operation of distributed energy resources.In [7], based on a multi-agent system, a coordinated optimal operation scheme for the distribution network is proposed and solved by the consistency algorithm.In [8,9], the ADMM algorithm is used to solve the distributed active and reactive power of new energy resources, and it achieves the optimal distribution of power.Moreover, the distributed energy resources is also used to control the system frequency.In [10], the coordinated frequency control strategy for DERs in a virtual power plant (VPP) is solved by the distributed subgradient-projection method while minimizing the generation cost.In [11], an optimization-based method is presented to determine the injection power of DERs, which the power is in proportion to their power rating.In [12], a coordinated frequency control scheme of VPP that enables a VPP concept by adaptively organizing the output power of an AS-PSH plant and energy storage systems according to their operating conditions is presented to improve the frequency characteristic.In [13], a two-layer control structure is designed to enhance the performance of frequency regulation and power optimization for microgrid.The primary control adopts the power-increment factor droop, and the secondary control uses a consensus protocol.It achieves the goals of frequency synchronization, frequency no-difference, and power optimization in a distributed manner.In [14], how to involve DERs in contingency frequency support through VPPs is investigated.An equivalent aggregation model for DERs is developed, and the frequency performance-to-cost map is built via the cost minimization optimization problem.In [15], a two-level coupling-based frequency control strategy for microgrids is developed.At the lower level, an adaptive dynamic compensation algorithm is designed to improve frequency characteristics.At the upper level, an adaptive distributed consensus algorithm is developed to address frequency restoration and active power sharing.In [16], based on the consensus algorithm, the distributed energy resources is used to participate in the primary frequency regulation of the power grid, and the frequency response coefficient is calculated adaptively to respond to the frequency changes rapidly when the fault occurs.In [17], based on the applied deep reinforcement learning algorithm, a multi-level control scheme of battery energy storage system (BESS) and electric vehicles in networked microgrids is proposed to participate in the energy management of active power and frequency control.The results show that the accuracy and efficiency of calculation are obviously improved.In [18], based on an intelligent probabilistic wavelet petri neuro-fuzzy inference algorithm, the power management and control platform is proposed to control the voltage/frequency metric of microgrids, considering renewable energy sources (RESs), BESS, and various uncertainties.In [19], considering the optimal participation of RESs, ESSs and the integrated demand response (IDR) programs execution, a hierarchical energy management approach is solved by the multi-agent deep reinforcement learning algorithm with the objective of minimizing the operating costs, environmental pollution costs, risk costs, and destructive effects of cyber-attacks.In [20], considering communication delay and input saturation, a distributed fractional-order predefined-time sliding mode controller is proposed to adjust the energy storage device, which has a fast response ability and strong regulation capability, to increase the transient stability of the power system.
The new power system with new energy resources as the main body in China requires that the new energy resources have the ability of frequency regulation.For the centralized new energy resources, the traditional centralized control is used to participate in frequency stability control.The control method can be categorized into two aspects.One approach is the response-driven frequency control, which mainly includes virtual inertia control and droop control.Another approach is the event-driven emergency frequency control, which includes emergency demand response and emergency power support.However, for distributed energy resources, the traditional centralized control method is not suitable.The distributed cooperative control method is widely used in power systems.Furthermore, the alternating direction method of multipliers, consensus algorithm, and distributed gradient descent algorithm provides effective methods for distributed energy resources to participate in frequency control, and these methods have been studied in frequency control of virtual power plants, microgrids, and distribution networks.While the above methods have a positive contribution to enhancing system frequency stability, there also remain some issues: 1) Control form of system frequency.Due to the large number of distributed energy resources, the potential capacity of frequency regulation is huge.However, the current most frequency control form, in which the output power responds to the change of the system frequency, may not be the best control strategy, since the primary frequency response speed of traditional synchronous units is relatively slow.Furthermore, because of the limited frequency regulation energy of these resources, the energy may not be fully utilized for the traditional frequency control form, and it possibly causes unnecessary control costs.Finally, most literature only focuses on the transient frequency control, but they do not consider the situation that the steadystate frequency does not meet the requirement due to insufficient reserve capacity.Therefore, considering the characteristics of flexible control and quick response of distributed energy resources, how to simultaneously improve the TF and SSF stability by using the limited frequency regulation energy after a large disturbance still remains a challenging problem.
2) Coordination among distributed energy resources.With the growth of the number and scale of the distributed energy resources, the distributed WTs and PVs will operate under different conditions.If the WTs and PVs participate in the frequency control after detecting frequency disturbances without any coordination, it is possible that they release more or less frequency regulation power to the power system.It will lead to a severe frequency excursion.Moreover, although the coordinated control of distributed energy resources has been achieved in some literature, the frequency characteristics of distributed energy resources are not fully considered, e.g., the detailed frequency response model of WT and PV.It does not obtain an optimal result.Finally, most of the frequency control schemes are the traditional centralized control.If this method is applied to distributed energy resources, the amount of computation will increase significantly.Therefore, considering the frequency response characteristics of distributed energy resources, how to coordinate the distributed energy resources and reasonably allocate the frequency regulation power among distributed energy resources is also a problem that needs to be solved.
So, a systemic frequency stability control strategy focusing on the coordination of distributed energy resources is required.
Based on the above issues, a hierarchical control scheme for the distributed energy resources is proposed to simultaneously improve the TF and SSF.The major contributions of this paper are as follows: 1) A modified multi-machine frequency response (MMFR) model, which is comprised of traditional thermal generators, traditional hydro generators, and distributed energy resources interfaced by power electronic converters is formulated to analyze the system frequency dynamic process.2) A hierarchical control framework is formulated to coordinate the frequency regulation power of multiple WTs and PVs clustering systems and WTs/PVs with different objectives after a large disturbance.3) For the system control level, considering the limits of the TF and SSF, an optimization model is constructed to determine the total optimal frequency regulation coefficients of multiple WTs and PVs clustering systems with minimum frequency regulation energy.4) For the WT/PV power distribution level, an optimization problem is constructed to determine the frequency regulation coefficients of WTs and PVs with minimum frequency regulation cost, considering the frequency response model of WT and PV.Furthermore, in order to improve the computational efficiency, the optimization problem is solved by ADMM.

Methods
The framework of the proposed hierarchical control strategy is illustrated in Fig. 1.
For the system control level, when a power disturbance is detected, the central controller receives the system operating information through WAMS, which includes the operating status of all traditional synchronous generators, the operating status of WTs and PVs, and corresponding parameters.According to the information, the disturbance power is estimated, and the frequency nadir (FN) and SSF are predicted [21].If the FN and SSF are less than the stability threshold value, the control scheme is activated.The optimal frequency regulation coefficients of WTs and PVs clustering systems are determined by the optimal control principle, with the aim of reducing total frequency regulation energy, considering the constraints of the TF and SSF.The optimal frequency regulation coefficient of the kth WT and PV clustering system is denoted as K WPk .
For the WT/PV control level, after receiving the K WPk from the system level, the frequency regulation coefficients of WT and PV are determined by MPC, and then the coefficients are sent to each WT and PV.The optimization target is to minimize the frequency regulation cost of the WT and PV.And the equality constraint is also satisfied for the K WPk .Besides, the optimization problem of the WT/PV level is solved by the ADMM algorithm, and it can reduce the computation burden of the WT and PV clustering system.The frequency regulation coefficients of the ith WT and jth PV in the kth clustering system are denoted as K wki and K pkj .

System level-frequency control
For the system level, the control objective is to ensure that the FN and SSF are within the frequency stability range after a large disturbance.The system level controller receives the operating information from the power system through WAMS and calculates the disturbance power.Meanwhile, the FN and SSF are predicted.If the FN and SSF are less than the stability threshold value, the system level control is activated.The optimal frequency regulation coefficients of WT and PV clustering systems are solved by the optimal control principle with the minimum frequency regulation energy.

Frequency response model with distributed energy resources
A modified multi-machine frequency response (MMFR) model is used to analyze and control the system frequency after a disturbance.This model includes the traditional thermal generator, traditional hydro generator, and distributed energy resources based on power electronic interface (WT and PV).The MMFR model with distributed energy resources is shown in Fig. 2.
1) Frequency characteristic of thermal generator.
The transfer function of the thermal governor is described as [22]: where T G is the response time of the thermal generator governor and R g is the frequency regulation coefficient.
The transfer function of the thermal turbine is described as [22]: (1) where F H is the proportion coefficient of the high-pressure turbine of the thermal generator; T C is the response time of the steam chest; T RH is the response time constant of reheat; and ΔP mg is the output power increment of the thermal generator.
2) Frequency characteristic of hydro generator.The transfer function of the hydro governor can be described as [22]: where T h is the response time of governor; R h is the frequency regulation coefficient.
The transfer function of the hydro turbine can be described as [22]: where R c is the comprehensive frequency regulation coefficient; T R is the reset time constant of the hydro generator governor; T w is the starting time constant of the hydro turbine; and ΔP mh is the output power increment of the hydro generator.
3) Frequency characteristic of distributed energy resources.In this paper, the deloaded operation is applied for WT and PV to ensure enough power reserve, which can provide the primary frequency response.Therefore, the relationship between the system frequency and WT output power can be expressed as: where K WPk is the frequency regulation coefficient of the kth WT and PV clustering system; ΔP WPk is the output power increment of the k th WT and PV clustering system.
4) Frequency characteristic of aggregated load.The characteristics of the aggregated load are described as [22]: where K L is the active power-frequency regulation coefficient of aggregated load.
The equivalent damping coefficient D can be expressed as: where D M is the damping coefficient of the equivalent generator.5) MMFR model.The MMFR model with distributed energy resources can be obtained as: where H is the system inertia constant; D is the equivalent damping coefficient; Δf is the system frequency deviation; and ΔP L is the disturbance power.6) System state-space equation.According to Eqs. ( 1)-( 10), the continuous state-space equation can be expressed as: (4) with.
x 1 = [ P vg , P cg , P mg , P vh , P wh , P mh , f ] T , u 1 = P WPk , e 1 = P L , y 1 = f , where x 1 is the state variable; u 1 is the control variable; e 1 is the disturbance variable; and y 1 is the output variable.A 1 , B 1 , E 1 , and C 1 are the corresponding matrices.
The SSF deviation can be calculated by the final value theorem, and it is obtained as: where

Objective function
Since there are control energy constraints on WTs and PVs participating in frequency regulation, it is hoped that the TF and SSF are controlled above the preset threshold while minimizing the control energy.Thus, the control problem is transformed into an optimization problem.The objective function can be written as: where t n is the ending time of optimization; γ k is the weight coefficient; and N is the number of the clustering system.
where ΔP WPk max is the maximum frequency regulation power provided by the kth clustering system.

Constraints
1) Power constraint of WT and PV in clustering system. ( In this paper, WT and PV operate with the deloaded mode and ensure a d% power margin.Therefore, the maximum frequency regulation power of the k th WT and PV clustering system can be obtained: where ΔP WTk max and ΔP PVk max are the maximum frequency regulation power of WT PV in the kth clustering system, respectively.
The kth WT and PV clustering system should follow the power limit: 2) TF constraint.
In the process of frequency change, the TF should be satisfied as: where f 0 is the standard frequency; f th is the TF threshold value; and ε is the stability margin.
The SSF should be satisfied as: where f ss th is the preset SSF threshold value and δ is the stability margin.

Frequency optimal control
According to Eqs. ( 11)-( 14) and ( 16)-( 18), the optimization model can be formulated as: Obviously, the above optimization problem is an optimal control problem with state vector constraint.L and f are defined: Then, introduce a new state variable x n (t) and define a function f n , the frequency constraint is expressed as: (15) P WPk max = d% • P WTk max + d% • P PVk max (16) According to the Pontryagin Minimum Principle, the Hamilton function is constructed as: where λ and β are Lagrange multipliers.
The canonical equation can be obtained as: Moreover, the boundary condition and transversal condition can be expressed as: The optimal frequency regulation coefficients can be obtained by solving the minimum of the Hamilton function: The solution of the optimal frequency regulation coefficients can be summarized as: However, it is difficult to solve the analytic expression of the optimal control K WPk .To address this problem, a numerical algorithm-Gaussian Pseudospectral method [23] is used.It approximates the state variables and control variables using a Lagrange interpolating polynomials, which are based on a series of Legendre-Gauss points.Then, the continuous optimal control problem is transformed into a nonlinear programming problem.The sequence quadratic programming algorithm can be used to quickly solve the problem. (23)

WT/PV level-power distribution
After obtaining the control coefficient K WPk from the system controller, the optimal droop coefficients of both WTs and PVs are solved by using the ADMM algorithm in the WT/PV level controller.

WT and PV model
1) WT model.
The mechanical power of WT can be expressed as [23]: where ρ is air density; R w is the blade length; v w is the wind speed; C p is the utilization coefficient of wind energy; λ is the tip speed ratio; β is the pitch angle; and ω r is the rotor speed of WT. when β = 0, C p can be fitted as a second-order polynomial function [24].
Substituting ( 22) into ( 21), the mechanical output power of WT can be obtained as: where, In this paper, the power reserve of WT is obtained by rotor speed-based deloading control.The deloading power of WT can be expressed as: where d w is the deloading percentage and C pmax is the maximum of C p .
The relationship between the system frequency and WT output power can be expressed as: where K w is the droop coefficient of WT and T wf is the response time constant of WT active power.
The electrical power can be expressed as: (29) The swing equation of WT is represented as: where H w is the inertia constant of WT.Equation ( 38) is linearized around the initial operating point and can be obtained as: Equation (36) can be written as: According to Eqs. ( 39) and ( 40), the state-space equation of the ith WT in a matrix form can be written as: with.
x Wi = [�ω ri , �P e2i ] T , u Wi = P wui , y Wi = [�ω ri , �P e2i ] T , where x w is the state variable of the WT model, u w is the control variable of the WT model, and y w is the output variable of the WT model.A w , B w , E w , and C w are the corresponding matrices.
2) PV model.In this paper, a two-stage PV system with frequency regulation capability is adopted [25].For the PV array, a simplified model can be expressed as [26]: where N pv is the number of PV modules; P m is the maximum power point of the PV array; S is the solar irradiance; S ref is the standard solar irradiance; b is the relevant constant of the battery; e is the base of the natural logarithm.
The increment power reference of the PV array can be expressed as: (37) P e = P e1 + P e2 where K p is the droop coefficient of PV.
For the DC/DC converter, two control loops are adopted [25].The inner loop is a PI controller which regulates the voltage of the PV array Δu pv .The outer power control loop is a PI compensator that regulates the output power of the PV array ΔP pv .Thus, the power control loop of the DC/DC converter can be simplified as shown in Fig. 3.
When the active power reference changes, the dynamic process of the DC/DC converter can be obtained as: In order to advantageously formulate the state-space equation, a variable ΔP mid is introduced, and it can be defined as: where k p and k i are the parameters of the outer loop controller; T f and T P are the response time of active power and inner loop controller, respectively; and i pv0 is the PV output current.
According to (45)-(47), the state-space equation of the ith PV in a matrix form can be obtained as: with.
x Pi = u pvi , P pvi , P midi T , u Pi = P pui , y Pi = P pvi, where x p is the state variable of the PV model, u p is the control variable of the PV model, and y p is the output variable of the PV model.A p , B p , E p , and C p are the corresponding matrices, respectively.
3) Model of WT and PV.
Based on the model of WT and PV, the state-space equation of N1 WT and N2 PV can be obtained as: with. x where x 2 , u 2 , and y 2 are the state variable, control variable, and output variable respectively.A 2 , B 2 , E 2 , and C 2 are the corresponding matrices, respectively.

Objective function
The objective function includes two parts.One part is to reduce the control cost of WT, and the other part is to reduce the control cost of PV.
For the WT, the objective function is written as: For the PV, the objective function is written as: The overall objective function can be obtained as: where α 1 , α 2 , and α 3 are the weight coefficients.The weight coefficients α 1 , α 2 , and α 3 can be calculated as: (49 where ΔP WTmax is the maximum power reserve provided by WT and ΔP PVmax is the maximum power reserve provided by PV.
The frequency regulation limit of WT can be expressed as: 2) Constraint of PV.
The frequency regulation limit of PV can be expressed as: 3) Total frequency regulation power constraint

MPC problem formulation
According to Eqs. ( 49), (52), and ( 54)-( 56), the MPC problem is formulated and the optimization model can be expressed as: The above MPC problem can be solved by the ADMM algorithm.

ADMM solving
The above MPC problem can be solved at the WT/PV level using centralized optimization algorithms.However, the number of distributed energy resources in the power system is large, and the problem of dimension disaster may occur in the calculation process of the traditional centralized control method.In addition, due to the long calculation time of traditional centralized control algorithms, it will produce the control delay, which eventually affects the control effect.To deal with these issues, the ADMM algorithm is used.ADMM is an algorithm that it can decompose the objective function of the global optimization problem into several sub-problems equitably, and then obtain (54) the solution of the global problem by solving each sub-problem in parallel and coordinating the solution of the sub-problems [27].Thus, solving the large-scale and distributed optimization problems can significantly improve the computation efficiency.
For the above MPC problem, the decomposed optimization problems and corresponding constraints are solved in parallel in WT/BESS controllers.The above MPC problem can be rewritten as: where x and z are the optimization variables; f(x) is the objective function; g(z) is the indicator function; A is the coefficient matrix; and m is the vector.
The augmented Lagrangian can be expressed as: where β is the Lagrange multiplier; ρ is the penalty parameter; and ρ > 0.
The steps of the detailed solution are as follows.
Step 2: Update z.The augmented Lagrangian with supplementary variable z is minimized by the WT/PV level controller and can be calculated as: Step 3: Update x.Firstly, z i+1 is sent to the WT and PV.Secondly, the sub-problem with constraints is individually solved by the local controller of WTs and PVs in parallel.The augmented Lagrangian with variable x is minimized and can be calculated as: where x Wk is the variable of the kth WT and x Pk is the variable of the kth PV.
Step 4: Update β. β can be calculated as: Step 5: Convergence condition.The original residuals and dual residuals are defined as: (58) The condition of the iteration ending can be expressed as: where ε o and ε d are constants for the feasibility tolerances of the primal and dual feasibility conditions, respectively.

Control process
The process of the control can be illustrated in Fig. 4.
1) After a contingency event is detected, the disturbance amount ΔP L can be estimated by the method [21].Then, the frequency nadir f na and steady-state frequency f ∞ are predicted by the MMFR model, and the maximum output power of WT and PV can be calculated.2) When f na < f th and f ∞ < f ∞ th , the frequency stability control scheme is activated.
(66) 3) For the system level, the optimization problem ( 19) is solved to obtain the optimal frequency regulation coefficients of WT and PV clustering systems, and the results are sent to the WT/PV level controllers.4) For the WT/PV level, the MPC-based optimization problem (57) is solved by the ADMM algorithm, and the frequency regulation coefficients of each WT and PV are determined in WT/PV controllers.

Results
In this section, the validity of the proposed control scheme is verified in the modified IEEE 10-machine 39-bus system.The IEEE 10-machine 39-bus system is used as the main network, and the distribution network (DN) is connected to the bus 4, 20, 24, 25, and 28, respectively.The total capacity of WTs and PVs in the distribution networks is 800 MW, 900 MW, 850 MW, 750 MW, and 800 MW, respectively, and the de-loading operation is used for all WTs and PVs.The penetration level of distributed energy resources is 45.24%.The structure of the improved IEEE 10-machine 39-bus system is shown in Fig. 5.In this paper, MATLAB is used to solve the optimization problem.The result is verified in the simulation software PSS/E.The measurement data is obtained by PSS/E to simulate the WAMS information.In addition, the central processing unit of the computer is an Intel Core i9-14900KF, and it has 64 GB of memory.Moreover, two other control schemes are also applied for WT and PV to test the control performance of the proposed approach.The two control schemes are described as follows: Control method 1 is presented in [28] as the conventional droop control strategy with a fixed frequency regulation coefficient.The method is implemented by controlling the injected additional active power from the power electronics interfaced energy resources, which is proportional to the frequency deviation.
Control method 2 is presented in [29] as the variable coefficient control, which one constant frequency regulation coefficient is used to control the TF, then the other constant frequency adjustment coefficient is used to control the SSF, and the switching of the two coefficients is realized by a fixed rate.

Case: generator fault
The generator G4 trips at 1 s, and it results in the active power loss of 820 MW.All distributed energy resources in the system participate in the frequency control.Moreover, we define the limit that the TF deviation and SSF deviation are maintained within 0.5 Hz and 0.2 Hz.

Prediction
In order to verify the performance of the proposed MMFR model, the SFR and PSS/E simulations are compared.The frequency response curves of MMFR, SFR, and PSS/E are shown in Fig. 6, and the key characteristics among MMFR model, SFR model, and PSS/E are compared in Table 1.
It is seen from Fig. 6 that the frequency response curve of the MMFR model is closer to the PSS/E model compared with the SFR model.Moreover, to be more intuitive, Table 1 shows the absolute error and relative error of key characteristics compared Fig. 6 The frequency response curve with PSS/E.It is seen that the error of the MMFR model is smaller than that of the SFR model.Therefore, the MMFR model has higher precision, and the key characteristics of the frequency, which includes the FN and SSF, can be accurately described after a severe disturbance.Besides, since both the FN and SSF are out of the frequency stability thresholds, the proposed control is activated.

Control performance
The frequency curves of the system and the curves of frequency regulation coefficients with three control methods are shown in Figs.7 and 8.
It is seen from Fig. 7 that the TF and SSF of the proposed method can be accurately restored above 49.5 Hz and 49.8 Hz, and the proposed method meets the requirements of frequency stability.To much better show the control effect, the error of the FN and SSF is calculated.The relative errors are 0.0062 Hz and 0.0048 Hz, and the absolute errors are 0.0125% and 0.0097%, respectively.The results show that the control accuracy of the TF and SSF of the proposed method is higher.Moreover, it is seen from Fig. 8 that the distributed energy resources are controlled to quickly inject the active power into the system after disturbance and repress the frequency decline successfully.The frequency nadir is controlled above 49.5 Hz.Then, the frequency regulation coefficient of the distributed energy resources begins to decline after about 2.15 s.When the coefficient is less than or equal to the coefficient of the steady-state frequency stability, the TF control is ended, and the constant coefficient is provided by the distributed energy resources to control the SSF after this.Therefore, the control objection of both TF and SSF are met simultaneously and the effectiveness and correctness of the proposed approach is verified.The frequency regulation coefficients of five distribution networks are shown in Fig. 9.It is also seen that the frequency regulation coefficients are assigned according to the total reserve capacity of the distributed energy resources in the distribution network, and the distribution network with strong reserve capacity bears more.In addition, the output power of PVs and WTs in the distribution network 1 and 3 are shown in Figs. 10  and 11.Similarly, the output power of PV/WT is determined by their reserve capacity, and the distributed energy resource with large reserve capacity bears more.

Consuming control energy
Compared with the two methods, the proposed approach needs less energy.For method 1, the distributed energy resources with a constant frequency regulation coefficient participate in frequency control, resulting in a significant increase in energy consumption.For method 2, although the variable frequency regulation coefficient is used to control the frequency, the TF is controlled by a constant frequency regulation coefficient, which leads to more energy usage.Although the three methods can recover the system frequency above the stability thresholds, there are obvious differences in control performance, as indicated in Table 2.It is seen that the proposed approach requires the least control energy.Furthermore, compared with the two approaches, the energy consumption of the proposed approach is decreased by 37.64% and 24.56%.It is indicated that the proposed approach can improve the utilization efficiency of control energy. 2 is relatively high.In addition, the total control cost is also calculated, as shown in Table 3.Compared with the two traditional approaches, the control cost of the proposed scheme is decreased by 30.83% and 19.64%.Thus, it can be concluded that the proposed control strategy reduces the control cost due to the coordinated control of the distributed energy resources.

Solving efficiency
The solving time of the optimal control and ADMM-based solution in each steptime simulation is shown in Fig. 13.The steptime T s is 0.1 s in this paper.The simulation period is set from the start to the end of transient frequency control (1.1 ~ 4 s), and a total of 30 steps of simulation are carried out.It is seen from Fig. 13a that the solution time of each steptime simulation is between 0.02 and 0.05 s, and the average solution time is 0.0363 s for the optimal control.Moreover, the ADMM-based solving time of distribution networks 2 and 5 are shown in Fig. 13b and c.For distribution network 2, the solution time is between 0.03 and 0.06 s, and the average solution time is 0.042 s.For the distribution network 5, the solution time is between 0.03 and 0.06 s, and the average solution time is 0.043 s.The total average solving time is 0.0793 s, and it is less than the steptime of 0.1 s.So the solving speed of the optimal control and ADMM-based solution can meet the requirements of real-time control.

Impact of non-linearity
In the primary frequency response process, the non-linearity factor of the synchronous generator is mainly the dead band.In the frequency stability control, the dead band is hardly handled, so it needs to be simplified.For the dead band of the governor, the linear approximation method is usually adopted.It is known that the actual governor response is slightly smaller than the simplified governor response.It means that the calculated power is slightly larger.However, the dead band is much less than the frequency deviation caused by a large disturbance.For example, the dead band of the governor is 0.033 Hz in China [30], while the frequency deviation caused by the generator's failure to trip is generally more than 0.5 Hz.So, the impact of the dead band can be ignored.Moreover, the calculated error caused by the dead band can be compensated by the frequency safety margin in Eqs.(17) and (18) to ensure control performance.

Impact of time
In order to ensure the performance of the proposed control, the total time including the evaluation time, the optimization time, and communication delay is analyzed.Firstly, the imbalance power, FN and SSF are estimated, and the evaluation time is < 20 ms [21].Secondly, the optimization solution time of the coordinated control scheme is < 100 ms.Finally, the communication delay is considered, and it is < 80 ms [31].Therefore, the total time is < 200 ms, which does not have too much impact on the control performance, and it can meet the requirements of real-time control.

Application scenarios
In this paper, the proposed method can control the TF and SSF above the frequency stability threshold.However, it is only applied to one of the scenarios, in which both the prediction value of FN and SSF are less than the preset frequency stability threshold, i.e., f na < f th and f ∞ < f ∞ th .For other scenarios, the proposed control method is still applicable.
When f na < f th and f ∞ > f ∞ th , it means that only transient frequency needs to be controlled.For the system level, the constraint of steady-state frequency is removed from the optimization model (19).For the WT/PV level, the optimization problem and solution method remain the same.
When f na > f th and f ∞ < f ∞ th , it means that only steady-state frequency needs to be controlled.For the system level, the constraint of transient frequency is removed from the optimization model (19).For the WT/PV level, the optimization problem and solution method remain the same.
When f na > f th and f ∞ > f ∞ th , it means that both the prediction value of FN and SSF are more than the preset frequency stability threshold, and it does not need to participate in frequency stability control.

Conclusions
In this paper, a hierarchical control strategy for the distributed energy resources is presented to simultaneously improve the TF and SSF.We can get the following conclusions.
1) The proposed MMFR model can precisely describe the key characteristics of system frequency after disturbance and has a good control performance.2) For the frequency stability control level, the optimal frequency regulation coefficients of multiple WT and PV clustering systems are determined by the optimal control principle.Compared with the traditional methods, the control energy is decreased by 37.64% and 24.56%, respectively.It shows that the control energy of the proposed method can be used efficiently.3) For the WT/PV power distribution level, considering the frequency response model of WT and PV, the frequency regulation coefficients of WTs and PVs are coordinated by MPC to reduce the frequency control cost, which it is solved by the ADMM algorithm.Compared with the traditional methods, the control cost is decreased by 30.83% and 19.64%.It shows that the proposed method is more economical in cost control.
The proposed control method significantly improves the frequency stability.Meanwhile, it reduces the control cost and achieves the optimal allocation of frequency regulation power.However, the proposed method neglects the influence of power loss and node voltage change on the system frequency.Moreover, it also does not consider the impact of cyber-attacks.The authors will continue to carry out in-depth research in the future.

Fig. 2
Fig. 2 MMFR model with distributed energy resources

Fig. 3
Fig. 3 Power control loop of converter

Fig. 8
Fig. 8 Total frequency regulation coefficients with three control methods

Fig. 11
Fig. 11 Output power of WTs in the distribution network 3

Fig. 12 Fig. 13
Fig. 12 Control cost with three control methods

Table 1
Comparison results for FN and SSF

Table 2
Energy of distributed energy resources participate in frequency control

Table 3
Cost of distributed energy resources participates in frequency control