Dynamic modeling of multiple microgrid clusters using regional demand response programs

: Preserving the frequency stability of multiple microgrid clusters is a serious challenge. This work presents a dynamic model of multiple microgrid clusters with di ﬀ erent types of distributed energy resources (DERs) and energy storage systems (ESSs) that was used to examine the load frequency control (LFC) of microgrids. The classical proportional integral derivative (PID) controllers were designed to tune the frequency of microgrids. Furthermore, an imperialist competitive algorithm (ICA) was proposed to investigate the frequency deviations of microgrids by considering renewable energy resources (RERs) and their load uncertainties. The simulation results conﬁrmed the performance of the optimized PID controllers under di ﬀ erent disturbances. Furthermore, the frequency control of the microgrids was evaluated by applying regional demand response programs (RDRPs). The simulation results showed that applying the RDRPs caused the damping of frequency ﬂuctuations.


Introduction
Recently, air pollution produced by fossil fuel power plants has been causing serious environmental problems [1].Due to environmental concerns, the tendency toward using renewable energy resources (RERs) for power generation is increasing [2].Generally, RERs increase the reliability of the system, improve the voltage profile, and decrease the dependence on fossil fuel [3][4][5][6].By increasing the influence of RERs, an imbalance between generation and demand occurs in power systems because of the stochastic nature of these resources [7].To overcome this problem, the use of energy storage systems (ESSs) can be a proper solution.ESSs can inject or store energy in power systems when it is required.Furthermore, they improve the reliability, stability, and power quality of the systems [8].The RERs have a small capacity.Therefore, they are linked to the grid at the distribution level.A distribution system that has small-scale energy resources is called a microgrid, which is also composed of loads, distributed energy resources (DERs), ESSs, and other infrastructure [9].Generally, microgrids connect to the main grid.However, situations such as errors, voltage fluctuations, and high oscillations of the main grid are causing them to separate from the main grid and they act in isolation [10].During the isolated operation of the main grid, the controllable power resources respond to power and load disturbances by injecting or absorbing power.Furthermore, they adjust the frequency and power of the microgrid.Due to the inherent delay and the generation rate constant (GRC) of power resources, the control instructions applied to the power resources cannot be implemented immediately.As such, frequency fluctuations will occur in the microgrids and the proper adjustment of control parameters can help the damping of frequency oscillations [11,12].Since the power resources are operating in parallel, to acquire a stable feedback loop, it is necessary to get the speed-droop characteristic (f-P) for the power resources to participate in the frequency control.
In addition, for the power-sharing between resources, it is necessary to obtain the reverse relation of frequency droop coefficient (R) by considering the nominal value of that power resource [13,14].Similarly, each power resource has a control signal for adjusting the output power in a manner that is proportional to the frequency oscillations, which is called frequency bias (β).Nowadays, many papers have presented the dynamic modeling of islanded microgrids [15][16][17][18][19][20].In Abazari et al. [16], a new load frequency control (LFC) model for an independent hybrid microgrid with RERs was presented and the dynamic models for the resources were introduced.In addition, their participation in frequency control studies was considered.In A-Hamed et al. [18], the performance of a microgrid in the presence of a new control method was evaluated.The proposed microgrid had various power generation resources, including photovoltaic, wind, diesel generator, fuel cell and electrolyzer, and battery resources.A model for the microgrid was implemented through Matlab/Simulink software and the integral proportional (PI) control scheme was implemented.In Fini and Golshan [20], a microgrid with controllable and uncontrollable resources was considered.The stability of the microgrid frequency was investigated by considering the virtual inertia of the ultracapacitor unit.
Due to the position of microgrids in the development of future energy substrates, maintaining the security of microgrid electrification, especially in critical situations, is a serious challenge.One way to maintain the security of microgrid electrification is by connecting microgrids.The multiple microgrid clusters are supporting each other and increasing the reliability of the whole system [21].In Toro and M-Nava [22], a frequency control strategy for multiple microgrid clusters that are connected by a tie-line was suggested.When load disturbances occur, the proportional integral derivative (PID) controllers are used to control the system frequency.On the other hand, the uncertainties in the generation of RERs make the PID controllers unable to perform the frequency control operation.Therefore, the use of intelligent techniques, such as a genetic algorithm and the imperialist competitive algorithm (ICA), to adjust the optimal parameters are suggested in References [23,24].One of the factors used to engender instability is by increasing the over-frequency in power systems.Increasing the over-frequency is damaging to system equipment and it can cause outages.Therefore, demand response programs can avoid these issues and enhance power systems security [25].References [26,27] showed that demand response programs could play an important role in the frequency regulation of the power system.In Moghadam et al. [28], an algorithm of distributed frequency control was studied through randomizing the frequency response of smart appliances.In this method, the controllers of the generation side are deactivated and the demand response is activated.
The participation of demand response programs in the frequency control of isolated microgrids and power systems was investigated in References [29,30].In Huang and Li [31], the effect of controllable loads on the frequency response of the system and load damping coefficient was examined.In Malik and Ravishankar [32], domestic loads were considered as a source of demand response during frequency control.In Soares et al. [33], the use of demand response programs to integrate the RERs and reduce the effects of changes in the generation of these resources on power systems was proposed.In C-Chien et al. [34], the collaboration of a spinning reserve and demand response programs in frequency restoration during contingencies of the system was suggested.
Despite the high participation of demand response programs in frequency adjustments, many operating concerns remain.On the one hand, overload shedding can cause unnecessary load shedding and consequently lead the system to excessive over-frequency.On the other hand, the light usage of demand response programs during system faults may reduce the positive impact of such a method.Therefore, knowing the location of the load shedding and the exact magnitude of the disturbances, particularly in multiple microgrids, is necessary [35,36].In Bayat et al. [37], a method Energies 2020, 13, 4050 3 of 20 was presented that used the changes of reactive power in buses to locate the changes of the active power consumption in controllable loads.In Gao and Ning [38], the wavelet analysis method was used to assess the location and size of the disturbances.In Babahajiani et al. [39], a method for determining the location of load disturbances and the magnitude of disturbances in interconnected power systems was suggested.,By monitoring the power flows of the tie-lines, the proposed method determines the magnitude of disturbance and the region that disturbance occurs in and applies the response to the demand program for the affected region.
In light of the aforementioned literature, the contributions of this work can be explained as follows: • The absence of a suitable dynamic model is one of the serious challenges for a study on LFC in multiple microgrid clusters, which restricts research activities in this field.Modeling can be done for a variety of goals.Thus, using a suitable model for a specific study goal is the preferred option.
The present paper provides a suitable dynamic model for multiple microgrid clusters based on tie-lines for studying the frequency and tie-line power control.

