A medium-term hybrid IGDT-Robust optimization model for optimal self scheduling of multi-carrier energy systems

,


Introduction
Nowadays, researchers tend to study multi-energy systems due to the synergic effects of the multi-carrier energy systems and their valuable merit in increasing the efficiency of the system [1].The cooperation and management of electric, thermal, and gas units to meet the customers' demands form the energy hub (EH) [2].The EH may consist of different sources such as wind turbine (WT), combined heat and power (CHP) plant, photovoltaic (PV), electric vehicles (EV), power to gas (P2G), as well as storage units of electric or thermal energy [3].The most crucial feature of EH is its capability in converting one type of energy to other types to meet the demands of other types of energy and assure the highest profit or reliability [4].Therefore, different timely topics such as EH planning and sizing, optimal operation of EH, demand response in the EH, the effect of uncertainty on the EH, etc., has become the primary focus of both the industrial and academic societies [5].

Literature review
Generally, one of the research directions is the short-term scheduling of EH.In Ref. [6], the authors have proposed a bi-level operating strategy of the EH to find the optimized day-ahead schedule at the higher level while applying a model predictive control to track the real-time schedule at the lower level.The higher-level optimization was directed to minimize the fuel cost and electrical load shedding for a 24-h horizon.The outcome of the higher-level trajectories was used as reference trajectories of the lower-level real-time model predictive control scheme, resulting in a shorter prediction horizon with better resolution.The authors in Ref. [7] investigated risk-constrained stochastic scheduling of EH using conditional value at risk method where uncertainties related to wind speed, solar irradiation, all demand types, and electricity market prices were managed by k-means data cluster analysis.In Ref. [8], the optimum operation of multiple EHs was investigated.
The EH consisted of WT, PV, CHP, P2G, and gas-fired generators to provide interdependencies of infrastructures for these energy sources.This article used the probability density function to consider the uncertainties of WT and PV production.A multiobjective optimization problem was applied to maximize social welfare while minimizing carbon dioxide emissions, and a genetic algorithm was deployed to handle the problem.Additionally, the Newton-Raphson method was applied to determine the flow of heat, gas, and electricity.The optimal operation of EH was considered in Ref. [9] taking into account the uncertainties related to electrical, heating, and cooling demands, as well as the wind power generation via a two-stage stochastic programming model.In this work, the Monte-Carlo simulation was used to generate initial scenarios, and the K-means method was applied to reduce the number of scenarios properly.The optimal scheduling of an EH was presented in Ref. [10] taking into account the uncertainties of WT and PV generation, natural gas, and electricity prices using the probability density function.In Ref. [11], an enhanced grey wolf optimization approach was used to handle optimal operation and scheduling problems of an EH via a scenario-based stochastic model aiming at maximizing its profit.In this problem, uncertainties related to PV generation and electricity price were taken into account by applying the Monte-Carlo simulations and k-means method on historical data.A stochastic programming approach was developed in Ref. [12] to handle the uncertainties of WT and PV generation in an EH considering different configurations for the N-1 contingency.The investigated literature in operation and management of EHs within a short-term time horizon shows that the majority of works deal with the uncertainties via stochastic programming [7,9,11,12] and probability density function [8,10], while their main purpose is to show how to manage the EH in the presence of uncertainties.

