Resilience Improvement of Distribution Networks Using a Two-Stage Stochastic Multi-Objective Programming via Microgrids Optimal Performance

Dealing with major power disruption during natural disasters is one of the most notable concerns in power systems. In this regard, the optimal application of microgrids as a potential solution in increasing and sustaining the distribution system resilience is considered. Furthermore, this purpose is pursued while maintaining the resilience of each DC microgrid connected to the distribution system, which is an essential and challenging issue. In the proposed method of this paper, a novel modeling strategy is formulated as a multi-period two-stage scenario-based stochastic mixed-integer linear programming (MPTSS-MILP) based on a multi-objective optimization problem (MOOP). In this framework, the operation associated with emergency and normal conditions, according to the influences of each situation on another one, is managed in multi-microgrids coordinately. In this regard, the technical constraints correlated to the operation of microgrids as well as the distribution system are satisfied simultaneously in specific to each condition which covers normal and critical operating under all uncertainty scenarios. Through introducing two evaluation criteria and also a resilience metric in microgrids and distribution systems, the efficiency of the proposed method is demonstrated. Meanwhile, innovative modeling is executed based on the subjective behavior of people affected by the disaster. The proposed method is implemented on a test system that involves a 34-bus distribution system with three distinct DC microgrids. In this regard, the impact of plug-in electric vehicles, as well as social behavior affected by severe events, demonstrate significant results according to the resilience criteria and resilience metric of microgrids and distribution system in three case studies based on the proposed approach.

INDEX TERMS Resilience, extreme events, microgrids, stochastic optimization, distribution system.  Severe natural and man-made disasters have increased over the past few decades. These events directly affect the electricity infrastructure, resulting in a lot of economic damage. Extreme weather and natural disasters (e.g., earthquakes, floods, and hurricanes) lead to economic, social and physical disruptions and cause significant inconvenience to people living in the affected areas [1], [2]. Power systems are established by the main aspects of reliability, so insured against events with low-effect but a high probability of occurrence. In the meantime, a new topic called resilience in power systems is raised, which represents the resistance of electrical infrastructure to incidents with high impact but low probability (HILP) events [3], [4]. The potential of an energy system to withstand turbulence and provide affordable energy services to consumers is a definition of resilience in reference [5]. Coming out of the shock as soon as possible and offering alternative options to meet the needs of energy services are features of a resilient system. The classification of strategies to increase the resilience of the electrical system is divided into methods based on planning and based on operations. In [6], [7] and [2] can be mentioned from the studies conducted in the field of increasing resilience of power system based on planning. For example, In [8], in the planning field, the authors presented a method for distribution system resilience enhancement via line hardening and DGs installation. To reduce the impact of disasters, operation-based resilience methods propose immediate solutions. In order to increase resilience in the field of the distribution network operation, issues based on microgrid (MG) formation and division of electricity networks into smaller MGs have been studied in [9], [10]. In contrast, some studies have focused on the integration of multiple MGs intending to increase resilience [11], [9]. In [12] state-ofthe-art methodologies and performance of networked microgrids (NMGs) have been investigated. In the operation-base method with enhancing power system resilience approach, various researches have been done. Distribution network reconfiguration-based methods to increase the resilience of the power grid are also among the problems that have been investigated extensively [4], [13]. In [14], the authors investigate a practical case study focusing on the following aims simultaneously against physical attacks: allocating the feeder routing problem and substation facilities, obtaining the models of conductor's installation, and the power lines hardening based on economic point of view. This study has been solved based on meta-heuristic optimization methods by considering grey wolf optimization (GWO) algorithm according to normal and emergency conditions, and total blackout restoration or local outage has been resulted. Reference [15] focuses specifically on the storm occurrence and dealing with its destructive effects on the distribution system. Therefore, three-stage framework includes disaster occurrence modeling, component vulnerability, and network resilience is proposed. Some hardening strategies such as vegetation management and upgrading of network poles are applied to distribution system resilience enhancement. Scheduling of energy storage units [16], [17] and use of mobile emergency resources [18] are the solutions to increase resilience in the distribution networks operation. Load restoration-based methods are another purpose of increasing resilience that has been used in [19]. MGs with DG are introduced as a viable solution in [20], [21] for local load survivability during contingencies due to their islanding schemes. In this regard, according to a hybrid AC/DC MG resilience enhancement, a resilience-driven approach is presented in emergency and normal operation mode. Among the most common methods mentioned, the use of MGs has become popular as a suitable source to increase resilience. MGs with distributed generation energy resources (DERS) in emergencies are a proper solution to increase resilience. Because they can join MGs or networks that do not have enough capacity to maintain their sensitive loads in critical conditions [9], [11], [22]. The authors in [23] believe that a RES-based MG is useful for the distribution system resilience. Reducing carbon emission is also mentioned as another benefit of using mentioned MGs. Furthermore, different categories of MG architecture besides the technical characteristics of energy storage devices and various communication channels are investigated in [23]. Among the applications of MGs, which are divided into four categories in reference [11]; we can refer to MGs that are used as a source to increase global resilience. Due to the fact that traditional networks such as conventional distribution networks do not have enough resources to maintain their sensitive loads in case of emergency, MGs with different resources can help such networks in these situations. By combining MGs with each other and also with the distribution feeder in critical situations, the resilience of the power system can be increased. Additionally, the integration of MGs in which the penetration of renewable energy is high and are also based on power electronics increases the inertia of their collection in emergencies. Because frequency deviation due to fluctuations in renewable energy resources is better managed and controlled [24]. To achieve this goal, both in the MG set and distribution network, the objective functions and constraints must be defined and considered under different operating conditions (emergency and normal). These should apply to both MGs and distribution feeders; Which has not been considered in this way in any of the previous references. The study on small-size decentralized MGs scheduling under normal and critical operating conditions is developed in [25]. In this paper, based on MGs' resilience performance an energy management system (EMS) is suggested. Consequently, the managing and coordination of multi-microgrids (MMGs) and distribution network in different operating situations is a very important and significant point. Due to the inherent uncertainty of resilience, in this situation, one of the main challenges is uncertainty about the duration of disconnection from the main network [26]. Uncertainty about the start of an unscheduled disaster as well as the duration of the emergency condition are important parameters in resilience programming, which is considered only in [16]. Because in every hour of the time horizon, the system conditions are different due to the high penetration of renewable energy sources (RESs). Considering this issue proves the efficiency of the proposed methodology and model. On the other hand, MGs management strategies to increase resilience include central, decentralized, and hierarchical methods [10], used in [27] and [28] as an outage management strategy, but in modeling the static behavior of the power system, the constraints related to reactive power and its equilibrium equations as well as the constraints related to voltage limits have been neglected.