•
In this study, the Monte Carlo simulation was used to generate scenarios regarding the uncertainties of RERs and loads by considering probability distribution functions when modeling them.However, due to the uncertainties in RERs and loads, the conventional controllers, such as PID controllers, will not be suitable for LFC.Hence, adjusting the control parameters was suggested by using the ICA to solve this problem.

•
Damping frequency fluctuations is one of the main issues of frequency control in power systems and microgrids.Due to its fast dynamics, the demand response program could be an attractive suggestion for damping system frequency fluctuations.Opposed to the other works, the regional demand response program is simpler and more accurate, and does not require complex and massive computational calculations.Furthermore, this method is applicable to almost all types of loads by considering the load's participation coefficient in each region during the frequency control.Therefore, the regional demand response was applied for each microgrid in this study to reduce the frequency oscillations.
The rest of this paper is organized as follows.Section 2 presents the case study.In Section 3, the modeling of the units is explained.In Section 4, the modeling of the microgrids is depicted.Section 5 presents the modeling of multiple microgrid clusters.The modeling of uncertainties is explained in Section 6.In Section 7, the modeling of regional demand response programs is discussed.The simulation results and discussions are provided in Section 8. Finally, conclusions are given in Section 9.

Case Study
Here, the multiple microgrid clusters, which are isolated from the main grid and are considered for this case study, are presented.The studied multiple microgrid clusters are shown in Figure 1.In this study, three microgrids with different specifications were connected to each other using a tie-line.The first microgrid consisted of a diesel engine generator, a wind turbine generator, a flywheel energy storage system, and a superconducting magnetic energy storage system.The second microgrid included a diesel engine generator, a battery energy storage system, and an ultracapacitor.The third microgrid consisted of a diesel engine generator, a fuel cell, a photovoltaic panel, and a battery energy storage system.In addition, the RERs and load uncertainties were considered in this case study.Furthermore, the regional demand response programs (RDRPs) were applied for each microgrid.