Binary variables associated with bilateral contracts
The medium-term and long-term horizons are the other focuses of researchers in EH studies.A chance-constrained optimization model for economic dispatch and operation of EH for a 1-year horizon was investigated in Ref. [13].In this work, the original nonlinear problem was recast into a linear programming problem, and the outcomes indicated that considering a well-tuned chanceconstrained model directly influenced the final solution in a costeffective manner.A MILP model was studied in Ref. [14] to consider a one-year time horizon with an hourly resolution to guarantee the optimal EH operation, including seasonal energy storage.The proposed methodology was then used to design an EH in a neighborhood in Switzerland, aiming at optimizing the total annual costs and carbon dioxide emissions.Under this investigation, no uncertainty issues were included.A stochastic programming approach was suggested in Ref. [15] for supplying the required energy of clients, taking into account the uncertainties of rival hub managers and electricity price.In the proposed approach, the hub manager is placed in the upper level, and the clients are at the lower level.Using the KarusheKuhneTucker optimality conditions and the strong duality theorem, the bilevel nonlinear stochastic program was recast into a linear MILP one.Additionally, the conditional value-at-risk (CVaR) indicated the effects of risk in the trading decisions of the hub manager.The authors in Ref. [16] proposed a tool for the medium-term time management of the EH to make a self-scheduling for a large consumer subject to wind speed and electricity price uncertainties.The problem was modeled based on a five-stage stochastic programming approach to minimize the total cost of EH in a risk-averse perspective via CVaR.The number of operating hours significantly increases within the medium-term and long term.Therefore, one of the main purposes in the medium-term is to define a framework to reduce the number of hours.The authors in Refs.[15,16] reduced the number of hours to half in order to decrease the computational time.In addition, some researches concentrated on the significance of the design aspects of the EH rather than uncertainties [14].
Also, the current trend concerns the application of robust optimization and the IGDT method in the EHs.The authors in Ref. [17] focused on the hybrid stochastic/IGDT optimization method to investigate the optimal scheduling of WT in EH.The considered uncertainties were wind power generation, electrical energy prices, and energy demands trying to enable the EH operator in selecting proper risk-seeker and risk-averse strategies using a mixed-integer nonlinear programming model.Ref. [18] discusses the riskconstrained self-scheduling of the EH plant based on IGDT to find the worst-case realization of wind power generation.It was then applied to investigate the opportunity and robustness of the models in order to assure the minimum target benefit of the riskaverse EH, if wind power generation is lower than the forecasted value.An IGDT-based robust scheduling strategy was proposed in Ref. [19] for optimal on-grid EH coordination under high penetration of wind power generation.This model aimed at minimizing the EH operational costs, wind power curtailment, and carbon emissions.The authors in Ref. [20] studied an optimal load dispatch model for a community EH, which reduced the total cost of the EH, including operation and emission costs of the system.In this work, the RO was deployed to tackle the uncertainties of the electricity prices while the electrical and thermal demand responses were included to shift the electrical and thermal loads from the peak intervals to the off-peak.In Ref. [21], the authors presented a robust-based mixed-integer linear programming model to handle the day-ahead scheduling of the EH integrated with fuel cell-based hydrogen storage system.The main scheduling was directed to convert energy to hydrogen during low electricity prices and converted back during higher electricity prices.Additionally, the RO method was applied to consider the uncertainty of electricity price.
The authors in Ref. [22] proposed an optimal operating model for managing multiple EHs with electrical and heat energy demands, taking into account the electrical and thermal demand response programs.The operations were directed to reduce the total cost of the EH, while the uncertainties of WT generation and electricity prices were considered, and the RO method was then proposed to deal with it.In Ref. [23], a decentralized robust optimization model was presented to take into account the multi-hub as multi-entities, and an Alternating Direction Method of Multipliers (ADMM) was used to solve the problem.Moreover, the RO model was deployed to tackle the electricity price uncertainty.The work in Ref. [24] proposed a min-max-min robust framework for the short-term operation of EH-based microgrids.The model was linearized, and then the column-and-constraint generation procedure was applied to decompose the framework into a master problem and a subproblem.The master problem tried to minimize the unit commitment cost, while the sub-problem indicated the dispatch cost associated with uncertainties via a max-min objective function.The investigated literature regarding the application of robust optimization in the operation of EH demonstrates that several papers focused on the importance of applying the IGDT or RO in the EH to reach the worst-case realization of uncertainties [18,19,23].Besides [17], commenced a new way toward hybrid optimization in the EH by deploying a combination of stochastic programming and IGDT approach.

