Robust Day-Ahead Dispatch for Integrated Power-Heat-Gas Microgrid considering Wind Power Uncertainty

Wind power generation has been widely deployed in the modern power system due to the issues of energy crisis and environment pollution. Meanwhile, the microgrid is gradually regarded as a feasible way to connect and accommodate the distributed wind power generations. Recently, more research studies also focus on incorporating various energy systems, for example, heat and gas into the microgrid in terms of satisfying diﬀerent types of load demands. However, the uncertainty of wind power signiﬁcantly impacts the economy of the integrated power-heat-gas microgrid. To deal with this issue, this paper presents a two-stage robust model to achieve the optimal day-ahead economic dispatch strategy considering the worst-case wind power scenarios. The ﬁrst stage makes the initial day-ahead dispatch decision before the observation of uncertain wind power. The additional adjustment action is made in the second stage once the wind power uncertainty is observed. Based on the duality theory and Big-M approach, the original second-stage problem can be dualized and linearized. Therefore, the column-and-constraint generation algorithm can be further implemented to achieve the optimal day-ahead economic dispatch strategy for the integrated power-heat-gas microgrid. The experimental results indicate the eﬀectiveness of the presented approach for achieving operation cost reduction and promoting wind power utilization. The robustness and the economy of the two-stage robust model can be balanced, of which the performances signiﬁcantly outperform those of the single-stage robust model and the deterministic model.