Indices & sets
The presence of electric vehicles (EV) in resilience studies can be observed in references [29], [30] and [31]. In [29], generation resources, battery energy storage systems (BESS), and EVs of MGs are programmed to be able to feed the most important MG loads in critical situations. But this scheduling is proposed for all EVs in the same approach and at the same time. In a MG that contains several EVs, not all EVs can be charged or discharged at the same time and following the same conditions and it does not resemble to be an appropriate solution. In [30], authors use BESS and EVs with a different approach to reducing load shedding amount in emergency states. But here, as in previous research, specific scheduling for each EV is not considered. In [31], in an energy management system to increase the resilience of MGs, a demand response program is modeled through the PHEV fleet. In this case, it is assumed that the PHEV fleet is available at different times of the day, but the driving pattern of electric vehicles and the times when they are not available in the parking lot are not considered.
Given that in emergencies and severe disasters, the people of that geographical area are additionally affected, modeling the social behavior of the community is effective in scheduling the power system, which has not been considered in any of the previous studies. The innovative aspects employed in this study in comparison with other recent studies are conferred in Table 1.
Based on mentioned points, in this paper, we seek to provide a way to manage MMGs under all conditions in a system to improve the resilience of the distribution system connected to them so that their resilience is maintained at an acceptable level. According to Figure 1, one of the objectives of this modeling is to have an economically efficient operation in normal system conditions. Although in emergency conditions, our objectives are focused on minimizing the energy curtailment of the distribution network and minimizing the load shedding of the MMGs, simultaneously. Since MMGs  are responsible to provide loads in critical situation, optimal operation in economic perspective is considered in second priority of objective function. In the following, the applied method is introduced in the paper as a resilience evaluation criterion.
It should be noted that in the real system, MGs have different dimensions with each other in terms of production and consumption resources. In this paper, following the mentioned point, the MGs in the studied system are completely different from each other, which poses a serious challenge in the planning of these resources. The method of coordination and preparation of MGs in each of the conditions in the system, according to the technical constraints is one of the prominent features of this research. For example, one of the most essential limitations is power flow equations which the DistFlow model is employed as the power flow model for the distribution system. All this has been done to improve and maintain the resilience merci of each MG and distribution network under uncertainty conditions through a stochastic approach. In fact, a two-step uncertainty is considered in this study. In the first step, the event occurrence time has uncertainty with a uniform distribution. In this probability distribution function, the probability of each event occurring is considered the same. Therefore, at every hour of the day-ahead scheduling time, the probability of a severe event is the same. The second step of uncertainty is the duration of the emergency condition. Due to the mentioned uncertainties, the modeling is based on a stochastic optimization method. One of the efficient methods for conducting stochastic problems is the scenario-based method that is utilized in this study This uncertainty is modeled by the stochastic scenario-based programming (SSP) method, which is a mixed-integer linear programming (MILP). Considering the importance of the reaction and social behavior of people in the area affected by the event, this behavior, which is reflected in the use of plug-in electric vehicles (PEV), has been modeled in this issue and its effects on the results have been studied. In scheduling a system that interacts closely with electricity consumers, their behavior and reaction in different situations are very important. Because in the programming, such reactions should be modeled and considered so that the proposed method is as close to reality as possible. Based on mentioned aims, the main contributions and assumptions of this paper can be listed as follows: • Proposing a new approach to achieving the goals related to MGs and distribution system in all normal and emergency operation according to the specific constraints of each operation condition under uncertain states with a multi-period two-stage scenario-based stochastic mixed-integer linear programming (MPTSS-MILP) based on a multi-objective optimization problem (MOOP).
• A novel scheduling based on two-stage multi-period programming to ensure adjustability and adaptivity in each scenario of the proposed approach is applied.
• A more comprehensive consideration in the interaction and coordination between several completely different MGs in terms of power generation sources, load profiles, electrical energy storage parameters, number of PEVs, and driving pattern to examine the impact and cooperation of each of them in different situations.
• Considering the different climate and environmental conditions for each MG, depending on their geographical location, alongside with considering various RESs, wind and solar as clean productions, that their generation at any hour depends on climate circumstances.
• Modeling a new and innovative assumption based on the reaction of people affected by the disaster in that area on using PEVs and the effect of this behavior on MGs and distribution system resilience.
An overview of the problem and the studies conducted with the mentioned issue in comparison with this paper is illustrated in Figure 2.
The remainder of the current paper is organized as follows ( Figure 3). In Section 2, the fundamental formulation, concepts and methodology associated with the proposed mathematical model have been thoroughly discussed. Section 3 particularly emphasizes the formulation including the objective functions and the constraints correlated with them according to the defined problem. Consequently, in Section 4, numerical results derived from the simulation of the suggested approach will be declared. Subsequently, the conclusions are outlined in Section 5.