The challenges and paper contributions
The literature shows valuable works related to the optimal operation of an EH.However, to the best of our knowledge, the electricity market prices and renewable energy resources have seldom been taken into account together in previous studies in the EH.In addition, no reported study has considered a hybrid IGDT-RO model to deal with the electricity market prices and renewable energy sources.Modeling both uncertainties results in non-linear terms in the IGDT model due to the presence of electricity price in the objective function [25].On the contrary, the IGDT approach provides a situation to consider the binary variables easily.On the other hand, modeling both uncertainties using robust optimization needs complicated tools such as Benders Decomposition or column-and-constraint generation procedure to solve the problem [24].In this situation, increasing the number of binary variables may result in intractability.The situation will be more complicated by considering medium or long-term problems that tackle with more variables.Therefore, for the sake of simplicity, researchers often neglect some binary variables such as the charging/discharging limitation of the energy storage.Therefore, this paper tries to capture the positive aspects of both methods by presenting a hybrid IGDT-RO model.The positive aspects are twofold; 1) The proposed hybrid model takes into account both uncertainties with the ability to consider a large number of binary variables while including a big number of integer or binary variables is always challenging for the decomposition methods such as the benders decomposition approach.2) Both uncertainties are included in the problem with the hybrid IGDT-RO method through a linear model, while it is not possible to present a linear model only by the IGDT approach.It means the proposed hybrid model is a linear framework that can host various binary variables.By doing so, a mixedinteger linear programming model is obtained, which can be solved using commercial software.Moreover, the literature review reveals that robust medium-term management has rarely been studied in previous works.Consequently, the problem also investigated the EH management in a medium-term time horizon.Plus, the proposed model is integrated with the P2G storage system in order to better deal with the operation of the natural gas network.
To summarize, the main contributions of the article are as follows.
Proposing a linear hybrid IGDT-robust optimization model for robust operation of EH to capture the positive features of both the IGDT and robust optimization models in handling uncertainties.Proposing a medium-term robust model for the optimal selfscheduling of the EH to enable the medium-term contracts in the problem for the first time.Although the medium-term contracts have been considered in the optimal self-scheduling of EH [15,16], the robustness of such contracts in EH has not been addressed so far.Integrating the proposed robust medium-term model with the P2G technology.Since the gas price is fixed during a day and it varies during a week, integrating the proposed model with the P2G system reduces the burden of the peak days when the prices of natural gas are high.

Paper organization
The rest of this paper is organized as follows.The main description of the problem is presented in Section 2. Section 3 gives the problem modeling.Section 4 provides the solution methodology.Numerical results are presented in Section 5, and finally, the conclusion is given in section 6.

Problem description
The EH is a concept to consider the multi-carrier energy systems.In this paper, a medium-term hybrid IGDT-RO is developed for the optimal self-scheduling of the EH.Four inputs from two energy types, including the wind turbine, DA electricity market, bilateral contracts, and the natural gas network, are taken into account in the EH, see Fig. 1.The outputs of the EH are the electricity, heat, and natural gas to satisfy the corresponding demands.The electricity demands are satisfied through the wind turbine, DA electricity market, bilateral contracts, electrical storage system (ESS), and a CHP unit.The electricity market can feed the P2G, ESS, and also the electricity demands directly, while the bilateral contracts and wind turbine generations are directly transmitted to the electricity demands.Among the three electricity resources, i.e., wind turbine, bilateral contracts, and electricity market, only the electricity market can be stored in the ESS.It is due to the fluctuation of the prices in the electricity market, while the other electricity resources have no price fluctuation.The CHP is the cogeneration unit that links the gas, electricity, and heat energy carriers.Meanwhile, the boiler is fed by natural gas and generates heat.The heat demands are met through the boiler and CHP units, while the natural gas demands are satisfied directly with the gas network and the P2G storage system.It is worth mentioning that the P2G is fed by electricity and discharged into the gas network.The EH operator is a large consumer that can directly participate in the electricity market.The boiler and CHP are also fed through the natural gas network.Due to the concerns regarding releasing the emissions and especially the CO 2 emissions, the cost of emissions related to electrical energy and natural gas are taken into account.The emission accompanies the electricity procured from the DA electricity market and the bilateral contracts as well as natural gas purchased from the gas network.In addition, the P2G storage system is fed by the DA electricity market and injects natural gas into the gas network.The DA electricity price and wind turbine electricity production are subject to uncertainties.The robust optimization-and IGDT-based models are two well-known methods to address the uncertainties.While there are advantages and disadvantages to these methods, this paper aims at capturing the positive features of the IGDT and RO models to achieve the worst-case realization of the uncertainties within a medium-term horizon.

CHP model
The electrical and thermal outputs of the CHP are mutually dependent.Their dependency is defined through an enclosed curve which is called a feasible operation region.Fig. 2 shows the feasible operation region represented by three lines that are described by a set of coordinates (P, Q).For instance, the coordinate of the point A is represented by (P A , Q A ). Eqs. ( 1 4) and ( 5) force the bounds for the amount of heat and electricity outputs and also assure staying at the right hand side of line AD.Considering all these equations simultaneously guarantees that the coordinates (P, Q) is definitely within the closed area of ABCD.