Modeling of Units
The modeling of the units, including DERs and ESSs, are presented below.The DERs consisted of a diesel engine generator, a fuel cell, a wind turbine generator, and a photovoltaic panel.The ESSs consisted of a flywheel energy storage system, a superconducting magnetic energy storage system, a battery energy storage system, and an ultracapacitor.

Modeling of the Diesel Engine Generator (DEG)
A DEG was used to compensate the microgrid required power by using a governor and a speeddrop system.The DEG conversion function is represented in Equation (1) [20]: where Tg and Tt are the governor and turbine time constants, respectively.Furthermore, RDEG is the speed regulation coefficient of the DEG.Generally, a power increment limitation block is considered for the DEG model.

Modeling of the Fuel Cell (FC)
The FC transforms the chemical energy stored in hydrogen into DC electrical energy.The FC begins to produce power by injecting hydrogen stored in its hydrogen tank.The model used for the FC was [40]: where TFC1 and TFC2 are the FC1 and FC2 feed time constants, respectively.The power generation increase rate limit was also applied for the FC model.

Modeling of Units
The modeling of the units, including DERs and ESSs, are presented below.The DERs consisted of a diesel engine generator, a fuel cell, a wind turbine generator, and a photovoltaic panel.The ESSs consisted of a flywheel energy storage system, a superconducting magnetic energy storage system, a battery energy storage system, and an ultracapacitor.

Modeling of the Diesel Engine Generator (DEG)
A DEG was used to compensate the microgrid required power by using a governor and a speed-drop system.The DEG conversion function is represented in Equation (1) [20]: where T g and T t are the governor and turbine time constants, respectively.Furthermore, R DEG is the speed regulation coefficient of the DEG.Generally, a power increment limitation block is considered for the DEG model.

Modeling of the Fuel Cell (FC)
The FC transforms the chemical energy stored in hydrogen into DC electrical energy.The FC begins to produce power by injecting hydrogen stored in its hydrogen tank.The model used for the FC was [40]: where T FC1 and T FC2 are the FC1 and FC2 feed time constants, respectively.The power generation increase rate limit was also applied for the FC model.

Modeling of the Wind Turbine (WT) Generator
In this study, the WT generator was not dispatchable.The WT generator was modeled as a specified random distribution.Furthermore, it was represented by the following transfer function [40]: where T WTG is the WT generator time constant.

Modeling of the Photovoltaic (PV) Panel
In this study, the PV panel was not dispatchable.The PV panel was modeled as a specified random distribution and it was modeled using the following transfer function [40]: where T PV is PV panel time constant.

Modeling of the Flywheel Energy Storage System (FESS)
The FESS plays an important role in a microgrid since it performs both the functions of releasing and storing energy at an appropriate time.Generally, a FESS is considered to provide mechanical energy.It is mostly identified using two various time constants: the electrical time constant and the mechanical time constant, which are identified by the electronic parts of the technology and the inertia of the system, respectively.The FESS generates active and reactive powers based on the angular speed of its rotational mass.Therefore, the speed control loop is added to the FESS model to provide an exact dynamic model of the FESS in the LFC scheme.As a mechanical rotational mass, the FESS has inertia.Due to the influence of the number of poles on the speed of this equipment, the number of poles should be applied to the mechanical part of this equipment [16].
The dynamic model of an FESS is represented in Figure 2, where T F1 , T F2 , and T F3 are the time constants of the measurement device, command device, and converter, respectively.F 1 is a reference angular speed.F 2 is the speed regulator proportional constant.P FESS is the number of poles of the FESS, and J FESS is the inertia of the rotational masses of the FESS.

Modeling of the Wind Turbine (WT) Generator
In this study, the WT generator was not dispatchable.The WT generator was modeled as a specified random distribution.Furthermore, it was represented by the following transfer function [40]: where TWTG is the WT generator time constant.