II. MODELING AND MANAGING APPROACH BASED ON RESILIENCE A. MODELING AND CONCEPTS
According to the features of MGs in disasters, three DC MGs are used in this issue, which is explained particularly in IV. Resources in these MGs include RESs, fossilfuel-based sources, energy storage systems (ESS), three-level priority loads, and PEVs. It is noteworthy that these MGs can exchange power with the upstream network according to the market price in normal situation. In addition, the amount of power transfer at the PCC point of MGs, as well as the total power injection of the upstream network to the studied system, is technically limited due to power flow distribution lines limitation. Further, weather conditions in each MG are considered according to the geographic area located in it. This is due to the high dependence of sustainable energy sources on climatic conditions. Furthermore, each of these MGs is connected to one of the buses of the distribution network and exchanges power with the upstream network through the distribution feeder. In the normal system operation, the aim is to schedule the resources of MGs optimally so that all the loads in them are supplied at the minimum cost. Now, if the system enters a critical operation mode, new conditions will prevail in the system. Therefore, the distribution system, along with the three MGs, are generally disconnected from the upstream network. Each MG in this situation is programmed to minimize the load shedding so that sensitive loads are necessarily supplied. In this regard, each MG can transfer its excess power to other MGs if possible due to technical and operational constraints.
In fact, with integrated management of MGs collection, it can be shown their cooperation with each other in critical situations. On the other hand, load supplying is more superior to minimizing operation costs in emergencies and is a priority for the objective functions. Meanwhile, the distribution feeder which has sensitive loads, attempts to provide its maximum priority loads. While the distribution network has no power generation source for its own. This assumption is intended to take into account the worst-case scenario for the system to ensure the efficiency of the proposed method. An important point in dealing with this problem is the combination of critical and emergency operating conditions in the system. According to the descriptions, we are confronted with a multi-objective optimization problem (MOPs), in a situation of emergency operation. These goals are pursued as resilience evaluation criteria in the distribution system and MGs. In these optimization problems, two or more objective functions are optimized simultaneously and are always in conflict with each other. An appropriate answer to one objective function may be an unsuitable solution for the other one. Therefore, for a multi-objective problem, it is difficult to find a solution that satisfies all the objective functions. Accordingly, MOP can be represented as the following general mathematical model [34], [35] and [36]: where k is assumed to be an optimized number of objective functions, θ and ϕ are a C C C-dimensional search space for the decision variable and k-dimensional vector space of objective functions respectively. Inequality and equality constraints are denoted by the equations g m ( ) and h n ( ) with the number of M M M and N N N, respectively. Analytical and numerical methods are generally two kinds of techniques to solve MOPs. In this manuscript, the analytical method is used because it involves accurate mathematical proof and can provide an exact solution. In this regard, there are two general approaches to deal with such issues: priori method and posteriori method [34] and [35]. According to the first approach, all objective functions are combined and become a single objective. One of these ways is to use the Weighted Sum Method (WSM) of objectives according to the importance of each of them.
where, weighted coefficients set is determined with w q . Since in the set of MGs, the precedence of the minimum load shedding is more than the operation cost, the WSM method has been used to achieve these goals. As a result, the objectives of MGs are aggregated into one objective VOLUME 9, 2021  according to their importance. On the other hand, the distribution network also seeks to achieve minimizing energy curtailment. For this purpose, the loads due to their importance, using weighted coefficients (3), are present in the objective function of load loss minimization. Now, given the conflict between the goals of MGs and distribution networks, a solution must be determined to solve this problem [37]. In contrast to the first method described, the latter approach, namely the posteriori method solves the MOP problem through the concept of non-dominance. In this method, several optimal points, which are the set of optimal Pareto-front solutions, are obtained that is displayed Figure 4.
The most appropriate solution from this set is selected based on the preferences of the determiner. To solve this part of the problem, the ε-Constraint method (ECM) has been used. This method is in the second category of the discussed methodologies by producing the Pareto front of optimal answers. In this method, one of the objective functions is selected as the optimized objective and the other ones are selected as the constraint for this optimization problem [34] and [35].
The set of upper limits of the objective functions in the constraints are denoted by the values of ε = [ε 1 , ε 2 , . . . , ε k ] T . The method used to obtain epsilon values in this manuscript is illustrated in Figure 5.
Due to the uncertainty related to the emergency condition duration and the occurrence time, a stochastic problem must be solved. To achieve the objectives under uncertainty, the two-stage scenario-based stochastic programming (TSSP) method is used. The general structure of SSP is as follows [38]: In the general form of SSP, there are two kinds of decision variables. The first category is Here-and-Now (HN) or First-Stage variables. Decisions about the values of these variables must be determined before applying uncertainty scenarios. In contrast, after determining the uncertainty scenario, some variables are decided, which are named Second-Stage or Wait-and-See (WS) variables. Equation (6) consists of two terms; the first part is related to the HN variables and the second part is related to the WS variables. Constraints (7) and (8) are related to the first-stage and second-stage variables, respectively. The important point of TSSP is the existence of a relationship between these two kinds of variables with each other in a constraint that is mentioned in equation (9). In this study, we have solved the problem using the classic form of the SSP model. By replacing the expression for the classical form instead of φ(H (x ω )), it becomes as follows [38]: where, π ω is the probability of any scenario's occurrence. Each of these scenarios is characterized according to the probability distribution function of an uncertainty parameter. So, in the MOP model, the objective functions are represented as (10). Using a combination of the above methods, a hybrid method for satisfying the objectives of MGs and distribution network is presented in Figure 6. In this case, two parameters, the starting time of a critical situation and the duration of its continuance, are uncertain. Given that the start time of the disaster is unknown, it is assumed that the probability of unscheduled occurrence in each study period is the same. Further, at any stochastic unscheduled event start time, there is another uncertainty regarding the duration of the emergency operation ( Figure 7).
The current circumstances in normal and critical operating conditions affect each other and this issue has been completely considered in this study. These requirements include technical constraints related to the elements in the system.   After obtaining the Pareto-front set, we must find a solution that creates a compromise between the two objective functions using a suitable method. There are many techniques to determine a compromise solution between Pareto front set of answers. The fuzzy satisfaction method or max (min (.)) method is a conventional way for selecting the best solution among Pareto optimal answers. In this method, if it is assumed that N objective functions exist to be minimized, the linear membership function for then th answer of each objective is defined as follows [39], [40]: fn k ≥ f max k k = 1, . . . , N ,n = 1, . . . , N P (11) Here f max k and , f min k are the maximum and minimum values of the k th objective function in the set of Pareto optimal solutions, respectively. µn k presents the degree of optimality of then th solution of the k th objective function. The membership function of solutionn can be calculated using the following equation: A conservative decision maker strives to maximize satisfaction among all objective functions. The final solution can be found as follows: Using this process, the optimal compromise solution between the objective functions of the MGs and the distribution network is obtained. At this point, the variables in the problem are in their optimal values. The output variables of the problem in both critical and normal operation modes are shown in Figure 8.