Demand response model
In this paper, electrical demand response is implemented to shift the peak hours energy consumption to the off-peak intervals.It reduces the cost of procuring electricity by the EH owner in the proposed model.It should be noted that the EH should keep the total electrical energy consumption constant while applying the demand response program.The electricity demand after applying the electrical demand response (EDR) is indicated in (6).The electricity energy shifting up and down are declared by (7)  ; dr E;down t !0 (10)

Uncertainty modeling
The uncertainty of the DA electricity market price is modeled by the RO approach.This uncertainty has to be maximized to consider the worst-case realization of the DA electricity price.The objective function in (11) aims at maximizing the deviation of the DA electricity price.The deviation is added to the expected DA electricity price in (12), while (13) declares that the deviation can be either positive or negative and the total value is limited by the uncertainty budget in (14).
The mentioned absolute values in ( 13) and ( 14) make the problem non-linear.Therefore, the uncertainty model needs to be linearized using the method presented in Appendix A. Then, the equivalent min problem of the objective function (11) is obtained in the following using the duality theory.Eqs.16e19 are the constraints of the dual problem.min x t !p DA t (16)

Main model
The aim of this paper is to minimize the total operation cost of the system in (20).The objective function consists of the costs related to the DA electricity market, purchased electricity from the bilateral contracts, purchased natural gas, P2G operation, and the CO 2 emission.Eqs.21e25 indicate the mentioned costs, respectively; the electricity, heat, and natural gas demands are presented in ( 26)- (28); the electricity purchased from the bilateral contracts are bounded by ( 29) and ( 30): the gas to heat conversion is performed by (31); and ( 32) and (33) indicate the generated electricity and heat through the CHP, respectively.The P2G technology produces natural gas through DA electricity.First, the electricity is converted to the equivalent natural gas via (34) and then it is charged in P2G storage system [21], while (35) stands for the energy balance in the P2G storage system.The state of charge of the P2G storage system is equal to the amount of charge in the previous hour plus the amount of charging in the current hour minus the discharging amount of natural gas to the gas network.The power injected to the P2G, the amount of discharging gas from the P2G storage system, and the state of the charge are limited by ( 35)-(38).The energy balance for the ESS is presented in (39).The ESS state of the charge and charging/discharging electricity energy to the ESS are bounded in (40)-(43), while (44) ensures that only one of the states of charging/discharging occurs at a time.Finally, (45) shows that the first part of the purchased electricity from the DA Fig. 2. Feasible operation region for the CHP unit [16].
electricity market is transferred into the output directly, and the second part is stored in the ESS and the third part is injected into the P2G system.max Subject to (1)e( 10), ( 12)e( 14), and, p BC bt !p BC b;min ,u BC bt (30) 4. Solution methodology

Solution steps
The volatile nature of the DA electricity price and wind turbine generation cause a challenge for the EH manager to make appropriate decisions to schedule the operation of the EH.Therefore, it is crucial to deploy a useful approach to tackle the uncertainties.The IGDT and robust optimization techniques are two commonly-used approaches to increase the robustness of uncertainty-based problems.These methods, besides their merits, have some shortcomings.The IGDT often converts the problem to a non-linear problem [25] resulting in difficulties for commercial solvers to find the optimal solution.Since the RO techniques are often applied to convex problems, increasing the binary variables in a problem makes a limitation in deploying the RO in the problems with various binary variables [24].This paper presents a hybrid IGDT-RO model aiming at taking advantage of both techniques in the process of maximizing the inner and outer layer of the (46) (as the uncertainty layers) and minimizing the middle layer (as the cost objective function).To tackle the uncertainties and the mentioned deficiencies, the Hybrid IGDT-RO is described within the following steps.
Step 1: linearizing the nonlinear terms of the proposed robust optimization model for the DA electricity price uncertainty presented in ( 12)-( 14), using the method introduced in the Appendix A.
Step 2: transforming the linearized version of the max problem of the DA electricity price mentioned in (A.1)-(A.5)to a min problem using duality theory, see (15) Step 4: applying the IGDT method to Step 3 to take into account the wind turbine generation, see (56)e(58).

