Resiliency-Oriented Optimization of Critical Parameters in Multi Inverter-Fed Distributed Generation Systems

: In the modern power grid, with the growing penetration of renewable and distributed energy systems, the use of parallel inverters has signiﬁcantly increased. It is essential to achieve stable parallel operation and reasonable power sharing between these parallel inverters. Droop controllers are commonly used to control the power sharing between parallel inverters in an inverter-based microgrid. In this paper, a small signal model of droop controllers with secondary loop control and an internal model-based voltage and current controller is proposed to improve the stability, resiliency, and power sharing of inverter-based distributed generation systems. The distributed generation system’s nonlinear dynamic equations are derived by incorporating the appropriate and accurate models of the network, load, phase locked loop and ﬁlters. The obtained model is then trimmed and linearized around its operating point to ﬁnd the distributed generation system’s state space representation. Moreover, we optimize the critical control parameters of the model, which are found using eigenvalue analysis, and Grey Wolf optimization technique. Through time-domain simulations, we show that the proposed method improves the system’s resiliency, stability, and power sharing characteristics.


Introduction
With increasing penetration of distributed renewable energy resources in power system, the grid's stability and operation is becoming more important. Due to the distributed characteristic of these generations, Micro-Grids (MGs) are formed as decentralized groups of generation units and loads. MGs may work in grid-connected or islanded configurations to provide efficient, resilient, reliable, and clean energy [1,2]. Unlike traditional utility grids, which have large synchronous generators, MGs generally consist of low inertia renewable and distributed generation systems. Therefore, they are more prone to instability due to contingencies or extreme events. Consequently, operation and control issues are extremely important, especially in an islanded mode.
The inverters of Distributed Generations (DGs) are mostly connected in parallel to the grid to improve power sharing and stable performances [3]. Hence, designing a robust controller for Inverter-based Distributed Generations (IDGs) plays a critical role to improve the stability and resiliency of MGs.
Sensitivity and eigenvalue analysis are two common methods used along with the Small Signal Stability (SSS) analysis to improve the stability and resiliency of the MGs. One benefit of the SSS is that it is highly dependent on the accuracy of the mathematical models of the devices in the MG. Different DG systems are discussed in [4] and different types of the inverters are described based on the frequency regions. Furthermore, it is shown that the power controller coefficients are within the low frequency region, which is generally known as the dominant mode. Authors [5] present a SSS study to optimize the parameters of the inverter's PI current controller using Particles Swarm Optimization (PSO); however, in [5], the SSS analysis was only performed on the PI controller and mathematical model of the other components was not considered. Authors [6] utilized PSO to directly tune the control parameters without any SSS analysis. This can significantly increase the computation burden since there are no constraints for the parameters during the optimization process. In [7], the problems with [5,6] were addressed and all the MGs' controllers were optimized after identifying the significant control parameters using SSS.
Droop control is one of the most common approaches introduced for improving MGs' stability, resiliency, and power sharing between parallel inverters. In [8] a virtual flux droop is used to achieve voltage and frequency stability and equal power sharing. Reference [9] proposes a modified droop controller to improve the transient response and stability margins. Droop coefficients have a direct effect on the SSS analysis. It is shown that large droop coefficients deteriorate the stability and resiliency which can lead to high power output in some of the sources that can cause failures and even blackouts. On the other hand, small droop coefficients can result in unequal power sharing. Network parameters play an important role in power sharing. A very small mismatch between different DG's network parameters could lead to large differences in the shared power between the inverters. Therefore, sensitivity analysis is often used to determine the droop coefficients to improve the stability, resiliency, and power sharing of the IDGs [5,7,[9][10][11][12][13].
In [14][15][16][17], SSS analysis was investigated under different load conditions. Furthermore, in [4,7,[9][10][11][14][15][16][17], small signal analysis was studied in a system which uses PI controller to control the voltage and current in the dq axes. However, the dq axes are not fully decoupled and variations of one of the axes might lead to transients on the other one.
To have a more stable dynamic performance, Internal Model Control (IMC) was used in IDGs using the state space [18]; while in [3] the developed model was validated using SSS analysis. The results from [18,19] show that IMC significantly improves the stability of the IDGs in comparison to conventional PI voltage and current controllers. Nonetheless, in the previously mentioned literature, the models of the filter and the Phase Locked Loop (PLL) are not considered in the stability analysis. Authors [20] present a SSS analysis, considering the filter and PLL models and using a PI controller. However, oscillatory modes are neglected in [20] due to the large amount of the damping resistor.
Further, SSS analysis of IDGs using conventional droop controller is investigated in [5,[9][10][11][12][13][14][15][16][17][18][19][20][21]. Although the conventional droop controller is generally incorporated to increase the reliability and resiliency of the system, it has two major drawbacks. The first problem is voltage and frequency changes due to a sudden load variation. The second limitation is improper load sharing when the output impedances of the inverters don't match. To solve these issues, a robust controller and a universal droop controller are suggested in [3,22], respectively. Furthermore, [21,23] show a SSS analysis for the universal droop controller. However, they do not investigate the effects of the critical control parameters on the resiliency of the system. Authors in [24] develop a droop controller with secondary control loop to achieve equal power sharing between parallel inverters regardless of the output impedances of the inverters. However, the effects of the filter elements, PLL, and network and load models are not considered.
In the past decade, the resiliency of modern power system has attracted a lot of attention. Resiliency is defined as the ability of the power system to withstand extreme events and recover rapidly with the minimum damage. Reference [25] investigates the effects of power-electronic interfaces, energy storage systems (ESSs) and distribution system architectures on the resiliency of MGs and improves the resiliency during extreme events. Authors [26] propose an advanced control strategy for power electronic converters and suggest an energy management scheme for MGs with photovoltaic and battery ESSs. They use model predictive control strategy to minimize the operating cost, maximize the profit and obtain resiliency of the MG system. Furthermore, [27] proposes a two-stage, coordinated power sharing strategy between battery ESSs and coupled MGs for overload management in MGs, using dynamic frequency control.
Heuristic algorithms have received a lot of attention in the recent years due to their high speed and efficiency in solving optimization problems. One of the promising heuristic methods especially in the power system and power electronic applications, is Grey Wolf Optimizer (GWO) [10]. In [28], a chaotic GWO is proposed to minimize the power losses and improve the voltage stability in order to size and site DGs in the system. Moreover, the chaotic GWO determines the optimal droop parameters. Reference [29] uses GWO to optimize the control parameters of the proportional resonant controllers and find the optimal design of the output filter of a grid-tied inverter. Furthermore, Djerioui et al. in [30] present a predictive torque control of a permanent magnet synchronous motor based on GWO for smooth operation of electric vehicles. Specifically, the authors try to minimize the oscillations at low-speed operation of the motor. In [31], the optimal allocation of DGs is investigated to minimize the power losses while meeting the real and reactive demands in a distribution network. Furthermore, [32] proposes an optimization problem to minimize the total cost of the two-terminal HVDC system by incorporating the GWO in the optimal power flow algorithm. In this paper, we propose a novel methodology to improve the stability, resiliency, and power sharing of IDGs using GWO algorithm. The contributions of this paper are as follows:

•
Detailed mathematical modeling and equations of an IDG containing droop control with secondary control loop and an IMC-based voltage and current controller is presented.

•
The state space equations of the system are derived by linearization using Linear Analysis Tool in MATLAB/SIMULINK. The system matrix is then introduced. • Critical control parameters of the system are identified after the eigenvalue and sensitivity analysis on the system matrix. After finding the stability domains of the critical parameters, GWO is used to find the optimal critical parameters for the controllers. • The DG system's stability and power sharing characteristics are verified in time domain simulations by applying different load disturbances at different time steps. • Stability of the system is maintained by fast restoration and the impacts of different disturbances are reduced in order to improve the resiliency of the DG System. • Finally, it is shown that the resiliency of the DG system is significantly improved by comparing different stability indices.
The remainder of this paper is organized as follows. Section 2 provides a detailed description of the mathematical modeling and equations of the MG components. In Section 3, the linearization and state space equations of the MG system are discussed in detail. In Section 4, eigenvalue, and sensitivity analysis of the MG system are investigated. Furthermore, critical control parameters of each component are recognized, and the parameter stability domains are presented. After which, the optimal control parameters are derived using GWO. Section 5 presents the simulation results and performance evaluation of the optimal parameters. Finally, conclusions and discussions of the future work are presented in Section 6.

Mathematical Modeling
In this work, we consider a microgrid with two distributed generation systems connected through a power line with impedance Z. The system is shown in Figure 1. Each DG is connected to a static load. Figure 2 shows the DGs' system block diagram including the LCL filter, voltage and current controller, PLL and power controller. The mathematical modelling of each component along with the dynamics of the loads and the network are discussed in the following sections. The equations are linearized around the operating point and the state space model of the system is derived later.

Droop Control with Secondary Control Loop
In this work, the power controller is a droop control with a secondary control loop. A low pass filter is used to remove the higher harmonics of the instantaneous powers P and Q . These instantaneous powers are calculated using inverter's output voltages and currents (v , i , v , i ) as expressed in (1) and (2). Figure 3 shows the secondary control loop design of the power controller used in the suggested MG [5]. The secondary control loop collects the power outputs P and Q of each inverter and then calculates the reference active/reactive powers according to the inverter's power ratings. In this work, the power ratings of the inverters are assumed to be equal (P = P ) and the studied MG is considered as a low voltage MG. In the low voltage MGs, the output reactance is lower than the output resistance of the inverter. Therefore, a

Droop Control with Secondary Control Loop
In this work, the power controller is a droop control with a secondary control loop. A low pass filter is used to remove the higher harmonics of the instantaneous powers P and Q . These instantaneous powers are calculated using inverter's output voltages and currents (v , i , v , i ) as expressed in (1) and (2). Figure 3 shows the secondary control loop design of the power controller used in the suggested MG [5]. The secondary control loop collects the power outputs P and Q of each inverter and then calculates the reference active/reactive powers according to the inverter's power ratings. In this work, the power ratings of the inverters are assumed to be equal (P = P ) and the studied MG is considered as a low voltage MG. In the low voltage MGs, the output reactance is lower than the output resistance of the inverter. Therefore, a

Droop Control with Secondary Control Loop
In this work, the power controller is a droop control with a secondary control loop. A low pass filter is used to remove the higher harmonics of the instantaneous powers P and Q. These instantaneous powers are calculated using inverter's output voltages and currents (v od , i id , v oq , i iq ) as expressed in (1) and (2). Figure 3 shows the secondary control loop design of the power controller used in the suggested MG [5]. The secondary control loop collects the power outputs P and Q of each inverter and then calculates the reference active/reactive powers according to the inverter's power ratings. In this work, the power ratings of the inverters are assumed to be equal (P 1r = P 2r ) and the studied MG is considered as a low voltage MG. In the low voltage MGs, the output reactance is lower than the output resistance of the inverter. Therefore, a P-f/Q-V droop can be well-implemented, which is commonly used for low voltage MGs. The block diagram of the droop controller is presented in Figure 4 and the corresponding equations are shown in (3) and (4) for each of the DGs: (i = 1, 2) Sustainability 2021, 13, x FOR PEER REVIEW 5 of 17 P-f/Q-V droop can be well-implemented, which is commonly used for low voltage MGs. The block diagram of the droop controller is presented in Figure 4 and the corresponding equations are shown in (3) and (4) for each of the DGs: (i = 1, 2)

Phase Locked Loop (PLL)
PLL is used to model the effects of load changes on the system's frequency. Figure 5 shows the PLL used for the small signal analysis of the DG system. This PLL is modeled in the dq-axex reference frame [7]. The input of the model is the d-axis voltage of the capacitor, C . The PLL tracks phase angle using a PI controller. The input of the PLL is the d-axis of the measured output voltage, and it is set to zero (v = 0). Furthermore, the voltage of the load is considered to be in q-axis (v ). Following are the equations related to the design of the corresponding PLL: P-f/Q-V droop can be well-implemented, which is commonly used for low voltage MGs. The block diagram of the droop controller is presented in Figure 4 and the corresponding equations are shown in (3) and (4) for each of the DGs: (i = 1, 2)

Phase Locked Loop (PLL)
PLL is used to model the effects of load changes on the system's frequency. Figure 5 shows the PLL used for the small signal analysis of the DG system. This PLL is modeled in the dq-axex reference frame [7]. The input of the model is the d-axis voltage of the capacitor, C . The PLL tracks phase angle using a PI controller. The input of the PLL is the d-axis of the measured output voltage, and it is set to zero (v = 0). Furthermore, the voltage of the load is considered to be in q-axis (v ). Following are the equations related to the design of the corresponding PLL:

Phase Locked Loop (PLL)
PLL is used to model the effects of load changes on the system's frequency. Figure 5 shows the PLL used for the small signal analysis of the DG system. This PLL is modeled in the dq-axex reference frame [7]. The input of the model is the d-axis voltage of the capacitor, C f . The PLL tracks phase angle using a PI controller. The input of the PLL is the d-axis of the measured output voltage, and it is set to zero (v od = 0). Furthermore, the voltage of the load is considered to be in q-axis (v oq ). Following are the equations related to the design of the corresponding PLL:   Figure 6 shows the IMC-based voltage and current controller which uses PI, PD, and PID controllers to control the output voltage and current of the filter along with the frequency of the system [3]. By using IMC, the system's output can be predicted through different control signals in the feed-forward loop. Another advantage of the IMC is that the set points can be changed using a feedback loop to balance out the different model mismatch and disturbances.   Figure 6 shows the IMC-based voltage and current controller which uses PI, PD, and PID controllers to control the output voltage and current of the filter along with the frequency of the system [3]. By using IMC, the system's output can be predicted through different control signals in the feed-forward loop. Another advantage of the IMC is that the set points can be changed using a feedback loop to balance out the different model mismatch and disturbances.   Figure 6 shows the IMC-based voltage and current controller which uses PI, PD, and PID controllers to control the output voltage and current of the filter along with the frequency of the system [3]. By using IMC, the system's output can be predicted through different control signals in the feed-forward loop. Another advantage of the IMC is that the set points can be changed using a feedback loop to balance out the different model mismatch and disturbances.

Voltage Controller
Mathematical modeling of the internal model-based voltage controller with filter dynamics is derived from [3]. The control parameters are described as below: k pv , k dv , k pv , and k iv stand for the coefficient of the first proportional controller, coefficient of the derivative controller, coefficient of the second proportional controller, and coefficient of the integral controller in the IMC-based voltage controller respectively.
As shown in Figure 6, the reference frequency and the output voltage from the power controller are used as set-points for the voltage controller. The difference between droop frequency and the voltage from the load and their set-points are calculated and sent to the IMC voltage controller. φ d and φ q are the states of the IMC voltage controller. (10)- (13) show the equations related to the voltage controller.

Current Controller
IMC current controller is designed and implemented using the same methods presented in [2,3]. As it is seen from Figure 6, the filtered currents (i ld , i lq ) are compared with the voltage controller output (i * ld , i * lq ) and then the error is sent to the current controller. The output of the integrator, γ d and γ q , are considered as the state variables. The corresponding differential equations and the control parameters of this controller are shown in (14)- (18). k pc , k ic , k dc k pc , and k ic show for the coefficient of the first proportional controller, coefficient of the derivative controller, coefficient of the integral controller, coefficient of the derivative controller, coefficient of the second proportional controller, and coefficient of the second integral controller in the IMC-based current controller respectively.

LCL Filter
The LCL filter is shown in Figure 2. This filter is used to eliminate the higher harmonics. r c and r f are parasitic resistances of the inductor. The resistance of the capacitor is also

Network and Load Model
As is seen in Figure 1, Z represents the line impedance between the two microgrids. This impedance has a direct effect on the power sharing between the two inverters. The transmission line is mathematically modeled as below: A static load including a combination of resistors and inductors is connected to each bus of the system. They are modeled using the following equations. For i = 1, 2 we have: Because the system frequency is constant, ω pll appears in the network and load equations. In the equations, the variables which have DQ as the subscript, are in the global reference frame. The first DG system is considered as the reference frame for the whole study. The new phase angles for the DGs considering the reference frame, are presented in (24).

Microgrid Model
The steady state operating point of the system is found by linearizing the nonlinear equations using linear analysis tool in MATLAB/SIMULINK. The DG system model is then linearized around the operating point and the state space representation is obtained as (25).
In this equation, u represents the input and consists of the following: In this microgrid every DG has 15 states. The load model has 4 states, and the network model has 2 states. Hence, the state space representation has 36 states in total. Matrix x in (25) represents the state variables. The bus voltages should also be considered in the system matrix (A) to accurately predict the load changes. In this respect, a large virtual resistor (R N ) is considered to have minimum effects on the dynamics of the system. By introducing R N , the input matrix u is determined in terms of the states using (28).
As previously discussed, DG1 is assumed as the reference. It is important to transform all the variables from their local reference to a global reference frame, because all the mathematical models were derived in each of the component's local frame. This transformation is carried out using (29). (29) δ represents the phase angle between the global and the local reference frame. This equation is also linearized and incorporated to the system matrix.

Control Parameter Optimization
The eigenvalues of matrix A in the state space representation shown in (25), are the system's poles. By tuning the control parameters of the system, the critical control parameters are found. As it is shown in Figure 7, the eigenvalue spectrum of matrix A can be divided into three frequency regions: high, medium, and low frequency. It is known that the control parameters related to the power controller are in the low frequency region. These parameters are also generally known as dominant mode. It is also demonstrated that the critical parameters concerning the IMC-based voltage and current controllers are in the medium and high frequency regions.
The values of λ v and λ c are estimated and then tuned (λ v = 10 −4 λ c = 10 −5 ). In Table 1, the control parameters, and parameters of the microgrid are shown. In mboxfigfig:sustainability-1219784-f008 the eigenvalue analysis is presented by changing the m q parameter of the droop controller. As the arrows show in Figure 8, by increasing the m q parameter from 10 −6 to 10 −2 , the poles of the DG system go more towards zero and also in some cases go on the right side of the real axis. This means that increasing the m q parameter will lead to more instability in the DG system. Therefore, we have to first consider a domain for m q to ensure the stability of the DG system and then aim to find an optimal value for that parameter. By repeating the same analysis for other control parameters, the domain of the critical parameters are obtained as follows: n p ∈ 10 −5 , 10 −3 ; m q ∈ 10 −5 , 10 −2 ; K e ∈ [ 1, 120]; k pv ∈ 6.72 × 10 −3 , 0.25 ; k iv ∈ [ 1, 105]; k pv ∈ [0.01, 12]; k pc ∈ [ 1272]. These domains will help to improve computational efficiency during the optimization process. It is noteworthy to mention that the other parameters do not have significant effects on the power sharing, stability, and resiliency of the system. Sustainability 2021, 13, x FOR PEER REVIEW 10 of 17 The values of λ and λ are estimated and then tuned (λ = 10 λ = 10 ). In Table 1, the control parameters, and parameters of the microgrid are shown. In Figure 8 the eigenvalue analysis is presented by changing the m parameter of the droop controller. As the arrows show in Figure 8, by increasing the m parameter from 10 to 10 , the poles of the DG system go more towards zero and also in some cases go on the right side of the real axis. This means that increasing the m parameter will lead to more instability in the DG system. Therefore, we have to first consider a domain for m to ensure the stability of the DG system and then aim to find an optimal value for that parameter. By repeating the same analysis for other control parameters, the domain of the critical parameters are obtained as follows: n ∈ [10 , 10 ]; m ∈ [10 , 10 ]; K ∈ [1, 120]; k ∈ [6.72 × 10 , 0.25]; k ∈ [1, 105]; k ∈ [0.01, 12]; k ∈ [1272]. These domains will help to improve computational efficiency during the optimization process. It is noteworthy to mention that the other parameters do not have significant effects on the power sharing, stability, and resiliency of the system.   In this paper, Grey Wolf Optimizer [10] is used to optimize the parameters of the controllers. This algorithm is inspired by the hierarchy and hunting methods of the grey wolves. The advantages of this method over other heuristic methods are that it is flexible and easy to implement, and it is applicable to different problems without any significant changes. Moreover, this method has a better exploration ability in comparison with other methods, especially in the power electronics applications.
Power sharing is an important aspect of the operation of parallel inverters. Therefore, in this problem, the objective is considered as a deviation from a reference reactive power (Q ref ). For the optimization, the objective presented in (30) is considered with the constraints shown in (31). A population size of 15 and an iteration number of 50 are assumed for the optimization. Figure 9 presents the objective function convergence with respect to the iterations. methods, especially in the power electronics applications.
Power sharing is an important aspect of the operation of parallel inverters. Therefore, in this problem, the objective is considered as a deviation from a reference reactive power (Q ). For the optimization, the objective presented in (30) is considered with the constraints shown in (31). A population size of 15 and an iteration number of 50 are assumed for the optimization. Figure 9 presents the objective function convergence with respect to the iterations.   methods, especially in the power electronics applications. Power sharing is an important aspect of the operation of parallel inverters. Therefore, in this problem, the objective is considered as a deviation from a reference reactive power ( Q ). For the optimization, the objective presented in (30) is considered with the constraints shown in (31). A population size of 15 and an iteration number of 50 are assumed for the optimization. Figure 9 presents the objective function convergence with respect to the iterations.  The derived optimal values for the critical parameters n p , m q , k e , k pv , k iv , k pv , k pc are: 0.008, 9.018 × 10 −5 , 120, 0.133, 1.320, 20.240 and 283.6611, respectively. The derived optimal values for the critical parameters n p , m q , K e , k pv , k iv , k pv , k pc are: 0.008, 9.018 × 10 −5 , 120, 0.133, 1.320, 20.240 and 283.6611, respectively.

Simulation Results
An active load of 5 kW is added to DG1 between 1 s and 1.5 s and another active load of 10 kw is connected to DG2 between 1.5 s and 2 s. Figure 10 shows the active and reactive power sharing between the two DGs using the proposed critical parameters in time-domain simulation. Equal power sharing between the DGs exists in steady-state condition. It is shown that from 1 to 1.5 s, DG1 is sharing more power due to the increase of load 1. Similarly, DG2 is sharing more power from 1.5 to 2 s when load 2 is increased. The stability of the frequency and the q-axis voltage of the DGs are also presented in Figure 10a,b, respectively. It is shown that both voltage and frequency are regulated within the acceptable range. Figure 10 also shows that the stability of voltage and frequency is maintained by fast restoration, thus proving the improvement of the resiliency of the system. In order to evaluate the performance of the proposed critical parameters of the system, two disturbances (a small and a large one) are applied. The results using the obtained critical parameters are shown in Table 2. The DGs' stability and resiliency performance are reported as rising time (t ), overshoot (M ) and settling time (t ). For comparison, the same analysis is executed for the control parameters derived from the conventional PIbased droop control method, which is implemented in [4] and primarily used in realworld applications. The results from the conventional method are presented in Table 3. For instance, the overshoot (M ) for both DGs in case of small disturbance in Table 2 are much less than the overshoot of the DGs with small disturbance in Table 3. We can see the same improvement for the t and t indices in Table 2 in comparison with Table 3. Hence, the effectiveness of the proposed parameters identification method is verified. Especially less settling time (t ) in the critical parameters case shows fast restoration after applying the disturbances, and less overshoot (M ) leads to less chance of damaging the equipment in the DG system while restoration. Consequently, the results show the improvement of the resiliency as well. In order to evaluate the performance of the proposed critical parameters of the system, two disturbances (a small and a large one) are applied. The results using the obtained critical parameters are shown in Table 2. The DGs' stability and resiliency performance are reported as rising time (t r ), overshoot (M p ) and settling time (t s ). For comparison, the same analysis is executed for the control parameters derived from the conventional PI-based droop control method, which is implemented in [4] and primarily used in real-world applications. The results from the conventional method are presented in Table 3. For instance, the overshoot (M p ) for both DGs in case of small disturbance in Table 2 are much less than the overshoot of the DGs with small disturbance in Table 3. We can see the same improvement for the t r and t s indices in Table 2 in comparison with Table 3. Hence, the effectiveness of the proposed parameters identification method is verified. Especially less settling time (t s ) in the critical parameters case shows fast restoration after applying the disturbances, and less overshoot (M p ) leads to less chance of damaging the equipment in the DG system while restoration. Consequently, the results show the improvement of the resiliency as well.

Conclusions and Future Work
In this work, the stability, resiliency, and power sharing characteristics of inverterbased distributed generation systems comprising of a droop control with a secondary control loop and an internal model-based voltage and current controllers are studied. Accurate mathematical modeling and equations of each of the mentioned components along with phase locked loop, and the network and load models are discussed. The equations are then modeled using MATLAB/SIMULINK. System Operating point is derived, and the equations are then linearized around this operating point using the Linear Analysis Tool. The state space equations of the system are presented, and the eigenvalue and sensitivity analysis of the inverter-based distributed generation systems are presented using the system matrix (A). The critical parameters are identified, and stability domains of these parameters are presented. Moreover, a grey wolf optimization algorithm is utilized to find the optimal control parameters and to improve the stability, resiliency, and power sharing in the system. Various load disturbances are applied in different time steps to show the effectiveness of the proposed algorithm in the time domain simulations. Finally, different stability indices are calculated to verify the stability and resiliency of the distributed generation systems. The results show that the proposed methodology improves the system's performance compared with the conventional droop controller.
In our future work, the performance of the proposed inverter-based distributed generation system will be compared with other control schemes. Furthermore, different optimization techniques will be investigated.

Data Availability Statement:
The data that support the findings of this study are available from the corresponding author upon reasonable request.

Conflicts of Interest:
The authors declare no conflict of interest. The coefficient for the proportional controller related to the voltage controller K iv

Nomenclature
The coefficient for the integral controller related to the voltage controller K dv The coefficient for the derivative controller related to the voltage controller K pc The coefficient for the proportional controller related to the current controller K ic The coefficient for the integral controller related to the current controller K dc The coefficient for the derivative controller related to the current controller K pd (s) Transfer function of the PD controller in the internal model-based controller K pi (s) Transfer function of the PI controller in the internal model-based controller K pid (s) Transfer function of the PID controller in the internal model-based controller λ v Tuning parameter for the internal model-based voltage controller λ c Tuning parameter for the internal model-based current controller T c Time period of the inverter output signal T pwm Time period of the inverter switching φ i States of the internal model-based voltage controller γ i States of the internal model-based current controller V bDi Input voltage in d-axis V bQi Input voltage in q-axis R N Virtual resistor m q A control parameter of the power controller used to control the voltage n p A control parameter of the power controller used to control the frequency K e One of the Control parameters of the power controller f sw Inverter's switching frequency