B. NEW APPROACH ON MODELING SOCIAL BEHAVIOR IN CRITICAL SITUATIONS
When a severe event occurs, people in that region are also psychologically affected. This event itself changes the behavior of people in society. One of these variations is the change in energy consumption. In this paper, a variation in the pattern of using electric vehicles by their owners at the time of the disaster is considered. This behavior pattern is observed immediately upon the start of an emergency by PEV owners. After the completion of these conditions, this behavior returns to the normal (Figure 9). The effect of this behavior on the results is also examined to ensure the efficiency of the proposed approach.
To more clarifying the paper's contributions advantages compared to other methods are illustrated in Figure 10.

III. PROBLEM FORMULATION
In this issue, which is considered in a day-ahead time horizon, we have normal and emergency operation conditions together. As mentioned in the previous section, the occurrence of events in the resilience debate has uncertainties. Therefore, stochastic programming is considered to deal with these uncertainties. Some decisions must be made without complete information about some random events (H&N) and then decisions must be made for second stage (W&S) variables [38]. Depending on these conditions, different goals are pursued. These objectives are defined in terms of the distribution network and the set of multi microgrids. These goals are conflict with each other in a stochastic emergency situation and must be achieved at the same time under uncertainty situation.
Therefore, multi-objective two-stage stochastic scenario base modeling (MO-TSSP) base on mixed-integer linear programming (MILP) is formulated in this paper. This formulation is provided for both normal and emergency operation modes which are interdependent on each other. Simultaneously, the objectives and constraints present in each situation are different. Some constraints also apply in all situations, regardless of the circumstances. Also, the relationships between MMGs and distribution network are very important and challenging, which is addressed in this section. According to the general description, details and formulations are provided below.

A. STOCHASTIC SCHEDULING IN NORMAL OPERATION MODE
In this case, the system is in normal operation mode. MGs can exchange power with the substation according to the market price. Power supply resources in MGs also work optimally from an economic point of view.

OF
The first three terms of objective function which is presented in Eq. III-A1, Pd s t,i,j × OC i,j , ONC i,j × Y t,i,j and OFC i,j ×Z t,i,j , represent the operating cost, start-up cost, and shut-down cost of each DG in each of the MGs. The unit commitment (UC) problem is considered in each MG. The commitment of units is considered as a here and now variables and regardless of scenarios, decisions are made before the scenario dependent variables are determined. Further, the fourth to sixth term, OCW i × Pw s t,i , OCS i × Psl s t,i and OCE i × Pb dch,s t,i , correspond to the operation cost of renewable wind and solar energy sources and ESSs, respectively. The buying price of electricity from the utility grid and the profit from the sale of electricity to the utility grid in normal operating conditions, corresponding to the seventh and eighth terms, Pg s t,i ×MPb t,i and Pgp s t,i ×MPs t,i ,. The last term, Pev dch,s t,i,ev × Cev i,ev , also refers to the cost of discharging all PEVs. These are calculated in all normal operating periods and all scenarios.