Information gap decision theory
The IGDT approach is deployed to model the robustness of the problem against wind generation.The IGDT is useful when there is a high level of uncertainty in the problem [26].Based on the decision maker's opinion, the uncertainties can be treated to have more cost/less profit or less cost/more profit.The IGDT considers these conflicts using two concepts, including robustness and opportunity functions [27,28].Since the robust version of the IGDT is utilized here to reach the worst-case realization for the wind power uncertainty, the opportunity function is not required to be considered.The uncertainty in the IGDT model is presented as follows [29,30].
where (47) declares that the value of uncertainty in the IGDT, r, falls within the bounds ð1 ÀaÞr and ð1 þ aÞr .
In the robustness function, the maximum uncertainty radius (a) is determined by dedicating a specific risk value to the base deterministic mode of the objective function.In other words, the objective function is robust against uncertain input parameters, and the amount of objective function is lower than the specific value within the obtained radius of the uncertainty.These decisions are made by risk-averse decision-makers to reach the worst-case realization.The mathematical representation of the robustness function is described as follows.max x a (48) where (48) maximizes the uncertainty radius with respect to the decision variable x.The sets of equality and inequality constraints are presented in (49) and (50).According to (51) and ( 52), for a minimization problem, the robustness target f s ðx; r Þ should be lower than a critical value, which is calculated by the risk level and the base objective function.The base objective function is obtained using the expected values as input for solving the problem.The uncertainty is bounded by (53) within either the negative or positive uncertainty radius, while (53) bounds the risk level.

Hybrid IGDT-RO model
In order to deploy the IGDT model for the wind power uncertainty, the uncertainty radius of the wind power generation is described as follows.10), ( 16)e( 19), ( 26)e(45) and, p w t ¼ ð1 À aÞ,p w t (57) The objective function of the Hybrid IGDT-RO is presented in (56).The worst-case for the wind power in a minimization problem means decreasing from the forecasted wind power value.In other words, the worst-case occurs for the wind power in (57).The constraint of the critical cost is given in (58).The left side includes two parts.The first part constitutes the bilateral contracts cost, natural gas and P2G operation costs, and the emission cost.The max second term is the dual objective of the DA electricity market cost, which imposes the worst-case realization of the DA electricity price, and it is obtained from Step 2. The right side declares the maximum deviation from the expected (based) cost.For more comprehension, Fig. 3 presents the flowchart stating the procedure of the proposed hybrid IGDT-RO framework, where the steps of the RO and IGDT are determined separately.

Case study
The proposed hybrid model is tested on the aforementioned EH presented in Fig. 1 containing three inputs and three outputs.The DA electricity market prices and wind turbine generations are subject to uncertainty.The electricity market prices in the mediumterm and long term are modeled by autoregressive integrated moving average (ARIMA) time series [31].The wind speed is also can be modeled by autoregressive and moving average (ARMA) time series [32].The data related to the electricity market prices and wind turbine characteristics and generations during four weeks, equal to 672 h, is obtained from Ref. [16].The natural gas prices are taken from Ref. [33].Fig. 4 shows the DA electricity market and natural gas prices.The electricity prices are changed hourly, daily, and weekly while the gas prices only vary daily.The electricity prices change from a minimum value of 12.15 $/MWh in hour 557 to a maximum value of 39.6 $/MWh in hour188, while the minimum and maximum values of the gas prices are 20 $/MWh and 29.16 $/MWh on days 3 and 9, respectively.
The wind generation is depicted in Fig. 5 that varies from 0 to 40 MW, and the average electricity production through the wind turbine within 672 h is 33.87 MWh.Fig. 6 shows the electricity, heat, and natural gas demands.The electricity and heat demands are adopted from Ref. [16].For creating the gas demands, the presented short-term demands in Ref. [21] has been extended to four weeks using the normal distribution function.The mean values are assumed like the short-term demands, and the standard deviation is considered 3%.The electricity demands range from 73.9 MW to 129 MW during the 672 h, while the natural gas demands change within a narrower interval of 74.2 MW and 106 MW.The minimum and maximum values of the heat demands are 22.27 and 119 MW.
Table 1 provides six candidates for the bilateral contracts, including the contract prices and minimum/maximum values of the electrical energy that could be purchased.Two first candidate contracts cover the whole period, and the other ones could be purchased for part of the periods.The maximum generation of the boiler is equal to the maximum heat demands with 85 % efficiency.Four (heat, electricity) coordinates of the feasible operation region of the CHP unit with respect to MW are (0,35), (25,25), (20,5), and   (0,10), considering 0.45 and 0.4 for the heat and electricity conversion, respectively.The other details of the technical specifics of the P2G and ESS are given in Table 2.The emission coefficients are taken from Ref. [22].To assess the proposed model properly, two studies are conducted.In the first study, the N p , risk level, and EDR are assumed to be 336, 0.1, and 10 %, respectively.The sensitivity analysis on the characteristics of the proposed hybrid model and the demand response is done in the second study.Finally, taking into account all the aforementioned data, the MILP model is carried out with CPLEX solver under GAMS environment version 34 [34].