Introduction
Currently, the generations of the intermittent energy resources (IER), especially the wind power generations, have been widely deployed in the modern power system due to energy crisis and environment pollution [1,2]. Meanwhile, the microgrid (MG) has been admitted as a feasible way to connect and accommodate the distributed wind power generations. In research studies [3][4][5], MG has been proved that it has remarkable ability to reduce the investment cost and the operation cost effectively. Recently, research studies [6][7][8][9][10][11][12][13] have also focused on incorporating various energy forms into MG to compose integrated energy-based MG, which is able to satisfy multiple energy demands simultaneously and has great potential to further accommodate IER.
However, there are two natural characteristics of wind challenge of the MG-based distributed wind power generation. e first one is that the uncertainty of wind power may impact the economy of the system operation because of the wind power curtailment [14]. e second one is that the wind power is relatively plentiful at night when the electricity load is relatively low [6]. erefore, to even out the wind power uncertainty and promote the wind power utilization, research studies claimed that the integration of additional energy equipments and technologies in MG is capable of storing or converting excessive wind power, thus effectively reducing the additional cost of wind power spillage [7][8][9][10][11][12][13]. For example, researchers in [7,8] pointed out that the cooperation of the energy storage system (ESS) and the heat storage system (HSS) successfully improves the flexibility of the combined heat and power (CHP) system and the wind power accommodation capability. Nielsen et al. [9] claimed that the deployment of the electric boiler (EB) effectively promotes the wind power utilization and reduces the operation cost by converting wind power to heat energy. Researchers in [10,11] reported that the employment of the power-to-gas (PtG) device converting unconsumed wind power to natural gas is highly capable of reducing the wind power spillage in the integrated power and gas system. Recently, studies have concentrated on the integrated power-heat-gas MG [12,13,15]. eir experimental results showed great significance with deploying ESS, HSS, PtG, and EB simultaneously in the power-heat-gas MG regarding the effective wind power utilization.
Although novel energy equipments and technologies are helpful for utilizing wind power, it should be realized that the uncertainty of wind power may affect the economic operation of MG due to imprecise prediction of wind [16][17][18][19]. Day-ahead ED determines the optimal allocation of the outputs of the controllable generators (CG) [19]. erefore, day-ahead economic dispatch (ED) [19,20] can be an effective solution to improve the economy of the integrated power-heat-gas MG involving wind power generations [21][22][23][24]. Among a wide variety of ED approaches, heuristic algorithms and mathematical optimization approaches have been regarded as the most effective solutions for achieving the optimal ED strategy. For example, researchers in [25,26] employed the differential evolution technique and the real-coded genetic algorithm to study the optimal ED problem for a CHP system, respectively. However, two main disadvantages of the heuristic algorithms cannot be simply neglected [27,28]. e first one is that the heuristic algorithms often suffer from local optimal issue. e second one is that the computational time is unstable. Consequently, the mathematical optimization approaches have been widely researched [21][22][23][24][25][26][27][28][29][30][31][32] due to its advantages of computational efficiency and global optimal solution guarantee.
Li et al. [21] presented a deterministic optimal ED model for a power-gas system. However, in their research, the uncertainty of wind power is not considered, which results in the accuracy loss issue of their optimization. Li et al. [13] further presented a stochastic optimization (SO) approach to achieve the optimal ED strategy for an integrated power-heat-gas MG involving the wind power uncertainty. However, Li et al. and Bertsimas et al. [22,23] pointed out that SO requires explicit and accurate distributions of the uncertainty which is an extreme difficulty in wind prediction. As a result, robust optimization (RO) that needs no accurate distributions of the uncertainty has been incrementally adopted in power system research studies to carry out the optimization tasks [16,[29][30][31][32]. Peng et al. [29] employed a single-stage robust optimization model to deal with the uncertainties of wind power and load demand. However, the research reported that the single-stage robust optimization tends to be overconservative. erefore, the two-stage robust optimization model with "min-max-min" structure has been studied in [16,[30][31][32] to deal with the uncertainty issues. For example, Sun et al. [30] employed a two-stage RO to investigate the unit commitment problem for a system with uncertain wind power. Liu et al. [31] introduced a twostage RO to achieve cooptimization of energy and reserve dispatch for an integrated power and gas system. Peng and Liu [16] presented a two-stage robust hydrothermal-wind ED model to improve the flexibility of the power system operation. According to the results of above researches although RO has successfully solved unit commitment [30], cooptimization of energy and reserve [31], energy management [32], and day-ahead ED problems [16], few research studies have investigated the potential of RO in enabling optimal day-ahead ED for the integrated powerheat-gas MG.
Motivated by the aforementioned research studies, this paper employs RO as the underlying approach to achieve the optimal day-ahead ED strategy for the integrated powerheat-gas MG involving the uncertain wind power. Firstly, a two-stage ED model is presented. In the first stage, the initial day-ahead decision is made before the realization of the wind power uncertainty. Once the uncertainty is observed, the additional corrective action is decided in the second stage to compensate the power imbalance. Secondly, regarding the nature characteristics of the wind power, an uncertainty set is constructed to describe the uncertainty of wind power. irdly, the strong duality theory and Big-M theory are employed to transform the original two-stage model into a mixed-integer linear problem (MILP). As a result, the original model is able to be solved directly by using the column-and-constraint generation (C&CG) algorithm. e rest of the paper is organized as follows: Section 2 constructs the two-stage model for the integrated powerheat-gas MG involving the uncertain wind power; Section 3 constructs the uncertainty set of wind power and presents the model solution; Section 4 shows and discusses the experiment results; Section 5 concludes the paper.

Two-Stage Day-Ahead ED for the Integrated
Power-Heat-Gas MG

Objective Function.
is paper focuses on achieving the optimal day-ahead ED strategy for integrated power-heatgas MG involving the uncertain wind power. Firstly, the structure of the researched MG is illustrated in Figure 1.
e two-stage robust model can be presented as the following equation: In equation (1), x represents the first-stage dispatch variables. e details of x are listed in Table 1. C da is the dayahead operation cost of MG, which is made before the observation of the uncertain wind power. e details of C da are shown in the following equation: where a FC and b FC are the coefficients of the fuel cell (FC); S FC represents the ON/OFF status of FC; λ gas is the price of the natural gas; λ buy Grid and λ sell Grid are the prices of the sold and purchased power with the grid; P sold Grid and P purchased Grid represent the sold and purchased power; and P firstwind and P forecasted represent the scheduled wind power in the first stage and the forecasted wind power, respectively.
C da consists of the operation cost of FC, the cost of the natural gas consumption, the benefit of the interaction with the grid, and the penalty cost of the wind power curtailment [33]. e first-stage strategy only considers the forecasted value of the wind power. However, the imprecise prediction may lead to power imbalance in the power-heat-gas MG. erefore, in the second stage, the decision C + da is the adjustment action after the observation of the wind power uncertainty.
Additionally, in equation (1), the second-stage max function searches for the worst-case scenarios which represents that the wind power is the most insufficient. e inner layer min function minimizes the adjustment cost under the discovered worst-case scenarios. u stands for the wind power. y is the second-stage adjustment variables. e details of y are listed in Table 2. C + da is shown in the following equation: where and λ down CG are the prices of the up/downregulation of CG; λ up gas and λ down gas are the prices of the up/downregulation of the natural gas supply; and P s wind is the uncertain wind power.
C + da consists of the additional adjustment cost of CG, the additional adjustment cost of natural gas supply, the adjustment cost of the exchange power with the grid, and the penalty cost of the wind curtailment once the uncertain wind is observed. In this paper, CG represents FC, microturbine, EB, and PtG.
A da (·) and A + da (·) in equation (1) represent the firststage constraints and adjustment constraints, respectively. e details of these constraints are discussed in Section 2.2.
where P dis ESS and P ch ESS are the discharging/charging power of ESS; P MT is the power output of the microturbine; P PtG and P EB are the electric power consumed by PtG and EB; Q dis HSS and Q ch HSS are the discharging/charging power of HSS; P load , Q load , and G load are the electric, heat, and gas load demands, respectively; Q MT and Q EB are the heat energy generated by the microturbine and EB; G gas is the total amount of supplied natural gas; G MT is the gas consumption of the microturbine; and G PtG is the gas generated by PtG.
In our model, the microturbine is a back-pressure CHP which generates electric power by consuming natural gas. Its converting ratio of the electric and heat power output is a constant value [34]. erefore, Q MT , P MT , and G MT can be presented in the following equation: where η MT is the constant ratio of converting the electric and heat power output; η MT is the converting efficiency of the gas consumption and the electricity output; similarly, G PtG and Q EB can also be calculated as follows: where η PtG and η EB represent the conversion efficiency of PtG and EB. In this paper, η MT , η gas , η PtG , and η EB are equal to 0.9, 0.85, 0.85, and 0.8, respectively [12,[34][35][36][37].
(2) CG Constraints. e constraints of CG include ramping rate and the limit of power output, as shown as follows: where P max CG,i represents the maximum power output of the i th CG; R up CG,i is the maximum ramping rate of the i th CG; and S CG,i represents the ON/OFF status of CG.
e storage system (SS) mainly includes ESS and HSS. In each hour t, the energy balance of SS can be constrained as follows: Moreover, the final stored energy should be equal to its initial stored energy, as shown in the following equation: (4) Interaction Constraints. e exchange power between MG and the power grid should be limited as follows: where P max Grid is the maximum exchange power with the grid and S purchased Grid is a binary variable, which represents the ON/ OFF status of the purchasing mode.
(5) Wind Power Constraint. e scheduled wind power in the first stage should not exceed the forecasted wind power [30], as shown in the following equation:

Adjustment Constraints
(1) Energy Balance Constraints. Once the uncertain wind power is observed, the second-stage decision should be made to compensate the power imbalance in MG. Power, heat, and gas energy balances should be satisfied after the adjustment action in the second stage as shown in the following equation: 4 Mathematical Problems in Engineering where P (2) Regulation Constraints where S up CG,i is a binary variable, which represents the ON/ OFF status of the i th CG. Equations (15)- (17) ensure that the ramping rate and the power output are within the permissible range.
Similarly, the gas power regulation of the gas supplier and the regulation of SS should be within the allowable range. And simultaneous up-and downregulation of the gas supplier and SS should be avoided, which is expressed as follows: where P ch,adjust SS,i and P dis,adjust SS are the charging and discharging power of the i th SS after the adjustment and E adjust SS,i is the stored energy after the adjustment. Notably, charging/ discharging mode remains the same status as in the firststage decision.

(3) Interaction Constraints
where S purchased Balance is a binary variable, which represents the selling/purchasing power status; P max Balance is the maximum adjustment power with the grid; the limits of the sold and purchased adjustment power should be constrained by equation (21). After the adjustment, the selling/purchasing mode with the grid keeps the same as the first-stage exchange status. e exchange power with the grid should not exceed the permissible range, which is expressed in equation (22).

(4) Wind Power Constraint.
After the uncertainty of the wind power is revealed in the second stage, the actual wind power should be less than the observed wind power, which is indicated as follows:

Uncertainty Set.
To guarantee the operation robustness, the formulation protects against all possible outputs of the wind power within an uncertainty set U. e budget of uncertainty Γ is also employed to adjust the degree of the conservativeness during the dispatch period T. e entire uncertainty set U adopted by this paper is presented in the following equation: where P error wind is the maximum forecasted error of the wind power. e value of Γ controls the conservativeness of the optimized result.

Compact Formulation. In terms of facilitating the presentation, the aforementioned objective function equations (1)-(3) and the constraints of the two-stage ED model equations (4)-(24) can be rewritten into a compact form:
where a, b, A, B, C, and D u are the constant matrices and d and e represent vectors.

Model Decomposition.
e presented model has a complicated min-max-min structure, which cannot be solved directly. However, C&CG has been proved to be an effective algorithm for solving the complex two-stage model with remarkable computational efficiency and convergence capability [38]. erefore, C&CG is employed to solve the presented two-stage ED model. Based on C&CG, the original model is firstly decomposed into the master problem (MP) and the subproblem (SP). e optimal first-stage decision x is decided in MP given the possible scenarios from SP, while SP searches for the possible worst-case scenarios given the first-stage decisions from MP. Afterwards, the problem is able to converge to the optimal solution during the iterations between MP and SP. erefore, MP can be rewritten as follows: Bx + Cy k ≤ e, where η is an auxiliary variable and y k and u k are the optimal adjustment decisions and the possible wind power scenarios at the k th iteration, respectively. By gradually identifying the possible wind power scenarios at each iteration, the updated constraint equations (27) 3.4. Transformation of SP. SP contains a max-min structure, which is difficult to be solved. erefore, strong duality theory is applied to transform SP into an equivalent singlestage problem SP1 as shown in the following equation: max u,π 1 ,π 2 e T π 1 − x * T B T π 1 + u T π 2 s.t.
where π 1 and π 2 are the dual variables for the constraints of SP equation (30) and u T π 2 is a nonconvex bilinear term, which also leads to the difficulties of directly solving SP1. However, u only contains the forecasted value u 0 , the upper limit value u + , or the lower limit value u − , which is not continuous. erefore, u T π 2 can be linearized by applying the Big-M method. Meanwhile, the uncertain budget is also introduced to adjust the uncertainty set. As a result, the transformed SP2 can be expressed by equation (32) which is able to be solved in a feasible way: where π + 2 and π − 2 represent the positive and negative values of π 2 , corresponding to the upper and lower bounds of the wind power and M is a large enough constant.

Model Solution Using C&CG.
e solution of the model using C&CG includes the following major steps given in Table 3.

Experimental Result
In this section, a number of evaluations are carried out to demonstrate the effectiveness of the presented ED optimization approach for an integrated power-heat-gas MG. e parameters of FC, the microturbine, EB, and the PtG devices are listed in Table 4. e parameters of ESS and HSS are listed in Table 5. e cost of gas and the adjustment cost of gas are 0.85 $/m 3 and 1.05 $/m 3 . e details of different load demands including power, heat, and gas are shown in Figure 2. e forecasted wind power and the price of the grid are shown in Figure 3. e maximum allowable fluctuation is set to 10%.

Experimental Results of Robust Optimization Regarding
Wind Power Utilization. In this section, to demonstrate the effectiveness of the presented approach for reducing wind power curtailment, the optimized results of the outputs of SS, CG, and the exchange power of the two-stage robust model are shown in Figures 4 and 5. Figure 4 indicates that when the wind power is plentiful during the 1 st -7 th and the 22 nd -24 th hours, the excessive wind power is converted into heat and gas energy by EB and PtG, respectively. erefore, the wind power utilization is considerably promoted. Figure 5 shows that ESS successfully improves the wind power accommodation when the wind power is plentiful during the 1 st -6 th and the 21 st -23 rd hours. e reason is that during these periods, ESS is able to keep charging to store the unconsumed wind power. Moreover, HSS keeps discharging heat power in the 23 rd -24 th hours. Consequently, the electricity and the heat outputs of the microturbine decrease at the same time, which creates an extra margin for wind power accommodation. Figure 5 indicates that MG sells power to the grid when the wind is plentiful, which improves the wind power utilization. Additionally, MG also sells power when the sold price is relatively high during the 9 th -11 th hours and 19 th -23 rd hours. erefore, the profit of MG can be significantly improved.

Impact of Uncertain Budget.
In this section, different values of the uncertain budget Γ � 0, 8, 16, and 24 are employed to study its impact on optimal day-ahead ED. In our model, Γ represents the allowable fluctuation times of the wind power. For example, if Γ equals 0, it indicates that in the period T, wind power should equal its forecasted wind power. If Γ equals 8, it means that the fluctuation times of the wind power should be no more than 8 times. As long as the optimal day-ahead ED strategy is achieved, another 500 scenarios are generated by Monte Carlo simulation as testing scenarios to evaluate the performance of the achieved optimal strategy. e optimized results based on different values of Γ are shown in Table 6.   Initialize lower and upper bounds UB � +∞ and LB � −∞. Set the iteration no. k � 0 and the convergence gap ε > 0.
Step 2 Solve MP. Derive an optimal solution (x * k , η k , y 1 , . . ., y k ), and update the lower bound LB � a T x k + η k .
Step 3 Solve SP3 with respect to the fixed (x * k , η, y 1 , . . ., y k ) from MP. Obtain the solution (u o k , y o k ) and update u k+1 � u o k and UB � min If UB − LB ≤ ε, return the solution x � x * k and the algorithm terminates. Otherwise, create a variable y k+1 and add constraint equations (27)- (29) to MP. And then update k � k+1 and go back to step 2.   Table 6 demonstrates that the day-ahead operation cost of the robust optimization keeps increasing along with the uncertain budget value rises from 0 to 24. As Γ is highly related to the number of the fluctuations, an increasing Γ leads to more fluctuation times of the wind power, which means wind power is more likely to become     Table 7, the total scheduled output of generators becomes higher to immunize against the wind power uncertainty. Although the increasing values of Γ guarantee the robustness of the model, the day-ahead operation cost becomes higher due to the increase in the scheduled output of the microturbine and FC. To sum up, larger Γ leads to worse economy but better robustness of the ED strategy. erefore, the uncertainty budget can be flexibly adjusted to balance the economy and robustness of the ED strategy in practical applications.

Impact of the Forecasted Error.
In this section, different maximum forecasted errors 5%, 10%, 15%, and 20% of the wind power are considered to reveal its impact on the dayahead ED. Also 500 scenarios are generated as testing scenarios to evaluate the economy of the achieved optimal strategy. e optimized results based on different forecasted errors are illustrated in Table 8.
As shown in Table 8, the day-ahead operation cost and the adjustment cost increase along with the increasing forecasted error. is is due to the fact that with the increasing forecasted error, the worst-case wind power becomes more insufficient. As a result, more resources and adjustment output are scheduled to handle the uncertain wind power in the robust optimization model, which causes the increasing costs.

Impact of the Confidence Interval.
In this section, confidence intervals 80%, 85%, 90%, and 95% of the wind power are considered to reveal its impact on the day-ahead ED. 500 scenarios are generated as testing scenarios to evaluate the effectiveness of the achieved optimal strategy. e optimized results based on different confidence intervals are illustrated in Table 9.
As shown in Table 9, the day-ahead operation cost and the adjustment cost decrease steadily as the confidence interval increases. is reveals the fact that with the increasing confidence interval, the uncertainty set narrows and the value of the worst-case wind power is closer to the forecasted wind power. erefore, fewer resources and adjustment output are needed to handle the uncertain wind power, which reduces the total operation cost.

Optimization Comparison with Different Approaches.
In this section, the deterministic model and the single-stage robust model are also implemented in terms of comparison. eir cost reduction optimization performances are compared to those of the presented two-stage robust model. Both the uncertainty budgets of the single-stage and the two-stage robust model are set to 8. e comparisons of different approaches are shown in Table 10.
According to Table 10, the cost of the two-stage robust model is lower than that of the single-stage model. e single-stage model schedules less wind power and more generation output than those of the two-stage robust model, as shown in Table 11. is is due to the fact that the singlestage model only considers one worst-case wind power scenario. However, the two-stage robust model is able to search for a branch of worst-case wind power scenarios by iterations between the first stage and the second stage. erefore, the two-stage model is less conservative than the single-stage model. In other words, the two-stage robust  model can better balance the robustness and economy than those of the single-stage model. It is also can be observed that the day-ahead operation cost of the two-stage robust model is higher than that of the deterministic model. e two-stage robust model considers the worst scenarios of the uncertain wind power. erefore, more resources and more operation cost are needed in the day-ahead strategy, but less adjustment cost is demanded, as shown in Figure 6. In this way, the total operation cost of the two-stage robust model is effectively reduced compared to that of the deterministic model. Figure 7 shows that during the two-stage robust model solution, the presented C&CG algorithm converges to the global optimal solution after only 4 iterations. is point significantly indicates the efficiency of the presented approach.

Conclusion
is paper presents a two-stage robust model in enabling optimal day-ahead ED for integrated power-heat-gas MG considering wind power. To be aware of the natural characteristics of wind, this paper employs an uncertainty set to describe the uncertainty of wind power. Additionally, to facilitate the problem solution, the original model with a complex min-max-min structure is decomposed into MP and SP with min and max-min structures, respectively, by the C&CG algorithm. Based on the duality theory and Big-M approach, the original SP can be successfully linearized.
erefore, the decomposed model can be directly solved by iterations between SP and MP. Experimental results demonstrate that the two-stage robust model can handle the worst-case wind power scenarios to achieve the optimal dayahead ED strategy. Compared to the single-stage robust model and the deterministic model, the presented approach can ensure a less-conservative and more-convincible solution. Moreover, the budget of uncertainty enables the     adjustment of the optimization performance for the dispatch scheme, so that the robustness and the economy of the ED strategy can be better balanced. In our future work, cooling energy supplying and uncertainty of load demands will be further explored.
Data Availability e data of the models and algorithms used to support the findings of this study are included within the article.