1) POWER BALANCE AND RELATED CONSTRAINTS
The power balance in MGs, as well as power exchange with the upstream network, are important equations due to the technical constraints of the distribution network.
Power balance constraint Eq. (15) for each MG guarantees that for each scenario and at each period of the scheduled time, the total power generated by microturbines and RESs, ESS charge/discharge power, PEVs charge/discharge power, and power exchanged with the main grid is equivalent to the summation of prioritized load.
The maximum power exchanged between each MG and the main grid at Eq. (16) and Eq. (17) is considered on the MGs side due maximum power flow through the PCC of each MG in each scenario and timeslot. Not buying and selling simultaneously is discussed in Eq. (18). Load supplying in MGs at each level is required according to Eq. (19).
The power balance constraint which is presented in Eq. (20) on the distribution system side ensures that the sum of the prioritized loads on each bus is equal to the power injected from the substation and the powers exchanged by each MG with the main grid. The maximum power injected from the main grid according to the power exchanges between the MGs and the distribution feeder with the main grid is presented in Eq. (21). Corresponding to the active power balance equation, Eq. (22) indicates the reactive power balance on the distribution system side.
Due to the requirement for the presence of power converters to connect DC MGs to the distribution feeder, the equations related to their technical constraints [41] are shown in Eqs. (23) and (24). The relation between the active and reactive power of the PCC point converters of each MG is shown here [42]. These converters operate in both directions. During transferring power from the distribution feeder to the MGs, they operate as a rectifier and vice versa, as an alternator [42]. The rated capacity constraint of the converters with respect to the overall apparent power is satisfied by Eq. (25) [43], [41]. It is notable that the quadratic expression in Eq. (25) is linearized using the method presented in [41].
Through normal operation, no load shedding in active and reactive power is observed on the distribution system in Eqs. (26) and (27) The active and reactive power flow constraints is presented in Eq. (28) and (29). The linear approximation of distribution system power flow has been used in several optimization problems [44], [18] and [45]. Hence, in this paper, the Dis-tFlow linear equations according to Eqs. (30) and (31) are used. The allowable limit for the bus voltage in the distribution system is according to Eq. (32). Additionally, due to the limitations of active and reactive power flow of distribution lines, constraints Eqs. (33) and III-A2 have been applied.

2) STOCHASTIC UNIT COMMITMENT (SUC) PROBLEM IN MG
The SUC problem is applied in each MG all over the conditions. The technical requirements of each MT located in each MG are recognized. The two main sections in SUC equations include ramp rates and minimum up / down time constraints. Increasing or decreasing the power of MTs is due to their ramp-up/down rate limit which is formulated in Eqs. (35)- (37). The minimum and maximum time-dependent operational constraints are also given in equations Eqs. (38)- (39). These values are not necessarily equal to the maximum and minimum production capacity of each MT. The on/off states of units and their start-up/shut-down state is considered according to minimum up /down time given in Eqs. (40)-III-A3 [46], [47].

3) ENERGY STORAGE SYSTEMS (ESS) IN MGs
Applications and advantages of ESSs in research include: contributing to the large-scale integration of renewable energy, reducing the cost of power grid dispatching, and meeting the needs of different types of load [48]. The scheduling of ESSs to increase power system resilience in the operation [16] and planning [6] field can also be mentioned. Given the presence of RESs and the advantages of ESSs as well as their effects on resilience, they will be presented in MGs.
The relationship between charging/ discharging power and state-of-charge (SOC) are modeled in Eq. (44) according to battery efficiency [49]. Constraint which is presented in Eq. (45) due to battery lifetime corresponds to permissible SOC limits. The maximum limits of charging/discharging power, as well as the prevention of simultaneous charging and discharging, are imposed in constraints Eqs. (46)-III-A4, respectively.

4) RENEWABLE ENERGIES
Recently, the utility of RESs, including wind turbine and photovoltaic units in MGs has increased significantly [50]. Consequently, MGs should cope with significant operational challenges. The following equations are used to simulate the productive power of RESs.
The maximum power generated by the wind turbines and photovoltaic system at a time period is defined in Eqs. (49) and (51). The amount of renewable energy scheduled in MGs is determined by Eqs. (50) and III-A5.

5) PLUG-IN ELECTRIC VEHICLE (PEV)
The popularity of electric vehicles with rechargeable batteries by the grid is increasing due to their low pollution. In addition to the ability to absorb power from the network (G2V), PEVs can offer benefits and a useful concept called the ability to inject power into the network (V2G). EV owners may obtain more benefits using the energy stored in their vehicles. The PEV battery can be charged as well as discharged according to the car owner's consent.
Simultaneous non-charging and discharging in PEV batteries are shown according to Eq. (53). It should be noted that this operation takes place during the hours of the day when the PEVs are not traveling according to its travel pattern, and otherwise, the charging and discharging operation that requires a connection to the network is not performed, which is the Eq. (54) expresses this issue. The energy stored in the PEV battery, the charge and discharge limits of the power, and the amount of energy that can be stored according to their lifetime, are given in Eqs. (55)-III-B [51]. The system enters its critical operation state following a severe disaster that results in the loss of the upstream network. In this situation, the MMGs and distribution network pursue their own aims. Besides, each MG cooperates with its adjacent MG as well as the distribution system according to the existing constraints.

1) FIRST OBJECTIVE FUNCTION AND CONSTRAINTS ON MMGs
In emergencies, the purposes of MGs are to provide sensitive loads, and if all their loads are satisfied, the economic operation will be present in the next goals. It should be noted that in this case, there is no power market between the MGs and the upstream network so is removed from the objective function. The first three-term of the objective function which is presented in Eq. (59) contains operation cost, start-up cost, and shut-down cost of MTs in each MG. The operation costs of RES and ESS are considered in the fourth to sixth terms in the objective function, respectively. The seventh term also refers to the cost of discharging all PEVs. Due to MG loads providing requirements and the importance of maintaining sensitive and critical loads during emergency operations, the penalty cost of loading shedding has come in the last term. The important matter here is that higher priority load shedding has a higher penalty cost than lower priority loads.

OF
The power balance constraint which is presented in Eq. (60) also changes, so that the load shedding, the power exchanged between the MGs and the power injected from each MG into the distribution system are present.
The power exchange limit between each of the MGs (injected and received power) is shown in Eqs. (63) and (64). The power injection constraint from each MG to the distribution system is obtained from Eq. (65) and (66) according to the amount of power exchanged between the MGs. The constraints on the linearization of these relationships are given in Eqs. (67)-(69). The equilibrium constraint of the exchanged powers between all MGs, according to converters, is satisfied by Eq. III-B2.