Modeling of the Photovoltaic (PV) Panel
In this study, the PV panel was not dispatchable.The PV panel was modeled as a specified random distribution and it was modeled using the following transfer function [40]: where TPV is PV panel time constant.

Modeling of the Flywheel Energy Storage System (FESS)
The FESS plays an important role in a microgrid since it performs both the functions of releasing and storing energy at an appropriate time.Generally, a FESS is considered to provide mechanical energy.It is mostly identified using two various time constants: the electrical time constant and the mechanical time constant, which are identified by the electronic parts of the technology and the inertia of the system, respectively.The FESS generates active and reactive powers based on the angular speed of its rotational mass.Therefore, the speed control loop is added to the FESS model to provide an exact dynamic model of the FESS in the LFC scheme.As a mechanical rotational mass, the FESS has inertia.Due to the influence of the number of poles on the speed of this equipment, the number of poles should be applied to the mechanical part of this equipment [16].
The dynamic model of an FESS is represented in Figure 2, where TF1, TF2, and TF3 are the time constants of the measurement device, command device, and converter, respectively.F1 is a reference angular speed.F2 is the speed regulator proportional constant.PFESS is the number of poles of the FESS, and JFESS is the inertia of the rotational masses of the FESS.

Modeling of the Superconducting Magnetic Energy Storage System (SMESS)
By injecting or absorbing power, an SMESS plays a significant role in increasing the frequency stability of the network.The exchange power of this unit is determined by the PI controller.This unit has a very small delay and it is modeled using a first-order conversion function that is described in Equation ( 5): where K SMESS is the SMESS gain and T SMESS is the SMESS time constant.Furthermore, the power generation increase and decrease rate limit are applied for the model [41].

Modeling of the Battery Energy Storage System (BESS)
A BESS is a basic voltage device that is used to prevent fluctuations in the power balance.The BESS conversion function is expressed in Equation ( 6): T BESS is the BESS time constant.Furthermore, R BESS is the speed regulation coefficient of the BESS.The ramp-up and ramp-down rate of the output power of the BESS were also considered in this model [20].

Modeling of the Ultracapacitor (UC)
The UC plays a significant role in increasing the virtual inertia of the network via power absorption or injection.The exchange power of the UC with the network is directly dependent on the derivation of frequency fluctuations and it is represented by the following transfer function [19]: where T UC is the UC time constant.

Modeling of the Microgrids
According to the presented dynamic models of each DER and ESS, and by considering the type of the resource of each microgrid, the dynamic models of the microgrids are shown in Figures 2-4.In Figure 2, the block diagram of the dynamic model of microgrid 1 is illustrated, where ∆f 1 is the frequency deviation of microgrid 1, β 1 is frequency bias of microgrid 1, D 1 is the load-damping coefficient of microgrid 1, M 1 is the total inertia of microgrid 1, ∆P L1 is the load change in microgrid 1, and ∆P tie,1 is the total tie-line power change between microgrid 1 and the other microgrids.In Figure 3, the block diagram of the dynamic model of microgrid 2 is given, where ∆f 2 is the frequency deviation of microgrid 2, β 2 is the frequency bias of microgrid 2, D 2 is the load-damping coefficient of microgrid 2, M 2 is the total inertia of microgrid 2, ∆P L2 is the load change in microgrid 2, and ∆P tie,2 is the total tie-line power change between microgrid 2 and the other microgrids.In Figure 4, the block diagram of the dynamic model of microgrid 3 is demonstrated, where ∆f 3 is the frequency deviation of microgrid 3, β 3 is the frequency bias of microgrid 3, D 3 is the load-damping coefficient of microgrid 3, M 3 is the total inertia of microgrid 3, ∆P L3 is the load change in microgrid 3, and ∆P tie,3 is the total tie-line power change between microgrid 3 and the other microgrids.The main parameters of microgrids 1-3 are given in Tables 1-3, respectively.The main parameters of microgrids 1-3 are given in Tables 1-3, respectively.The main parameters of microgrids 1-3 are given in Tables 1-3, respectively.In the first microgrid, the WT generator did not participate in the frequency control.However, DEG, FESS, and SMESS participated in the frequency control.In the second microgrid, DEG and BESS were responsible for the power generation and the UC was used to inject or absorb power quickly to increase the virtual inertia of the microgrid.In addition, in the third microgrid, the PV generator did not participate in the frequency control.Nevertheless, the DEG, FC, and BESS participated in the frequency control.In each microgrid, two PID controllers were used for the LFC of the microgrids.