Discussion
First, the results of the dispatches and purchased energy carriers are presented.For to sake of clarity and reducing the line congestion, the results of the second week are shown in the figures.The analysis of the other weeks is similar.The results are presented in three types of solutions.The expected solution gives the results without considering the uncertainties.Having the results of the expected solution provides an opportunity to compare the results with and without considering the uncertainties.In our analysis, the term Robust means solving the problem considering the robust optimization only for the DA electricity price uncertainty, while the term Hybrid means considering both uncertainties, i.e., robust optimization for the DA electricity price and IGDT model for dealing with wind power uncertainty.Fig. 7 demonstrates the robust uncertainty sources and the expected values.The DA electricity prices increase in the robust solution to reach the worst-case realization, while the wind generation decreases for the same goal.Fig. 8 depicts the purchased electricity energy in three solutions.More electricity is purchased from the DA electricity market in the expected mode.It is because the DA electricity price is lower in the expected situation compared to the robust and hybrid solutions.Since the generated electricity by wind turbine reduces in the hybrid model, the amount of electricity purchased from the DA market in the hybrid model is more than the robust model.Fig. 6.Electricity, natural gas and heat demands [16,21].Therefore, the purchased electricity from the DA electricity market raises to compensate for the lack of wind generation in the hybrid solution.
The results of gas entering the EH are given in Fig. 9, where it declares that more natural gas is procured in the robust and hybrid solutions.That is, at least one part of the shortage of electricity purchased from the DA electricity market in the robust solution is compensated through the natural gas using the CHP unit.It can be verified via Fig. 10-a.For instance, for the expected state, from hours 169 to 240, no electricity is generated by the CHP, while a significant amount of electrical energy is produced in the robust and hybrid states during the same period.As mentioned in Fig. 8, more electricity is purchased from the DA electricity market in the expected solution during the range of hours 169e240.It is due to the lower expected electricity price compared to the robust electricity prices.It means, since more electricity is purchased from the DA market, there is no need to generate electricity by the CHP in the expected solution during this range.On the other hand, since less electricity is purchased from the DA market in the robust and hybrid solutions (compared to the expected solution), the lack of required electricity is compensated through the CHP during the time interval 169e240.Meanwhile, the amount of electrical energy generated in the hybrid model shows a higher amount in comparison with the robust model.As mentioned before, there is a reduction in wind turbine generations in the hybrid model.Therefore, other devices have to compensate for the lack of electrical energy to meet the electricity balance equation.It should be noted that the heat generated through the CHP has a similar pattern to the CHP electricity output due to the mutual dependency of the heat and electricity in the CHP unit.
The results of procuring electricity from the bilateral contracts in the hybrid solution are provided in Fig. 11.No electricity energy has been purchased from contracts 1 and 6 due to the higher prices of these contracts.On the contrary, the transaction of contract 5, as the cheapest contract, is fully conducted.Less electrical energy is purchased from contract 4 compared to contract 3. It is due to the fact that their validation times are different.Contract 3 is for the first two weeks and contract 4 is valid for the last two weeks.The charging/discharging state of the ESS is given in Fig. 12 where the negative part shows the discharging state and positive values stand for the charging state.For a period of one week, it indicates seven charging and seven discharging states.The charging occurs in offpeak hours when the electricity market price is low and the discharges take place in peak hours.For example, it charges between hours 194 to 197 and discharges from hours 207 to 212, which are corresponding to the off-peak and peak hours electricity prices, see Fig. 4. It is worth mentioning that since charging in off-peak intervals and discharging in peak intervals are always beneficial, the pattern of the expected, robust, or hybrid models do not experience a significant difference.
It can be seen from Fig. 13 that more charging/discharging natural gas occurs in the expected state.The input charging natural gas is obtained using the purchased electricity from the DA market.In both the robust and hybrid models, the DA electricity prices increase to reach the worst-case realization.Therefore, when the electricity price goes up, less natural gas is charged into the P2G storage system.The other feature which can be captured from the figure is that the P2G storage system is charged in the off-peak hours of electricity prices and discharges in natural gas peak intervals.For instance, it charges from hours 194 to 197, which is equal to the off-peak electricity market prices, and discharges from hours 198 to 215, which is the peak hours of the natural gas prices.It means that the P2G system can help to supply a part of the natural gas entering the EH.For better evaluation, Fig. 14 shows the effect of P2G on the reduction of the total gas entered the EH.Obviously, in some hours, the amount of purchased natural gas decreases significantly, especially in the times between 146 and 220 and also 384 to 435, which are the peak intervals of the natural gas prices in Fig. 4.
For more comprehension, the average values of energy dispatching of the CHP and boiler, the ESS and P2G storage system as well as average values of procuring energy from the natural gas network and purchasing electricity from the bilateral contracts and DA electricity market are provided in Table 3 in all periods.This table verifies the mentioned results for one sample week.For example, the amount of natural gas purchased from the network and gas entered the CHP experience an increase in the RO and hybrid model.Similar to the purchased gas, the procured electricity from the bilateral contract raises from 1.57 MWh to 3.2 MWh and 3.3 MWh for the robust and hybrid approaches, respectively.On the other hand, the amount of electrical energy purchased from the DA electricity market shows a decrease of about 21.14 MWh in the RO solution and an increase of about 9.65 MWh in the hybrid model to  compensate for the lack of wind turbine generation.Meanwhile, the operation mode does not affect the charging/discharging pattern of the ESS and consequently, it does not change the amount of charging/discharging values via different solution approaches.Finally, the operation costs are provided in Table 4. Achieving the worst-case realization considering the DA electricity price uncertainty increases the gas network operation cost from 2,366,226 to 2,532,198 $ while it decreases the electricity market operation cost from 858,839 to 630,263 $ in the robust solution.It does not have a remarkable effect on the reduction of the emission cost due to replacing the DA electricity with a source with similar emissions in the robust solution.Since wind power (as clean energy) reduces in the hybrid mode, the emission cost goes up.By comparing the robust optimization solution and hybrid solution, which includes both sources of the uncertainties, it can be found that the cost increases due to the decrements of the wind turbine free-cost generation.Totally, the proposed hybrid IGDT-RO model shows the total cost of 3,799,229$, which is the highest cost among all since it considers the worst-case realization for both the DA electricity price and wind turbine generation.