2) SECOND OBJECTIVE FUNCTION AND CONSTRAINTS ON DISTRIBUTION SYSTEM
At starting the emergency operation, the distribution network attempts to minimize the energy curtailment, which is presented in Eq. (71). It has three types of loads with different levels. This means that on the distribution feeder, three goals are pursued so that they are superior to each other in order of importance.

OF LSh = s∈ sn
The power balance constraints of active and reactive power in the distribution network are satisfied by Eqs. (72) and (73).
About the active power, the sum of loads with different priorities is equal to the sum of the power received from the MGs and the load shedding in the distribution network. The same explanation is valid of reactive power. The loss of power exchange with the upstream grid is seen in Eq. (74).
The ratio of power curtailed in each bus due to maintenance of load power factor is considered in Eq. (75). The limitation of active and reactive load curtailment according to the load priority is also mentioned in Eqs. (76) and (77), respectively.
Eqs. (78) and (79) are technical constraints for converters according to power passing through them. Eq. (25) is also acknowledged here. The power flow is determined by Eqs. (80) and III-C in the conditions of a critical operation. The OPF constraints in this situation are the same as Eqs. (30)-III-A2. Resilience metrics are employed to evaluate the resilience in different conditions of a power system or MGs. Metrics is a suitable tool for comparing the resilience of a system in different situations when severe circumstances occur [52]. Grid resilience is an evolving research field hence hitherto the prescribed resilience metrics are new and no standard method has been defined for them. Thus, all resilience metrics are under comprehensive investigation [53]. Eq. (82) is used to calculate the resilience metric of the distribution system and MGs with regards to the amount of load supplied at the critical situation and the total load at that time [29] and [30]. Load priority from the highest to the lowest emphasis in the distribution system and MGs are considered in this formulation by Eq. (83). To calculate a normalized metric utilizing the maximum possible value of resilience metric, Eqs. (84)

IV. CASE STUDIES A. TEST SYSTEM AND MAIN ASSUMPTIONS
In order to evaluate the performance and efficacy of the proposed method, the suggested model is applied to a modified version of the IEEE 34-node test system, equipped with three MGs that are connected to specific buses [27]. These DC MGs are connected to the distribution system with bi-directional converters [43]. Subsequently, each of the MGs according to their available generation, sustain and enhance the resilience of other MGs and distribution system. The studied MGs have different sources. The first MG includes two kinds of microturbines (MTs), a photovoltaic (PV) unit, a wind turbine (WT), ESS, and prioritized load. The second MG includes six various MTs, a WT, ESS, and prioritized load. The third MG includes a MT, a WT, a PV unit, ESS, and prioritized load. Same as in practice, the parameters, capacities, and characteristics of each electrical device are different in mentioned MGs. On the other hand, with this assumption, the dependence of power generation sources is  diverse in each MG. The Contributions of each DG in total generation of MGs is given in Table 2.
Additionally, the load profiles and peak load of each MG are considered variously. It is assumed that the distribution network has no power generation sources owned. This assumption in the critical situation due to occurring a severe disaster, is supposed the worst-case scenario for distribution system as well as MGs. The mentioned test system is depicted in Figure 11 and the normalized profiles of load and RESs in MGs are exhibited in Figure 12.
Because the sources and elements employed in the three MGs are different from each other, their detailed technical information is classified in Table 3 based on their references, in order to avoid entering numerous information.
A summary of the essential information of the elements used in the MG is manifested in Table 4.
Moreover, the emergency operation durations follow a normal distribution with mean of 5 hours and standard deviation of one hour [16]. The probabilities associated with each scenario is indicated on the corresponding area, as shown in Figure 13.
Furthermore, the various aspects of the proposed framework are investigated in the following cases: Case I: A stochastic severe event causes the upstream grid outage. In this case as a base case, the distribution system and the MGs connected to it are not able to exchange power with the upstream grid in emergency condition. Moreover, PEVs are not considered here.
Case II: Same as Case I, except that the effect of PEVs is considered.
Case III: The reaction of electricity consumers and associated industries to electricity is affected by disasters in the electricity grid. This means that in addition to the technical issues in system (for example lines power flow and bus voltage limitations, and other power generations constraints), the social behaviors of the people affected by the disaster are also an important factor in managing the system. Therefore, to investigate this important factor, its effect on the behavior of electric vehicle owners illustrated in Figure 9 during emergency operation of the system is studied.