Modeling of the Multiple Microgrid Clusters
In this study, three microgrids were connected to each other with tie-lines.The exchanged tie-line power between two microgrids depends on the difference in frequency deviations of those microgrids and the tie-line synchronizing torque coefficient.Each of microgrids was connected to each other by two tie-lines and the total power of those two tie-lines was the exchange tie-line power of that microgrid.The dynamic model of multiple microgrid clusters is shown in Figure 5, where T 12 , T 13 , T 21 , T 23 , T 31 , and T 32 are the tie-line synchronizing torque coefficients, which were set to Energies 2020, 13, 4050 8 of 19 In the first microgrid, the WT generator did not participate in the frequency control.However, DEG, FESS, and SMESS participated in the frequency control.In the second microgrid, DEG and BESS were responsible for the power generation and the UC was used to inject or absorb power quickly to increase the virtual inertia of the microgrid.In addition, in the third microgrid, the PV generator did not participate in the frequency control.Nevertheless, the DEG, FC, and BESS participated in the frequency control.In each microgrid, two PID controllers were used for the LFC of the microgrids.

Modeling of the Multiple Microgrid Clusters
In this study, three microgrids were connected to each other with tie-lines.The exchanged tieline power between two microgrids depends on the difference in frequency deviations of those microgrids and the tie-line synchronizing torque coefficient.Each of microgrids was connected to each other by two tie-lines and the total power of those two tie-lines was the exchange tie-line power of that microgrid.The dynamic model of multiple microgrid clusters is shown in Figure 5, where T12, T13, T21, T23, T31, and T32 are the tie-line synchronizing torque coefficients, which were set to T12 = T13 = T21 = T23 = T31 = T32 = 1.

Modeling of the Uncertainties
In power systems, there are many uncertainties and the definite parameters cannot indicate the exact status of these systems.Therefore, probabilistic analysis was considered as an attractive solution.In this study, RERs and load uncertainties were modeled with probability distribution functions.

Probabilistic Model of the Load
Weather conditions and time are two factors for the deterministic component of the load change.
Due to the variable nature of the load, the random behavior of the load was suitable for modeling.

Modeling of the Uncertainties
In power systems, there are many uncertainties and the definite parameters cannot indicate the exact status of these systems.Therefore, probabilistic analysis was considered as an attractive solution.In this study, RERs and load uncertainties were modeled with probability distribution functions.

Probabilistic Model of the Load
Weather conditions and time are two factors for the deterministic component of the load change.Due to the variable nature of the load, the random behavior of the load was suitable for modeling.This modeling could be calculated using measured data analysis.The load behavior was modeled using a normal distribution function.The normal distribution function is defined as following [42]: where P L is the load power, µ is the mean, and σ is the standard deviation.

Probabilistic Model of the WT Power
The power produced by the WT generator depends on the wind speed.The wind speed changes regularly, which shows the significance of a probability model.The Weibull distribution function was used in modeling the wind speed.The Weibull distribution function is as follows [43]: where V is the wind speed, and y and x are the scale and shape parameters of the Weibull distribution function, respectively.After producing the wind speed samples, the real power generation by the WT generator can be acquired from the following equation: where P WT (V) is the power generated at wind speed V by the WT generator.P r,WT is the rated power of the WT generator; V ci , V co , and V r are the low cut speed, high cut speed, and rated speed of the WT generator, respectively.

Probabilistic Model of the PV Power
Solar radiation is one important parameter for power generation in the PV panel.This parameter is variable at any time of day.Thus, the location of the PV panel plays an essential role in generating electricity.The irradiance was modeled using the beta distribution function.The PV panel was tested under the standard test condition (STC).After radiation sampling, the final power output of the PV panel was calculated as follows [44]: where P PV (R) is the power generated by the PV panel; P r,PV is the rated power of the PV panel; and R, R C , and R STC are the solar radiation, certain radiation point, and the solar radiation under the STC, respectively.

