Welfare Maximization-Based Distributed Demand Response for Islanded Multi-Microgrid Networks Using Diffusion Strategy

Integration of demand response programs in microgrids can be beneficial for both the microgrid owners and the consumers. The demand response programs are generally triggered by market price signals to reduce the peak load demand. However, during islanded mode, due to the absence of connection with the utility grid, the market price signals are not available. Therefore, in this study, we have proposed a distributed demand response program for an islanded multi-microgrid network, which is not triggered by market price signals. The proposed distributed demand response program is based on welfare maximization of the network. Based on the welfare function of individual microgrids, the optimal power is allocated to the microgrids of the network in two steps. In the first step, the total surplus power and shortage power of the network is determined in a distributed way by using the local surplus/shortage information of each microgrid, which is computed after local optimization. In the second step, the total surplus of the network is allocated to the microgrids having shortage power based on their welfare functions. Finally, the allocated power amount and the initial shortage amount in the microgrid is used to determine the amount of load to be curtailed. Diffusion strategy is used in both the first and the second steps and the performance of the proposed method is compared with the widely used consensus method. Simulation results have proved the effectiveness of the proposed method for realizing distributed demand response for islanded microgrid networks.


Introduction
Demand response (DR) programs encourage consumers to reduce their electricity demands when required. This helps energy consumers to adjust their energy consumption in response to the variations in electricity prices and incentives offered by utilities, resulting in bill reduction. DR helps in regulating the load at peak hours and during system contingencies [1,2]. Therefore, it can help to avoid the need for building new power plants. Further, it can help in enhancing the utilization of renewables by adjusting loads and improving system reliability [3]. Thus, DR can provide economic assistance to both consumers and suppliers [3,4]. The DR problem has been an active research area in recent years and various methods have been suggested for its implementation [1,5,6]. 1) In existing studies, distributed DR is implemented on grid-connected microgrids which requires market price signals and a central aggregator. However, market prices are not available for islanded microgrids. Thus, existing methods cannot be applied for islanded microgrids. Therefore, we have proposed a distributed DR which is not dependent on market price and can be implemented for islanded microgrids. 2) In existing studies, distributed DR is applied to only single microgrids. We have proposed a distributed DR for a multi-microgrid system, which is more suitable for islanded microgrids to reduce load shedding by sharing power. 3) Diffusion strategy is utilized in the proposed method; hence, the convergence time is reduced in comparison with conventionally used distributed methods, thus making it more suitable for re-optimization during system contingencies.

System Model
An islanded network of I microgrids is considered in this study, where I is any finite number. The microgrids of the network are interconnected, i.e., they can exchange power with other microgrids in the network. However, being an islanded network, the microgrids cannot trade power with the utility grid. Each microgrid contains dispatchable generators (DGs), a battery energy storage system (BESS), renewable energy generators (RDGs), and different priority loads. Figure 1 shows an illustration of the proposed microgrid network.
Each microgrid in the network is equipped with a local energy management system (EMS), named as microgrid EMS (MG-EMS). The MG-EMS is responsible for collecting all data from its local components, deciding optimal schedules, and informing each component of the results. After local optimization, the MG-EMSs of the network share their surplus and shortage information to determine optimal power sharing among the microgrids. Power sharing among the microgrids at each interval is based on the welfare of individual microgrids in that particular interval. Details of the welfare function are provided in the following section. (distributed DR) are also compared with the results of the consensus-based welfare maximization method. The main contributions for this study are as follows. 1) In existing studies, distributed DR is implemented on grid-connected microgrids which requires market price signals and a central aggregator. However, market prices are not available for islanded microgrids. Thus, existing methods cannot be applied for islanded microgrids. Therefore, we have proposed a distributed DR which is not dependent on market price and can be implemented for islanded microgrids. 2) In existing studies, distributed DR is applied to only single microgrids. We have proposed a distributed DR for a multi-microgrid system, which is more suitable for islanded microgrids to reduce load shedding by sharing power. 3) Diffusion strategy is utilized in the proposed method; hence, the convergence time is reduced in comparison with conventionally used distributed methods, thus making it more suitable for reoptimization during system contingencies.

System Model
An islanded network of I microgrids is considered in this study, where I is any finite number. The microgrids of the network are interconnected, i.e., they can exchange power with other microgrids in the network. However, being an islanded network, the microgrids cannot trade power with the utility grid. Each microgrid contains dispatchable generators (DGs), a battery energy storage system (BESS), renewable energy generators (RDGs), and different priority loads. Figure 1 shows an illustration of the proposed microgrid network.
Each microgrid in the network is equipped with a local energy management system (EMS), named as microgrid EMS (MG-EMS). The MG-EMS is responsible for collecting all data from its local components, deciding optimal schedules, and informing each component of the results. After local optimization, the MG-EMSs of the network share their surplus and shortage information to determine optimal power sharing among the microgrids. Power sharing among the microgrids at each interval is based on the welfare of individual microgrids in that particular interval. Details of the welfare function are provided in the following section.