B. RESULTS AND ANALYSIS 1) BASE CASE (CASE I)
As mentioned in section IV-A, the event occurrence time is uncertain and therefore there will be 24 scenarios in case studies corresponding to 24 hours. There is also another uncertainty in each of those scenarios regarding the duration of the emergency operation. In addition, two objective functions are followed in each case. The result is exhibited as a Pareto Front  curve in Figure 14. The operation cost of MMGs (OCMG) and distribution system energy curtailment (DEC) are considered as the resilience evaluation criteria in this paper. The dependency of OCMG and DEC on event occurrence time is evidence in Figure 14. One of the important reasons for this issue is the load level of the distribution feeder and each MG, as well as their generation capacity, particularly RESs in critical time horizon. It should be noted that the operation cost of MMGs during critical times is significantly dependent on the penalty cost for load shedding in MGs. Therefore, in some hours, it is observed a high slope increase in operation cost, which indicates the possibility of load shedding in MGs. It can be asserted that the contingencies that occur around 18:00 are the worst possible situation. Figure 15 presents the compromised point which is chosen from the Pareto Front curve solutions. It can be deduced that during times which average sufficient sources of renewable generations are available if a disaster occurs, distribution feeder experiences less energy curtailment.   For example, from Figure 15 it can be perceived that in events that occur during 10:00-13:00, the average curtailed energy of the distribution network has a low value. In contrast, this value has a higher average in events that occur during the interval 15:00-20:00. In the events that occur after 20:00, due to the reducing load level of MGs and distribution system, it is verified reduction of curtailed energy in distribution feeder and the operation cost of MGs.
According to the presented issues, it can be affirmed that the OCMG and DEC as evaluation criteria provide a general overview of the system. Defining and considering an evaluation criterion for assessing resilience is the first step in developing resilience metrics. Therefore, the metrics based on resilience feature are used to examine the results more exactly. This metric quantifies the resilience of the distribution network ( Figure 16) and each of the MGs (Figure 17) in all uncertainty scenarios and at all critical operating intervals of the system.
The minimum value of the resilience metric is considered to be 0.75, which is equivalent to the total loss of the third level load in the distribution feeder and MGs. It means that the first and second level loads are definitely supplied. It should be noted that for 0.025 changes in the resilience metric, the load shedding changes by 10%. If interruption starts during 11:00-14:00 and lasts less than four hours, the resilience metric has a desirable value and in some scenarios is at maximum value in Figure 16. During these hours, the appropriate sources of renewable energy are available in MGs to deliver into the distribution system. It is worth mentioning that each MG deliver power to the distribution system due to the excess of available power. On the other hand, due to the downswing load level of MGs and distribution system, the resilience metric has a desired value due to the events that occur during 21:00-24:00. Another subject is related to scenarios that take the system in a critical state for a more lasting period. During this period, unexpectedly the average of this metric gradually takes on higher values, which is due to the decreasing load profiles of MGs and distribution system. This chart can be illustrated for each of the MGs. MGs, in addition to maintaining their resilience metric at a certain sufficient level, also sustain the distribution system resilience metric value. It can be inferred from Figure 17.b that the resilience metric of MG2 in most of the event occurring time and scenarios have values close to or equal to the maximum value. However, the resilience metric of first and third MGs decreases further in some scenarios and occurring times. These hours are near and beyond the sunset. At the same time, the MGs are experiencing peak loads. Additionally, due to the dependence of MG3 generations on WT, during 14:00-20:00, when the available WT power is diminishing, MG3 experiences a worse situation. By analyzing the power generation resources of MGs, it can be concluded that MG2 has more accessible and dispatchable resources than the other two MGs and is always able to generate power without dependence on environmental and climatic conditions. Following Figure 17 it can be confirmed that the minimum average of resilience metrics in the scenario sets in MG1, MG2, and MG3 is equal to 0.859, 0.978, and.787 respectively. Consequently, MG2 has a desirable qualification in the worse critical situations.

2) CASE II
The effect of electric vehicles on the resilience of the distribution system and MGs is investigated. The number of PEVs and the behavior pattern of electric vehicle owners in the three MGs are considered different from each other, as follows (Table 5): Presently, to investigate the impact of the PEVs presence, the results of the second case are compared to the first case.
For comparing the resilience evaluation criteria in the two cases, results are depicted in Figure 18. According to this  figure, the impact of PEVs on resilience criteria varies during the different occurring times of events. In the hours that PEVs are in the parking lot (1:00-8:00) and have not started their trip, the event occurrence increases the OCMG, but it has no significant effect on the DEC.
In events that occur during 17: 00-23: 00 in the presence of PEVs, the OCMG criteria decreases. Compared to the first case, discharge of PEVs supply the lost load of the MGs at the critical time, and therefore the penalty cost for the load shedding is reduced or equated to zero. Therefore, the DEC criteria decrease. This means that, due to discharging the PEVs, the portion of MGs available power is also injected into the distribution system, which is a positive point.
To adequately compare the results, compromise solutions are obtained as shown in Figure 19.
In Figure 19 occurring events after 21:00 increase the DEC criteria. During these periods, PEVs start charging due to providing the necessary energy to travel by the time they leave the parking lot. Consequently, the MG transferring power to the distribution system is reduced because of supplying the PEVs energy required. It can be verified from this figure that the most decrease in OCMG and DEC criteria is equal to 75% and 44%, respectively, which event occurs at 20:00. Moreover, the most increase in MGOC and DEC criteria is equal to 84% at 1:00 and 40% at 24:00, respectively.
To investigate the resilience metric, Figure 20 is presented. These figures reveal case I and II difference in resilience metrics.   Figure 20(b-d), it is concluded that in the scenarios which the MG resilience metric has increased, the distribution system resilience metric also experiences a significant increase in value. During events occurring 15:00-20:00, the distribution system resilience metric is increased generally. Before 7:00 and after 21:00, if a disaster occurs, in most scenarios are observed a decrease in the distribution system resilience metric. However, no change in the MGs resilience metric is observed at the same hours. Because the excess power that the MGs delivered to the distribution system is allocated to charging PEVs on those hours. For example, if an event occurs at 1:00, on average, 40.57%, 38.71%, and 59% of the charging in MG1, MG2, and MG3 respectively, is conducted in the emergency condition. Another notable issue is the extent and effectiveness of PEVs in resilience metrics of each MG. On the other hand, if the disaster starts almost at 18:00-23:00, MGs resilience metrics enhance desirable due to discharge of PEVs. For instance, when the event occurs at 20:00, on average, 93.2%, 95%, and 100% of discharge the PEV batteries perform during an emergency situation in MG1, MG2, and MG3 respectively.