Sensitivity analysis
This section performs a sensitivity analysis on the main  characteristics of the proposed hybrid IGDT-RO model and EDR to properly assess the model.Fig. 15 shows the results of evaluating the effectiveness of the risk level (b) on the uncertainty radius and total operation cost.The higher risk level means the higher values for both the uncertainty radius and the total cost.For example, a risk level equal to 0.04 means the objective function is allowed to increase by 4 %.In this situation, the uncertainty radius reaches 0.238.In other words, the wind generation decreases by 23.8 %, due to Eq. (57) to reach the 4 % risk level for the objective function.
Figs. 16 and 17 demonstrate the sensitivity analysis on the number of hours allowed to affect the DA electricity prices.The former one (Fig. 16) shows the variation of the total discharging gas   from the P2G storage system, total generated electricity through the CHP, and the total purchased electricity from the bilateral contracts.The latter figure (Fig. 17) depicts changing natural gas operation cost, the DA electricity procuring cost, and the total operating cost.Fig. 16-a shows that by increasing the N p , the total discharging gas from the P2G system decreases.It is because the number of affected DA electricity market prices increases, and the DA electricity market feeds the P2G storage system.Therefore, it is not an economic decision to charge the P2G system from an expensive source.On the contrary, according to the Fig. 16-b and Fig. 16-c, the total procured electricity from the CHP and the bilateral contracts increases.In order to justify this, it should be noted that the natural gas and bilateral contracts are more reasonable sources compared to the DA electricity market, especially when the DA electricity market prices increase in the worstcase.It is because the bilateral contracts and natural gas prices are lower compared to the DA electricity prices in the worst-case.Meanwhile, the trend of the total cost of the natural gas network operation in Fig. 17 confirms increasing the gas entered the CHP.Moreover, the decrement in the total cost of interaction with the DA electricity market by increasing the N p means that the proposed model prefers to procure the required electricity from cheaper sources such as the CHP and the bilateral contracts, instead of buying electricity from the DA electricity market.Fig. 17-c shows that on the whole, a higher total cost is imposed on the system by increasing the number of hours subject to uncertainty (N p ) in the DA electricity prices.It is noteworthy to mention that the curves in the mentioned figures are flattened after hour 500.This happens since the number of defined uncertainty budgets limits the amount of number that they can affect the problem.As mentioned in (13), there is a bound for the deviation of each hour, while the total deviations are limited by the uncertainty budget.In other words, the total amounts of the deviations reach the defined budget.
Table 5 provides the effect of changing the EDR on the natural gas and DA electricity energy purchasing costs.The EDR ranges from 0 to 15 %.As can be seen from Table 5, the EDR is significantly effective in reducing the DA electricity purchasing cost from 813,533 to 790,969$, which is equal to a 2.88 % reduction in the DA electricity purchasing cost.It also reduces the total operation cost from 3,825,441 $ to 3,786,332$, however, its effect on the natural gas purchasing costs is not remarkable.