Welafare Function
Each microgrid of the network is considered to be an individual entity with different load profiles and different load priorities. Based on the load profiles and the ratio of critical loads in each microgrid, a welfare function can be defined for each microgrid. The welfare function implies the satisfaction level of a microgrid for consuming an amount of energy. The welfare function is a non-decreasing concave function [23], as illustrated in Equations (1) and (2). It implies that the microgrid is always interested in consuming more power and the level of satisfaction gradually gets saturated. Equation (1) describes the marginal benefit, which is the derivative of the welfare function and is a decreasing function, as described in Equation (2).
In literature, several welfare functions have been proposed, but due to the suitability of quadratic welfare functions to mimic the satisfaction level of consumers, they are widely used [8,23]. Due to these merits, a quadratic welfare function is used in this study also, as given by Equation (3). In Equation (3), α is a predetermined parameter, x is the amount of energy consumption, and w is load priority determining parameter. With the help of w, we can distinguish energy consumers in accordance with the importance of their loads. Figure 2a shows the examples of graphical representations of welfare functions for α = 0.3 and different values of w (w1 = 3, w2 = 12, w3 = 7), corresponding to three different microgrids. For a fixed energy consumption x, larger w implies a higher welfare. Since w2 > w3 > w1, for the same power consumption, the welfare of MG2 is highest and the welfare of MG1 is lowest, as illustrated in Figure 2a. An example of the marginal benefit function for w = 1 is shown in Figure 2b, which represents the marginal benefit of a microgrid for different power consumptions.

Welafare Function
Each microgrid of the network is considered to be an individual entity with different load profiles and different load priorities. Based on the load profiles and the ratio of critical loads in each microgrid, a welfare function can be defined for each microgrid. The welfare function implies the satisfaction level of a microgrid for consuming an amount of energy. The welfare function is a nondecreasing concave function [23], as illustrated in Equations (1) and (2). It implies that the microgrid is always interested in consuming more power and the level of satisfaction gradually gets saturated. Equation (1) describes the marginal benefit, which is the derivative of the welfare function and is a decreasing function, as described in Equation (2).
In literature, several welfare functions have been proposed, but due to the suitability of quadratic welfare functions to mimic the satisfaction level of consumers, they are widely used [8,23]. Due to these merits, a quadratic welfare function is used in this study also, as given by Equation (3).
In Equation (3),  is a predetermined parameter, x is the amount of energy consumption, and w is load priority determining parameter. With the help of w, we can distinguish energy consumers in accordance with the importance of their loads.  and different values of w (w1 = 3, w2 = 12, w3 = 7), corresponding to three different microgrids. For a fixed energy consumption x, larger w implies a higher welfare. Since w2 > w3 > w1, for the same power consumption, the welfare of MG2 is highest and the welfare of MG1 is lowest, as illustrated in Figure  2a. An example of the marginal benefit function for w = 1 is shown in Figure 2b, which represents the marginal benefit of a microgrid for different power consumptions.

Diffusion Strategy
Several methods are available in the literature for distributed operation of microgrids [24][25][26]. Diffusion strategy, as described in Equation (4), has the benefit of fast convergence as compared to other distributed operation methods. This is due to the utilization of a gradient term in an intermediate step which is dispersed to neighboring nodes, allowing information to diffuse more rapidly to the neighborhood of any given node [24]. A sample of a communication topology with i microgrids is shown in Figure 3; here, each MG-EMS is referred to as a node, which can only share information with its neighboring nodes. As shown in Figure 3, node j can only share information with a group of nodes Nj, within its communication reach (highlighted in Figure 3). The network's

Diffusion Strategy
Several methods are available in the literature for distributed operation of microgrids [24][25][26]. Diffusion strategy, as described in Equation (4), has the benefit of fast convergence as compared to other distributed operation methods. This is due to the utilization of a gradient term in an intermediate step which is dispersed to neighboring nodes, allowing information to diffuse more rapidly to the neighborhood of any given node [24]. A sample of a communication topology with i microgrids is shown in Figure 3; here, each MG-EMS is referred to as a node, which can only share information with its neighboring nodes. As shown in Figure 3, node j can only share information with a group of nodes Nj, within its communication reach (highlighted in Figure 3). The network's communication topology is represented by adjacent matrix A, determined through the Metropolis rule [24], as given in Equation (5). Di f f usion Strategy : where φ i,k−1 is the intermediate state of node i, ∇s i,t (φ i,k−1 ) is the stochastic gradient for node i, Q i,k is the state of node i, and u > 0 is a non-negative updating parameter.
where a ij is the entry of adjacent matrix A and n i and n j are the number of neighboring nodes of agent i and j.  [24], as given in Equation (5).
where aij is the entry of adjacent matrix A and ni and nj are the number of neighboring nodes of agent i and j.

Problem Formulation
In this section, the mathematical model of the proposed distributed DR scheme for an islanded MMG network is presented. The proposed model is formulated for a scheduling horizon of T with a time interval of t, which can be any uniform interval of time. The step-by-step operation of a microgrid in the network is shown in Figure 4. The microgrid operation is divided into two steps, namely local optimization and distributed DR. In the first step, the generation capabilities of DGs, RDGs, and BESS, along with the load profile of each microgrid are taken as input data. After receiving all data, each MG-EMS performs local optimization and ensures optimal operation of its components. Since each microgrid is being operated in an islanded mode, power exchange with the external grid is not possible. In cases where the available power is less than the power demand, load shedding is carried out to maintain the stable operation and power balance of the network. Local optimization is a parallel process for all the microgrids in the network. After local optimization, each MG-EMS has