Scenario Generation and Reduction
The uncertainty of the input parameters, including the generation of RERs and the load demands, was considered by applying the scenario generation.For the uncertainty parameters, a Monte Carlo simulation was used to generate scenarios.Based on Hemmati [45], probability distribution sets were used to model the uncertainties.To reduce the complexity of the problems, a scenario reduction method was utilized.The scenario reduction is a method used for reducing the implementation time of the simulation [46].

Imperialist Competitive Algorithm
Since there are many variables and control signals, we needed to use optimization algorithms to adjust the control parameters.For these algorithms, the specific objective functions must be defined.The objective functions should be optimized according to the constraints to determine the optimal value of the control parameters.In this study, the ICA was used as described below.
The ICA is an evolutionary optimization method.Similar to other evolutionary methods, the ICA produces some random initial population, which is named "country."It is separated into two types, imperialists and colonies, which together form empires.The stages of the suggested algorithm are described as follows: the generation of initial empires, the assimilation, the revolution, the exchange between the best colony and imperialist, the imperialistic competition, and the elimination of a powerless empire [47].

Modeling of Regional Demand Response Programs
Based on Babahajiani et al. [39], an algorithm was utilized for applying the RDRPs in frequency control.The proposed solution involved monitoring the deviations of the tie-line flows, computed the magnitude of the disturbances, and synchronously specified the microgrid where disturbances occurred to apply the demand response programs exactly to the involved microgrid.The proposed algorithm is given as follows: At first, the second derivative of the tie-line flows of all microgrids was computed.Some of the demand response programs were voluntary and contract-based.Therefore, when applying the impact of this limitation, the participation factor 0 < γ < 1 in Equation ( 12) was considered.The participation factor determines how much load could participate in the demand response programs at the time of disturbance.γ = 0 means that no loads participate in frequency control and γ = 1 means that all active loads participate in frequency control.The Federal Energy Regulatory Commission (FERC) in the United States has announced that the residential controlled loads, such as refrigerators, water heaters, air conditioners, and heat pumps, make up about 20% of the total electricity consumption in the United States [48].Furthermore, in Bao and Li [49], it has been declared that electric heaters use about 11% of the total electricity consumption.In this study, for evaluating the proposed algorithm, it was supposed that in each microgrid, 30% of loads were controllable by the demand response programs.The proposed solution in this work was considered as a direct load control (DLC) that was usually voluntary and contract-based.Figure 6 displays how the RDRPs contribute to the frequency control, where ∆P tie,i is the total tie-line power change between microgrid i and the other microgrids, γ i is the participation factor of demand response programs of microgrid i, τ i is the demand response programs time delay of microgrid i, and RDRP i is the calculated load for the demand response programs of microgrid i.

Simulation Results and Discussions
The frequency control system of multiple microgrid clusters was implemented in Matlab R2016b/Simulink software.To investigate the proposed methods, three scenarios were examined:


Scenario 1: In the first scenario, the performance of the multiple microgrid clusters was studied with non-optimal parameter values. Scenario 2: In the second scenario, the performance of the multiple microgrid clusters was examined with the optimal control parameters. Scenario 3: In the third scenario, the performance of the multiple microgrid clusters was investigated by applying RDRPs.
These scenarios show the effectiveness of the selected solutions at overcoming frequency deviations and power tie-line deviations.In all the scenarios, the RERs and loads uncertainties were considered in the same profile.In the first microgrid, the production of the WT generator had uncertainty and it was modeled using the Weibull distribution function.In the second microgrid, the power load had uncertainty and it was modeled using a normal distribution function.In the third microgrid, the production of the PV panel had uncertainty and it was modeled using the beta distribution function.By applying uncertainties, the output powers of the WT generator and the PV panel were generated and are shown in Figures 7 and 8, respectively.In addition, the power load uncertainty is shown in Figure 9.