3) CASE III
Comparing the optimal points of base case and case III in Figure 21, it is concluded that the reaction of electric vehicle owners to the occurred disaster plays one of the principal roles in the results. Because it increases operation costs, much of which depends on the load shedding penalty cost in MGs during emergency operation. In addition, it will cause part of the distribution load to be curtailed. According to different disaster start times, this behavior has various effects on obtained results. For example, if an event occurs at 17:00 and 18:00, the resilience criteria will be worse than at 8:00 or 9:00. Event happening time is an important factor for examining the impact of electric vehicle owners' behavior and social reaction. The positive strength of the PEVs presence on the   system resilience criterion becomes the weakness of the system within some scenarios in case III. During 13: 00-23: 00, the OCMG and DEC criteria have risen distinctly.
By considering the compromise answers, the most increase in the objective functions can be recognized. The worst time of the disaster occurrence is at 24:00, as the DEC and OCMG have increased 2.15 times and approximately 5 times, respectively. This is one of the hours when the presence of PEVs in the parking is the most overlapping with emergency operation condition of the system. The battery of PEVs start charging immediately after the event occurs. PEVs continue to be charged until their batteries are completely charged or the emergency condition ceases. This behavior provokes the OCMG to increase approximately 4 times compared to case II. This value can be considered as the penalty cost that MGs have to undergo for this reaction in an emergency situation. The DEC, which is related to the distribution system resilience evaluation criteria, in this case, has increased by 15.9% at 24:00 compared to case II. Figure 22 is achieved by obtaining compromise solutions for this case.
The DEC, which is related to the distribution system resilience evaluation criteria, in this case, has increased by 15.9% at 24:00 compared to case II. The general resilience evaluating criteria of MGs and distribution system were scrutinized in this part. In the following, resilience metrics are investigated to analyze the results more accurately. A comparison of Figure 20 and Figure 23 shows that the PEV owners' reaction to the disaster and the charging pattern of vehicles in the critical operation have reduced resilience metric of MGs. This proves that the behavior of energy consumers has a significant impact even on MGs that have adequate power generation. As can be observed, the impact of this social behavior on the second MG is more noticeable because the number of electric vehicles in this MG is more than the others. The third MG, is less affected by this event than the other MG. Because MG3 has both fewer PEVs and less PEV battery capacity. Despite these descriptions, the allowable limits for resilience metrics in the distribution system and all three MGs are satisfied and indicate high efficiency of the proposed approach. It should be noted that the approach for applying demand response in the studies are considered in normal circumstances.

V. CONCLUSION
In this study, using the optimal performance of existing dissimilar DC microgrids, maintaining and increasing the resilience of the distribution system connected to them has been concluded while maintaining the MGs resilience at the desired value. Throughout this purpose, under the uncertainties, both normal and emergency operations of the system are investigated. To achieve these goals, the proposed method is used, which is as a multi-period two-stage scenario-based stochastic mixed-integer linear programming (MPTSS-MILP) based on a multi-objective optimization problem (MOOP). Then, using the fuzzy satisfying method based on the conservative approach, compromise solutions are obtained among the non-dominated results obtained from a Pareto-Front curve.
Considering the Distflow model for power flow, voltage and line flow limitation in all situations, as well as the optimal management of power injected in contingencies by MGs to the distribution system, are conspicuous features considered at the distribution system level. Furthermore, considering the stochastic unit commitment of dispatchable generations, exploitation of RESs, and management of ESSs as well as electric vehicles under all conditions is the prominent considered at the MGs level.
Due to the objectives and constraints of the problem, this investigation was analyzed in three case studies. In the first case, which is considered as a base case, the study was conducted without the presence of PEVs. 24 scenarios are considered for a possible 24 hours of the severe event start time. In addition, within each assumed scenario for disaster start time, 7 scenarios are considered, which indicate the duration of the emergency condition. Therefore, the worst disaster start time according to OCMG and DEC as resilience evaluation criteria was determined. It was concluded that if the emergency condition started at 18:00 and lasted an average of 5 hours, from the resilience point of view worst condition would be experienced for microgrids and distribution networks. It is noteworthy that the distribution system resilience metric in this situation has a desirable value and is 17% better than the minimum allowable resilience metric. In the second case study, the effect of PEVs according to the owners' driving pattern and battery life of PEVs on the distribution network and MGs resilience was investigated. The results confirmed that the effect of PEVs on DEC and OCMG was highly dependent on the disaster occurrence time. According to the DEC, 23:00 and 24:00 and according to the OCMG, 1:00-and 2:00 are the worst disaster start time. On the other hand, if an emergency condition starts at 20:00 and 21:00, the DEC and OCMG decrease 44% and 75% respectively due to the operation of PEVs. The distribution network's resilience metric in this case also improves by 6.7% compared to the base case. Finally, in the third case study, the impact of the social behavior affected by disaster on the results was comprehensively discussed. The social behavior effect on the results confirmed that the beneficial performance of PEVs at critical conditions could turn into a weakness for the system. At 16:00, OCMG and DEC experience the greatest increase of 54% and 5.1% compared to the first case, respectively.
Also, the resilience metric of the first to third MGs and the distribution network had decreased by 14.19%, 15.7%, 13.8%, and 6.7%, respectively, compared to base. It can be said that the vastest effect of the social reaction is on microgrids rather than distribution networks. However, these values satisfied the resilience metric minimum allowable value. Results confirm the flexibility of the proposed model and it can satisfy the minimum requirements even in the worst cases of crisis.