Conclusion
In this paper, a hybrid IGDT-RO framework has been proposed to investigate the robust operation of an EH within a medium-term time horizon for a large consumer.Since the problem deals with the uncertainties of the wind generations and DA electricity prices, the developed model maximizes the DA electricity prices deviation using the RO approach, and then the IGDT maximizes the wind generation uncertainty radius when the total operation cost of the EH has been incorporated in the IGDT part.The main results of the paper are as follows.
The robust electricity prices rise most of the time, and consequently, the procured electricity from the DA electricity market is reduced.On the other hand, the electricity purchased from the bilateral contracts and generated through the CHP unit increases to compensate for the lack of purchased electricity from the DA electricity market.Therefore, the cost of procuring energy from the bilateral contracts and natural gas network increases in comparison with the expected solution.On the whole, the total operation cost of EH increases by 8.6 % in hybrid solution compared to the expected one.The emission costs go up in the hybrid IGDT-RO model due to the decreasing generation of the emission-free wind turbine generator.The P2G storage system purchases electricity from the DA market in the off-peak intervals and discharges into the natural gas network during the peak intervals.The ESS always charges in off-peak hours of the DA electricity price and discharges in the peak hours without dependency on the operation mode.The sensitivity analysis reveals that the higher risk level results in more operation costs and more uncertainty radius for the wind turbine generation.Changing the uncertainty budget of the DA electricity price shows that the problem is robust against different levels of the DA electricity price uncertainty.Sensitivity analysis on the EDR demonstrates that the EDR reduces the DA electricity market operation cost significantly, in   which the total cost of purchasing electricity from the DA market reduces by 2.77 %, in case 15 % EDR is included.
As the prospect of future works, the proposed hybrid IGDT-RO can be further developed to reach the worst-case realization of various kinds of uncertainties in the operation of power system problems with a large number of binary variables due to the operation conditions of assets such as energy storage systems, distributed generators, ramp up and ramp down of generators, etc.

Declaration of competing interest
The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper.
)e(5) indicate the constraints of the CHP feasible operation region.Eq. (1) indicates the coordinates (P, Q) stay below line AB.(2) and (3) show the coordinates upper the BC and CD lines.( 1 to 4, the hybrid IGDT-RO model is obtained as follows.

Fig. 9 .
Fig. 9. Natural gas purchased from the gas network.

Fig. 16 .
Fig. 16.Effect of the N p on the generation of the resources.

Fig. 17 .
Fig. 17.Effect of the N p costs.
Conversion coefficient of the boiler unit g e,CHP , g h,CHP Electrical/thermal coefficient of the CHP unit g P2G Dual variables of the inner objective.a x t ; d t ; d t ; s t ; 6 (19)). max

Table 1
Candidate bilateral contracts data.

Table 2
Details of the ESS and P2G systems.

Table 3
Generation/dispatching mean values of the resources (MWh).

Table 4
Operation costs in different solutions ($).
Fig. 15.Effect of the risk level on uncertainty radius and total cost.

Table 5
Effect of DR on the operation costs of natural gas network and the electricity market.