Simulation Results and Discussions
The frequency control system of multiple microgrid clusters was implemented in Matlab R2016b/Simulink software.To investigate the proposed methods, three scenarios were examined: • Scenario 1: In the first scenario, the performance of the multiple microgrid clusters was studied with non-optimal parameter values.• Scenario 2: In the second scenario, the performance of the multiple microgrid clusters was examined with the optimal control parameters.• Scenario 3: In the third scenario, the performance of the multiple microgrid clusters was investigated by applying RDRPs.
These scenarios show the effectiveness of the selected solutions at overcoming frequency deviations and power tie-line deviations.In all the scenarios, the RERs and loads uncertainties were considered in the same profile.In the first microgrid, the production of the WT generator had uncertainty and it was modeled using the Weibull distribution function.In the second microgrid, the power load had uncertainty and it was modeled using a normal distribution function.In the third microgrid, the production of the PV panel had uncertainty and it was modeled using the beta distribution function.By applying uncertainties, the output powers of the WT generator and the PV panel were generated and are shown in Figures 7 and 8, respectively.In addition, the power load uncertainty is shown in Figure 9.

Simulation Results and Discussions
The frequency control system of multiple microgrid clusters was implemented in Matlab R2016b/Simulink software.To investigate the proposed methods, three scenarios were examined: • Scenario 1: In the first scenario, the performance of the multiple microgrid clusters was studied with non-optimal parameter values.• Scenario 2: In the second scenario, the performance of the multiple microgrid clusters was examined with the optimal control parameters.• Scenario 3: In the third scenario, the performance of the multiple microgrid clusters was investigated by applying RDRPs.
These scenarios show the effectiveness of the selected solutions at overcoming frequency deviations and power tie-line deviations.In all the scenarios, the RERs and loads uncertainties were considered in the same profile.In the first microgrid, the production of the WT generator had uncertainty and it was modeled using the Weibull distribution function.In the second microgrid, the power load had uncertainty and it was modeled using a normal distribution function.In the third microgrid, the production of the PV panel had uncertainty and it was modeled using the beta distribution function.By applying uncertainties, the output powers of the WT generator and the PV panel were generated and are shown in Figures 7 and 8, respectively.In addition, the power load uncertainty is shown in Figure 9.   Furthermore, the magnitude of the load demand in the first microgrid at t = 50 (s) was equal to 0.05 per unit and the magnitude of the load demand in the third microgrid at t = 100 (s) was equal to 0.1 per unit.Figure 10 shows the load power changes of microgrids 1 and 3.   Furthermore, the magnitude of the load demand in the first microgrid at t = 50 (s) was equal to 0.05 per unit and the magnitude of the load demand in the third microgrid at t = 100 (s) was equal to 0.1 per unit.Figure 10 shows the load power changes of microgrids 1 and 3.  Furthermore, the magnitude of the load demand in the first microgrid at t = 50 (s) was equal to 0.05 per unit and the magnitude of the load demand in the third microgrid at t = 100 (s) was equal to 0.1 per unit.Figure 10 shows the load power changes of microgrids 1 and 3.  Furthermore, the magnitude of the load demand in the first microgrid at t = 50 (s) was equal to 0.05 per unit and the magnitude of the load demand in the third microgrid at t = 100 (s) was equal to 0.1 per unit.Figure 10 shows the load power changes of microgrids 1 and 3.By applying the uncertainties, the frequency deviations the microgrids, non-optimal and optimal parameters, and RDRPs were simulated and the results of these simulations are shown in Figures 11-13.
Energies 2020, 13, 4050 13 of 19 By applying the uncertainties, the frequency deviations of the microgrids, non-optimal and optimal parameters, and RDRPs were simulated and the results of these simulations are shown in Figures 11-13.By applying uncertainties, the exchange tie-line power deviations between the microgrids, nonoptimal and optimal parameters, and RDRPs were simulated and are shown in Figures 14-16.By applying the uncertainties, the frequency deviations of the microgrids, non-optimal and optimal parameters, and RDRPs were simulated and the results of these simulations are shown in Figures 11-13.By applying uncertainties, the exchange tie-line power deviations between the microgrids, nonoptimal and optimal parameters, and RDRPs were simulated and are shown in Figures 14-16  By applying the uncertainties, the frequency deviations of the microgrids, non-optimal and optimal parameters, and RDRPs were simulated and the results of these simulations are shown in Figures 11-13.By applying uncertainties, the exchange tie-line power deviations between the microgrids, nonoptimal and optimal parameters, and RDRPs were simulated and are shown in Figures 14-16  By applying uncertainties, the exchange tie-line power deviations between the microgrids, non-optimal and optimal parameters, and RDRPs were simulated and are shown in Figures 14-16.In the first scenario, the PID controllers were used for the LFC of the microgrids and the nonoptimal parameters were acquired using a trial and error method.The non-optimal values are given   In the first scenario, the PID controllers were used for the LFC of the microgrids and the nonoptimal parameters were acquired using a trial and error method.The non-optimal values are given    In the first scenario, the PID controllers were used for the LFC of the microgrids and the nonoptimal parameters were acquired using a trial and error method.The non-optimal values are given In the first scenario, the PID controllers were used for the LFC of the microgrids and the non-optimal parameters were acquired using a trial and error method.The non-optimal values are given in Table 4, where K p is a proportional gain parameter, K i is an integral gain parameter, and K d is a differential gain parameter, and the subscript numbers refer to the relevant microgrids.