Problem Formulation
In this section, the mathematical model of the proposed distributed DR scheme for an islanded MMG network is presented. The proposed model is formulated for a scheduling horizon of T with a time interval of t, which can be any uniform interval of time. The step-by-step operation of a microgrid in the network is shown in Figure 4. The microgrid operation is divided into two steps, namely local optimization and distributed DR. In the first step, the generation capabilities of DGs, RDGs, and BESS, along with the load profile of each microgrid are taken as input data. After receiving all data, each MG-EMS performs local optimization and ensures optimal operation of its components. Since each microgrid is being operated in an islanded mode, power exchange with the external grid is not possible. In cases where the available power is less than the power demand, load shedding is carried out to maintain the stable operation and power balance of the network. Local optimization is a parallel process for all the microgrids in the network. After local optimization, each MG-EMS has the information of local shortage power and surplus power available in its microgrids. Depending upon the importance of the amount of shortage power, MG-EMS decides a welfare function.
After local optimization by the MG-EMS, prior to the second step, MG-EMSs share necessary information with each other and wait until all necessary information is shared among the microgrids in the network. This step is carried out to overcome the issues which might occur due to communication delays in the network. The second step, i.e., distributed DR, is to determine the optimal power allocation and the load curtailment. In order to determine optimal power sharing, the information of total shortage and surplus power available in the network is required. Using Algorithm 1, we can determine the total available surplus and shortage power in the network in a distributed way without using any central entity. The total surplus is optimally shared to microgrids with shortages using Algorithm 2. After local optimization by the MG-EMS, prior to the second step, MG-EMSs share necessary information with each other and wait until all necessary information is shared among the microgrids in the network. This step is carried out to overcome the issues which might occur due to communication delays in the network. The second step, i.e., distributed DR, is to determine the optimal power allocation and the load curtailment. In order to determine optimal power sharing, the information of total shortage and surplus power available in the network is required. Using Algorithm 1, we can determine the total available surplus and shortage power in the network in a distributed way without using any central entity. The total surplus is optimally shared to microgrids with shortages using Algorithm 2.
The proposed operation method has been developed for an islanded MMG network. During the islanding of the microgrid, the MG-EMS should switch from normal operation algorithm to islanded operation algorithm. Various algorithms are available in literature for distributed DR in gridconnected mode [9,10,15,[20][21][22][23]. Microgrids can utilize any of the algorithms available in the literature for distributed DR in grid-connected mode, while the proposed method can be utilized during islanded mode. The existing literature lacks research on distributed DR schemes for islanded microgrids; therefore, in this study, we have developed an operation algorithm for islanded microgrid networks. The detailed mathematical modeling of local optimization and distributed DR is explained in the following sections.

Local Optmization
The objective of the proposed optimization model is to minimize the total operational cost of each microgrid, as shown in Equation (6). The first line of Equation (6) corresponds to all costs related to DGs, i.e., generation cost, startup cost, and shutdown cost, respectively. The first term of the second line is the penalty cost when load shedding is carried out. The second term of the second line is the incentive price for the microgrid to show their surplus. The constraints for DG are given by Equations (7)-(10). Equation (7) represents the operation bounds of DG at time t. On-off status of DGs is determined by Equation (8). The start-up and shutdown statuses of DG are determined using Equations (9) and (10), respectively. The BESS model can be represented by Equations (11)- (14). The maximum and minimum bounds for charging and discharging BESS units are given by Equations (11) and (12), respectively. At each interval, the state of charge (SOC) is computed according to  The proposed operation method has been developed for an islanded MMG network. During the islanding of the microgrid, the MG-EMS should switch from normal operation algorithm to islanded operation algorithm. Various algorithms are available in literature for distributed DR in grid-connected mode [9,10,15,[20][21][22][23]. Microgrids can utilize any of the algorithms available in the literature for distributed DR in grid-connected mode, while the proposed method can be utilized during islanded mode. The existing literature lacks research on distributed DR schemes for islanded microgrids; therefore, in this study, we have developed an operation algorithm for islanded microgrid networks. The detailed mathematical modeling of local optimization and distributed DR is explained in the following sections.

Local Optmization
The objective of the proposed optimization model is to minimize the total operational cost of each microgrid, as shown in Equation (6). The first line of Equation (6) corresponds to all costs related to DGs, i.e., generation cost, startup cost, and shutdown cost, respectively. The first term of the second line is the penalty cost when load shedding is carried out. The second term of the second line is the incentive price for the microgrid to show their surplus. The constraints for DG are given by Equations (7)-(10). Equation (7) represents the operation bounds of DG at time t. On-off status of DGs is determined by Equation (8). The start-up and shutdown statuses of DG are determined using Equations (9) and (10), respectively. The BESS model can be represented by Equations (11)- (14). The maximum and minimum bounds for charging and discharging BESS units are given by Equations (11) and (12), respectively. At each interval, the state of charge (SOC) is computed according to Equation (13). The constraint for maximum and minimum bounds for SOC at time t, is given in Equation (14). The constraint for power balance for a single microgrid is given by Equation (15). The first two terms of (15) are the DG and RDG power generation, respectively, the third and fourth terms are the discharging and charging power for BESS, the fifth and sixth terms are the local surplus power and shortage power in the microgrid, while the last term is the electric load of microgrid at interval t. Equation (15) shows that the power generated by DGs, RDGs, BESS discharging, and amount of shortage power should be balanced with BESS charging power and load amount.

Distributed Demand Response
After all of the MG-EMSs have completed local optimization for their microgrids and the necessary initial information has been received and shared, the next is step is distributed DR. The proposed distributed DR method comprises two steps, namely information sharing and optimal power allocation. Both steps are carried out in a distributed way using diffusion strategy, as described in detail in the following sections.

Information Sharing
In the proposed distributed DR strategy, the first step is the exchange of information between all the MG-EMSs in the network to determine the total available surplus and deficit power in the network. The detailed mechanism of information sharing is explained in Algorithm 1. In order to determine the total surplus in the network, at each iteration the node i updates its current state s sur k−1,i to a new state s sur k,i using the local stochastic gradient at the iteration. The local stochastic gradient can be calculated from the difference of the intermediate state φ sur at this iteration [24]. Similarly, the total shortage in the network can be determined. By the end of Algorithm1, the amount of surplus power and shortage power in the network is known. According to the proposed approach, each MG-EMS participates in determining the total surplus and shortage power. Each MG-EMS in the network, shares information in a distributed way until all MG-EMSs reach a consensus, which is used to compute the total surplus/shortage power of the network. The proposed approach requires no central entity to determine total surplus/shortage power of the network.
compute the total surplus/shortage power of the network. The proposed approach requires no central entity to determine total surplus/shortage power of the network.

Optimal Power Allocation
In the proposed distributed DR strategy, the second step is to maximize the welfare of the network by optimal power sharing among microgrids. To maximize the welfare of the distributed  compute the total surplus/shortage power of the network. The proposed approach requires no central entity to determine total surplus/shortage power of the network.

Optimal Power Allocation
In the proposed distributed DR strategy, the second step is to maximize the welfare of the network by optimal power sharing among microgrids. To maximize the welfare of the distributed : end for 10: Update error 11: end while 12: end if

Optimal Power Allocation
In the proposed distributed DR strategy, the second step is to maximize the welfare of the network by optimal power sharing among microgrids. To maximize the welfare of the distributed network, Energies 2019, 12, 3701 9 of 18 the global welfare objective function is described in Equation (16); here, W i,t is the welfare function of each microgrid. Therefore, the goal of the objective function (17) is to maximize the global welfare of the network and is obtained by summing all the individual welfare functions of each microgrid.
Subject to: In contrast to the first step, i.e., information sharing, the gradient of the welfare function is used as the stochastic gradient (∇W i (φ k−1,i )). The detailed mechanism is explained in Algorithm 2. In order to ensure the optimal solution lies inside a feasible region as described in Equation (18), penalty function (19) is used. With the addition of the penalty function, the objective function described in (17) is redefined as in (20): Here, η is a parameter which represents the importance of the penalty function and z is a tuning scalar. Penalty function (19) is a squared deviation from the constrained function described in (18), tuned by the scalar z.

Numerical Simulations
In order to realize the effectiveness of the proposed distributed DR scheme, an islanded MMGs network of five microgrids was considered in this study, as illustrated in Figure 1. Each microgrid includes a DG, RGDs (wind turbine or photovoltaic array), a BESS, and electric load. The load within each microgrid is prioritized depending upon its importance, and its priority is represented using a welfare function. All microgrids in the network are interconnected and can share power with each other. In addition to the power exchange, MG-EMSs can also share information in a distributed way. The proposed method is formulated for a 24 h scheduling horizon with a time step of 1 h. For each time step t, the algorithm takes some iterations to reach a convergence point. In the simulation results, iteration numbers are the repetition of the same algorithm until it reaches a consensus, which is the optimal value. The simulations are carried out in a MATLAB environment by integrating it with IBM ILOG CPLEX [27] on a laptop (Intel Core i7-7500U CPU @ 2.70GHz, 8GB RAM). Figure 5 shows the simulation setup schematic and represents the individual roles of MATLAB and IBM ILOG CPLEX, for microgrid operation.

Input Data
DGs have the main role in optimizing the operation cost of microgrids. The parameters associated with DGs considered for each microgrid are given in Table 1. A BESS can store excess power during the intervals when demand is low and discharge it in peak load intervals to reduce load shedding. The parameters related to each microgrid's BESS are tabulated in Table 2. The penalty costs associated with load shedding at MG1, MG2, MG3, MG4, and MG5 are 500, 430, 350, 300, 260, and 240, respectively in Korean Won (KRW). The forecasted values of load and output power of the RDGs, taken as input, are shown in Figure 6a,b, respectively. The input data, i.e., DG parameters, BESS parameters, and forecasted values of load and renewable generation, are taken from [28]. Results of local optimization based on the input data mentioned above are discussed in the next section.

Input Data
DGs have the main role in optimizing the operation cost of microgrids. The parameters associated with DGs considered for each microgrid are given in Table 1. A BESS can store excess power during the intervals when demand is low and discharge it in peak load intervals to reduce load shedding. The parameters related to each microgrid's BESS are tabulated in Table 2. The penalty costs associated with load shedding at MG1, MG2, MG3, MG4, and MG5 are 500, 430, 350, 300, 260, and 240, respectively in Korean Won (KRW). The forecasted values of load and output power of the RDGs, taken as input, are shown in Figure 6a,b, respectively. The input data, i.e., DG parameters, BESS parameters, and forecasted values of load and renewable generation, are taken from [28]. Results of local optimization based on the input data mentioned above are discussed in the next section.

Input Data
DGs have the main role in optimizing the operation cost of microgrids. The parameters associated with DGs considered for each microgrid are given in Table 1. A BESS can store excess power during the intervals when demand is low and discharge it in peak load intervals to reduce load shedding. The parameters related to each microgrid's BESS are tabulated in Table 2. The penalty costs associated with load shedding at MG1, MG2, MG3, MG4, and MG5 are 500, 430, 350, 300, 260, and 240, respectively in Korean Won (KRW). The forecasted values of load and output power of the RDGs, taken as input, are shown in Figure 6a,b, respectively. The input data, i.e., DG parameters, BESS parameters, and forecasted values of load and renewable generation, are taken from [28]. Results of local optimization based on the input data mentioned above are discussed in the next section.

Local Optimization Results
Results of local optimization for each microgrid are shown in Figure 7. In Figure 7a-c respectively, it can be seen that the MG1, MG2, and MG3 have shortages during peak load intervals, whereas, due to the high generation capabilities of MG4 and MG5, they do not have shortage power throughout the day, as shown in Figure 7d,e. Microgrid-wise surplus power and shortage power at each interval of time are shown in Table 3. The negative values indicate the shortage power, whereas the positive values indicate surplus power in Table 3. Similarly, in Figure 7, battery charging power is taken as negative and battery discharging power is taken as positive.

Local Optimization Results
Results of local optimization for each microgrid are shown in Figure 7. In Figure 7a, 7b, and 7c respectively, it can be seen that the MG1, MG2, and MG3 have shortages during peak load intervals, whereas, due to the high generation capabilities of MG4 and MG5, they do not have shortage power throughout the day, as shown in Figures 7d,e. Microgrid-wise surplus power and shortage power at each interval of time are shown in Table 3. The negative values indicate the shortage power, whereas the positive values indicate surplus power in Table 3. Similarly, in Figure 7, battery charging power is taken as negative and battery discharging power is taken as positive.  After local optimization, MG-EMS decides the value of the w parameter depending upon the priority of shortage power. A w parameter value is selected between 0 and 100. The available surplus at each interval can be shared to microgrids with shortages in an optimal way to decide load curtailment (demand response) in the network. The results for distributed DR are discussed in detail in the next section.

Distributed Demand Response Results
In this study, we have considered two different cases in order to illustrate the effectiveness of the proposed distributed DR method. In the first case (case 1), the difference between the shortages in microgrids is not significant, whereas in the second case (case 2), shortage power in one microgrid is significantly lower as compared to other microgrids in the network. For each case (case 1 and case 2), the total information of surplus and shortage power in the network is determined first, followed by the optimal power allocation to microgrids with shortage power. Meanwhile, during the intervals when the total surplus of the network is greater or equal to the shortage in the network, the allocated  After local optimization, MG-EMS decides the value of the w parameter depending upon the priority of shortage power. A w parameter value is selected between 0 and 100. The available surplus at each interval can be shared to microgrids with shortages in an optimal way to decide load curtailment (demand response) in the network. The results for distributed DR are discussed in detail in the next section.

Distributed Demand Response Results
In this study, we have considered two different cases in order to illustrate the effectiveness of the proposed distributed DR method. In the first case (case 1), the difference between the shortages in microgrids is not significant, whereas in the second case (case 2), shortage power in one microgrid is significantly lower as compared to other microgrids in the network. For each case (case 1 and case 2), the total information of surplus and shortage power in the network is determined first, followed by the optimal power allocation to microgrids with shortage power. Meanwhile, during the intervals when the total surplus of the network is greater or equal to the shortage in the network, the allocated power to each microgrid is equal to the shortage power in that microgrid. The results of information sharing and optimal power allocation for case 1 and case 2 are discussed in the next section.

Case 1: Information Sharing
By applying Algorithm 1, the total available surplus power and shortage power in the network are determined for interval 10, as shown in Figure 8. Interval 10 has been taken as a representative interval of case 1. The same method is implemented for other intervals. In order to determine the total shortage power in the network, each MG-EMS shares information with its neighbors. The shortage power information of all the MG-EMSs converge to a single value, i.e., 63.4 kW, which is the average of the shortage power in the network (Figure 8a). In a similar way, the total surplus power in the network is determined for the same interval, shown in Figure 8b, where the information converges to 38.8 kW.
Energies 2019, 12, x FOR PEER REVIEW 12 of 19 power to each microgrid is equal to the shortage power in that microgrid. The results of information sharing and optimal power allocation for case 1 and case 2 are discussed in the next section.

Case 1: Information Sharing
By applying Algorithm 1, the total available surplus power and shortage power in the network are determined for interval 10, as shown in Figure 8. Interval 10 has been taken as a representative interval of case 1. The same method is implemented for other intervals. In order to determine the total shortage power in the network, each MG-EMS shares information with its neighbors. The shortage power information of all the MG-EMSs converge to a single value, i.e., 63.4 kW, which is the average of the shortage power in the network (Figure 8a). In a similar way, the total surplus power in the network is determined for the same interval, shown in Figure 8b, where the information converges to 38.8 kW.

Case 1: Optimal Power Allocation
By using Algorithm 2, we can maximize the welfare of the MMG network, by determining the optimal power allocated to each microgrid with shortage power. The optimal power allocated to microgrids depends upon the welfare function of each microgrid. The results for optimal power allocation for interval 10 are shown in Figure 9 and Table 4. In this case, w2 > w3 > w1; therefore, the power allocated to MG2 is the highest and the power allocated to MG1 is the lowest. In this case, none of the microgrids received the allocated power equal to their initial shortage, even though the value of w2 is much higher than w1 and w3. This is due to the limited amount of surplus power available in the network.

Case 1: Optimal Power Allocation
By using Algorithm 2, we can maximize the welfare of the MMG network, by determining the optimal power allocated to each microgrid with shortage power. The optimal power allocated to microgrids depends upon the welfare function of each microgrid. The results for optimal power allocation for interval 10 are shown in Figure 9 and Table 4. In this case, w2 > w3 > w1; therefore, the power allocated to MG2 is the highest and the power allocated to MG1 is the lowest. In this case, none of the microgrids received the allocated power equal to their initial shortage, even though the value of w2 is much higher than w1 and w3. This is due to the limited amount of surplus power available in the network. power to each microgrid is equal to the shortage power in that microgrid. The results of information sharing and optimal power allocation for case 1 and case 2 are discussed in the next section.

Case 1: Information Sharing
By applying Algorithm 1, the total available surplus power and shortage power in the network are determined for interval 10, as shown in Figure 8. Interval 10 has been taken as a representative interval of case 1. The same method is implemented for other intervals. In order to determine the total shortage power in the network, each MG-EMS shares information with its neighbors. The shortage power information of all the MG-EMSs converge to a single value, i.e., 63.4 kW, which is the average of the shortage power in the network (Figure 8a). In a similar way, the total surplus power in the network is determined for the same interval, shown in Figure 8b, where the information converges to 38.8 kW.

Case 1: Optimal Power Allocation
By using Algorithm 2, we can maximize the welfare of the MMG network, by determining the optimal power allocated to each microgrid with shortage power. The optimal power allocated to microgrids depends upon the welfare function of each microgrid. The results for optimal power allocation for interval 10 are shown in Figure 9 and Table 4. In this case, w2 > w3 > w1; therefore, the power allocated to MG2 is the highest and the power allocated to MG1 is the lowest. In this case, none of the microgrids received the allocated power equal to their initial shortage, even though the value of w2 is much higher than w1 and w3. This is due to the limited amount of surplus power available in the network.   Table 6. If s is correct.
Commented Table 6 from Figure 9. Results for optimal power allocation (case 1).

Case 2: Information Sharing
As with case 1, Algorithm 1 is used to determine the total available surplus power and shortage power in the network for interval 17. Interval 17 has been taken as the representative interval of case 2. As shown in Figure 10a, the shortage power information of all the MG-EMSs converges to a single value, i.e., 48.8 kW, which is the average of shortage power in the network. In a similar way, the total surplus power in the network is determined for the same interval as shown in Figure 10b, where the information converges to 18.6 kW. As with case 1, Algorithm 1 is used to determine the total available surplus power and shortage power in the network for interval 17. Interval 17 has been taken as the representative interval of case 2. As shown in Figure 10a, the shortage power information of all the MG-EMSs converges to a single value, i.e., 48.8 kW, which is the average of shortage power in the network. In a similar way, the total surplus power in the network is determined for the same interval as shown in Figure 10b, where the information converges to 18.6 kW.

Case 2: Optimal Power Allocation
The results for optimal power allocation for case 2 are shown in Table 5 and Figure 11. In this case, w1 > w2 > w3; therefore, the power allocated to MG1 is more than MG2. Despite of low value of w for MG3, it gets allocated power equal to its shortage power. This is because of significantly low shortage power at MG3 and comparatively low difference between values of w. The proposed method decides the power allocation according to the importance of load, rather than deciding power allocation based on size of the microgrid. Figure 11. Results for optimal power allocation (case 2).

Effect of Distributed Demand Response on the Network
The optimal power allocated to all microgrids in the network during each interval is given in Table 6. As mentioned earlier, MG4 and MG5 have no shortages so the power allocated to them is zero in each interval. The overall effect of the proposed distributed DR method can be observed from Figure 12. The load curtailed is the amount of load shed in the whole network and load served is the amount of total load which receives power. Available power is the maximum generation potential of the network and total load is the original load of the network. Microgrid-wise load curtailment and

Case 2: Optimal Power Allocation
The results for optimal power allocation for case 2 are shown in Table 5 and Figure 11. In this case, w1 > w2 > w3; therefore, the power allocated to MG1 is more than MG2. Despite of low value of w for MG3, it gets allocated power equal to its shortage power. This is because of significantly low shortage power at MG3 and comparatively low difference between values of w. The proposed method decides the power allocation according to the importance of load, rather than deciding power allocation based on size of the microgrid. As with case 1, Algorithm 1 is used to determine the total available surplus power and shortage power in the network for interval 17. Interval 17 has been taken as the representative interval of case 2. As shown in Figure 10a, the shortage power information of all the MG-EMSs converges to a single value, i.e., 48.8 kW, which is the average of shortage power in the network. In a similar way, the total surplus power in the network is determined for the same interval as shown in Figure 10b, where the information converges to 18.6 kW.

Case 2: Optimal Power Allocation
The results for optimal power allocation for case 2 are shown in Table 5 and Figure 11. In this case, w1 > w2 > w3; therefore, the power allocated to MG1 is more than MG2. Despite of low value of w for MG3, it gets allocated power equal to its shortage power. This is because of significantly low shortage power at MG3 and comparatively low difference between values of w. The proposed method decides the power allocation according to the importance of load, rather than deciding power allocation based on size of the microgrid. Figure 11. Results for optimal power allocation (case 2).

Effect of Distributed Demand Response on the Network
The optimal power allocated to all microgrids in the network during each interval is given in Table 6. As mentioned earlier, MG4 and MG5 have no shortages so the power allocated to them is zero in each interval. The overall effect of the proposed distributed DR method can be observed from Figure 12. The load curtailed is the amount of load shed in the whole network and load served is the amount of total load which receives power. Available power is the maximum generation potential of the network and total load is the original load of the network. Microgrid-wise load curtailment and

Effect of Distributed Demand Response on the Network
The optimal power allocated to all microgrids in the network during each interval is given in Table 6. As mentioned earlier, MG4 and MG5 have no shortages so the power allocated to them is zero in each interval. The overall effect of the proposed distributed DR method can be observed from Figure 12. The load curtailed is the amount of load shed in the whole network and load served is the amount of total load which receives power. Available power is the maximum generation potential of the network and total load is the original load of the network. Microgrid-wise load curtailment and power allocation can be seen in Figure 13. MG4 and MG5 did not have shortage; therefore, load curtailment is only carried out in MG1, MG2, and MG3. The amount of load to be curtailed in each interval depends upon the initial shortage power in the microgrid and power allocated to that microgrid. For instance, at interval 11, the initial shortage power at MG2 and MG3 is 180 kW and 64 kW, respectively, whereas surplus power in MG1 MG4, and MG5 is 7.36 kW, 101 kW, and 11 kW, respectively. It can be seen that the total shortage power is more than the total surplus power in the network. Power allocated to MG2 and MG3 is 64 kW and 55.36 kW, respectively. Since w2 > w3, the power allocated to MG2 is more than the power allocated to MG3. The amount of load curtailment in MG2 and MG3 is 116 kW and 8.64 kW, respectively.  power allocation can be seen in Figure 13. MG4 and MG5 did not have shortage; therefore, load curtailment is only carried out in MG1, MG2, and MG3. The amount of load to be curtailed in each interval depends upon the initial shortage power in the microgrid and power allocated to that microgrid. For instance, at interval 11, the initial shortage power at MG2 and MG3 is 180 kW and 64 kW, respectively, whereas surplus power in MG1 MG4, and MG5 is 7.36 kW, 101 kW, and 11 kW, respectively. It can be seen that the total shortage power is more than the total surplus power in the network. Power allocated to MG2 and MG3 is 64 kW and 55.36 kW, respectively. Since w2 > w3, the power allocated to MG2 is more than the power allocated to MG3. The amount of load curtailment in MG2 and MG3 is 116 kW and 8.64 kW, respectively.   power allocation can be seen in Figure 13. MG4 and MG5 did not have shortage; therefore, load curtailment is only carried out in MG1, MG2, and MG3. The amount of load to be curtailed in each interval depends upon the initial shortage power in the microgrid and power allocated to that microgrid. For instance, at interval 11, the initial shortage power at MG2 and MG3 is 180 kW and 64 kW, respectively, whereas surplus power in MG1 MG4, and MG5 is 7.36 kW, 101 kW, and 11 kW, respectively. It can be seen that the total shortage power is more than the total surplus power in the network. Power allocated to MG2 and MG3 is 64 kW and 55.36 kW, respectively. Since w2 > w3, the power allocated to MG2 is more than the power allocated to MG3. The amount of load curtailment in MG2 and MG3 is 116 kW and 8.64 kW, respectively.

Convergence Analysis
In order to analyze the fast convergence of the proposed DR program, we also performed a convergence comparison between the consensus algorithm and the proposed method based on diffusion strategy. Consensus algorithm and its applications have been extensively studied in the microgrid network [29][30][31][32]. The optimal power allocation determined by using the proposed method based on diffusion strategy, for case 1 and case 2, are shown in Figure 14a,b, respectively. The same results were obtained from the consensus algorithm; the results of optimal power allocation for case 1 and case 2 using the consensus algorithm are shown in Figure 14a,b, respectively. Result comparison shows that the proposed method converges faster compared to the consensus. Table 7 shows a summary of the comparison between the consensus algorithm and the proposed method, for both the cases. The reduction of iterations for case 1 and case 2 are 88.8% and 89.01%, respectively.

Convergence Analysis
In order to analyze the fast convergence of the proposed DR program, we also performed a convergence comparison between the consensus algorithm and the proposed method based on diffusion strategy. Consensus algorithm and its applications have been extensively studied in the microgrid network [29][30][31][32]. The optimal power allocation determined by using the proposed method based on diffusion strategy, for case 1 and case 2, are shown in Figure 14a,b, respectively. The same results were obtained from the consensus algorithm; the results of optimal power allocation for case 1 and case 2 using the consensus algorithm are shown in Figure 14a,b, respectively. Result comparison shows that the proposed method converges faster compared to the consensus. Table 7 shows a summary of the comparison between the consensus algorithm and the proposed method, for both the cases. The reduction of iterations for case 1 and case 2 are 88.8% and 89.01%, respectively.

Implementation Complexity and Practical Implementation
Conventionally, DR programs are executed in a centralized way; these algorithms require a powerful central controller to collect the global information and process massive amounts of data. Moreover, centralized approaches are costly and make systems more susceptible to single point failure, which may result in an outage in the whole system [33,34]. In order to overcome these issues, distributed/decentralized methods are considered as a possible solution. In the distributed approach, instead of a central controller, each unit in the system is controlled by a local controller; therefore, the computational burden is divided among different local controllers. The system configuration of a distributed network makes the system more robust to communication failure [33,34]. In this study we have utilized a fully distributed approach based on diffusion strategy. The proposed method inherits the features of the distributed approach as mentioned above, which makes the proposed approach less complex to implement and more robust regarding communication failure, as compared to the conventional DR programs.
In this study we have used diffusion strategy, which require few iterations to determine the optimal values. Time taken by the proposed method to determine the total surplus power, shortage power, and optimal power sharing in the network for case 1 and case 2 is 0.058 s and 0.043 s, respectively. If the proposed approach is to be implemented in a hierarchical control approach, where usually 15 min is used as a time step, the proposed approach will easily fit in the 15 min time step, as the time taken by the proposed approach for a single time step is enough to fit in a 15 min time step.

Implementation Complexity and Practical Implementation
Conventionally, DR programs are executed in a centralized way; these algorithms require a powerful central controller to collect the global information and process massive amounts of data. Moreover, centralized approaches are costly and make systems more susceptible to single point failure, which may result in an outage in the whole system [33,34]. In order to overcome these issues, distributed/decentralized methods are considered as a possible solution. In the distributed approach, instead of a central controller, each unit in the system is controlled by a local controller; therefore, the computational burden is divided among different local controllers. The system configuration of a distributed network makes the system more robust to communication failure [33,34]. In this study we have utilized a fully distributed approach based on diffusion strategy. The proposed method inherits the features of the distributed approach as mentioned above, which makes the proposed approach less complex to implement and more robust regarding communication failure, as compared to the conventional DR programs.
In this study we have used diffusion strategy, which require few iterations to determine the optimal values. Time taken by the proposed method to determine the total surplus power, shortage power, and optimal power sharing in the network for case 1 and case 2 is 0.058 s and 0.043 s, respectively. If the proposed approach is to be implemented in a hierarchical control approach, where usually 15 min is used as a time step, the proposed approach will easily fit in the 15 min time step, as the time taken by the proposed approach for a single time step is enough to fit in a 15 min time step.

Conclusions
A novel distributed DR program for an islanded MMG network has been proposed using diffusion strategy. The proposed method is based on welfare maximization, which aims to maximize the welfare of the network by maximizing the welfare of each microgrid in the network. Welfare maximization of the network is achieved by deciding the optimal power sharing among the microgrids without utilizing any central entity. Power allocation and load curtailment is decided based on the importance of load, rather than the size of the microgrid. Moreover, the proposed method has a fast convergence speed as compared to the consensus algorithm-based welfare maximization. The number of iterations were reduced by almost 89% as compared to the consensus algorithm.

Conflicts of Interest:
The authors declare no conflicts of interest.

Indices t
Index to time, running from t to T g Set of generations, running from g to G i Set of microgrids, running from i to I

Variables and Parameters
W g,t (·) Welfare function C DG g,t Per-unit generation cost of DG unit g at t C SU g Start-up cost of DG unit g at t C SD g Shut-down cost of DG unit g at t C P Penalty cost for load shedding C inc Incentive price for showing surplus power P RDG t Forecasted output of RGD at t P B_Cap Capacity of BESS P Load Forecasted load at t P min g Minimum generation limit of DG unit g P max On/Off status of DG unit g at t y g,t Start-up status of DG unit g at t z g,t Shut-down status of DG unit g at t P min