Figure 5 .
Figure 5. Dynamic model of the multiple microgrid clusters.

Figure 5 .
Figure 5. Dynamic model of the multiple microgrid clusters.

Figure 11 .
Figure 11.Comparison of the frequency deviations of microgrid 1 in the three scenarios.

Figure 12 .
Figure 12.Comparison of the frequency deviations of microgrid 2 in the three scenarios.

Figure 13 .
Figure 13.Comparison of the frequency deviations of microgrid 3 in the three scenarios.

Figure 11 .
Figure 11.Comparison of the frequency deviations of microgrid 1 in the three scenarios.

Figure 11 .
Figure 11.Comparison of the frequency deviations of microgrid 1 in the three scenarios.

Figure 12 .
Figure 12.Comparison of the frequency deviations of microgrid 2 in the three scenarios.

Figure 13 .
Figure 13.Comparison of the frequency deviations of microgrid 3 in the three scenarios. .

Figure 12 .
Figure 12.Comparison of the frequency deviations of microgrid 2 in the three scenarios.

Figure 11 .
Figure 11.Comparison of the frequency deviations of microgrid 1 in the three scenarios.

Figure 12 .
Figure 12.Comparison of the frequency deviations of microgrid 2 in the three scenarios.

Figure 13 .
Figure 13.Comparison of the frequency deviations of microgrid 3 in the three scenarios. .

Figure 13 .
Figure 13.Comparison of the frequency deviations of microgrid 3 in the three scenarios.

Figure 14 .
Figure 14.Comparison of the exchange tie-line power deviations between microgrids 1 and 2 in the three scenarios.

Figure 15 .
Figure 15.Comparison of the exchange tie-line power deviations between microgrids 2 and 3 in the three scenarios.

Figure 16 .
Figure 16.Comparison of the exchange tie-line power deviations between microgrids 3 and 1 in the three scenarios.

Figure 14 . 19 Figure 14 .
Figure 14.Comparison of the exchange tie-line power deviations between microgrids 1 and 2 in the three scenarios.

Figure 15 .
Figure 15.Comparison of the exchange tie-line power deviations between microgrids 2 and 3 in the three scenarios.

Figure 16 .
Figure 16.Comparison of the exchange tie-line power deviations between microgrids 3 and 1 in the three scenarios.

Figure 15 .
Figure 15.Comparison of the exchange tie-line power deviations between microgrids 2 and 3 in the three scenarios.

Figure 14 .
Figure 14.Comparison of the exchange tie-line power deviations between microgrids 1 and 2 in the three scenarios.

Figure 15 .
Figure 15.Comparison of the exchange tie-line power deviations between microgrids 2 and 3 in the three scenarios.

Figure 16 .
Figure 16.Comparison of the exchange tie-line power deviations between microgrids 3 and 1 in the three scenarios.

Figure 16 .
Figure 16.Comparison of the exchange tie-line power deviations between microgrids 3 and 1 in the three scenarios.