A CVaR-Robust-Based Multi-Objective Optimization Model for Energy Hub Considering Uncertainty and E-Fuel Energy Storage in Energy and Reserve Markets

The increasing demand for energy carriers has expanded the use of energy hubs that employ distributed demand response programs to improve power system reliability and efficiency. Moreover, the unstable behavior of renewable resources, as well as the indeterminate electrical and thermal demands, create major problems for energy hub operation. Inspired by this, this paper presents a day-ahead scheduling framework for energy hubs (EH) in energy and reserve markets considering two main objectives of economy and pollution emission. The studied energy hub consists of a novel hybrid energy storage facility based on a fuel cell, wind power, photovoltaic energy, and a particular fuel cell unit in the presence of elastic demand. This multi-component system participates in energy and reserve markets as a single entity to optimize energy hub operation. The proposed method also models the uncertainty of wind speed, photovoltaic irradiance, and load using the Mont-Carlo method. The energy hub risk level is analyzed using the conditional value at risk (CVaR) approach to increase the EH operation and efficiency. The proposed multi-objective energy hub model is solved using the MINLP method in General Algebraic Modeling System (GAMS) to minimize operation cost and pollution emission. Finally, to prove the effectiveness of adding a new E-fuel energy storage system and considering uncertainties on energy hub operation, the proposed method is compared with other reported models.

Efficiency of the turbine ρ Density of water g gravity k Step size H p pumping head Q p pump discharge η p efficiency of RPT κ PV production price of PV κ e t price of supplied power via upstream grid κ gas t natural gas price κ wind production price of wind turbine κ wa water price κ DRP demand response price κ TDRP thermal demand response price κ WDRP water demand response price κ GDRP Gas demand response price κ TSS Price of produced thermal by TSS κ E−fuel Price of produced power by E-fuel κ ESS Price of produced power by ESS κ DG Price of produced power by DG σ, ω and γ d cost factor of DG π, θ & ϑ emission rates of DG CO is equal to 1 if the CHP unit is active otherwise is equal to 0 are equal to 1 if the related zone is active otherwise, are equal to 0.

A. MOTIVATION
The energy hub is a new technology with high efficiency and reliability [1], [2]. An energy hub is a subsystem of the power system in which conversion, production, storage operations, and consumption of various energy carriers are performed simultaneously, unlike other power system structures [3]. In other words, an energy hub is a unit containing inputs, outputs, conversation, and storage with the use of multiple energy carriers [4]. The main energy hub structure includes three parts, made of input units such as natural gas, distributed generation, wind turbine, photovoltaic, fuel cell, and converter units such as boiler and CHP, as well as a different type of storages to ensure optimized operation of the energy hub system [5]. The main advantages of the energy hub system can be classified as minimizing operating costs and pollution emission, enhancing reliability, supplying demand response, and assisting energy management [6].
• Minimizing operating costs and emission pollution; VOLUME 9, 2021 • Enhancing reliability, supplying demand response, and energy management. Exploiting new methods to achieve optimal performance in EH scheduling is of utmost importance. One such method is the demand response program [7]. Due to the scarcity of non-renewable units and to reduce the pollution emission caused by these resources, employing renewable energy resources such as wind turbines and solar system is a good alternative in the management of an energy hub integrated with the energy and reserve market [8], [9].

B. LITERATURE REVIEW
The usual power systems have a simple structure so that electricity is produced in central power plants and transported with a long transmission line to supply demand. However, it has some disadvantages such as great operational cost, low efficiency, high emission problem, and some control problem. Also, in this structure, all types of equipment are operated independently of each other [10]. The researchers introduced a new structure called an energy hub to solve the mentioned problems, [3]. In the EH, economic, technical, and environmental models in multi-carrier energy systems have been investigated.
One of the first and most essential steps in designing an energy hub is determining the needed equipment and the optimal size of them, etc. In Ref. [11], the optimal size of the equipment in an energy hub structure for robust scheduling were analyzed. only non-renewable sources such as gas and electricity without storage sources are considered. Different kinds of renewable sources were applied to provide demand response and enhancing the EH operation by ignoring the pollution emission problem in [12]. To optimize the scheduling of the energy hub, a novel demand response program was represented in [13]. On the other hand, integrating the renewable sources such as wind turbine, solar system, etc. with hub energy have a significant impact on the hub energy performance. However, due to the uncertain nature of these sources, this impact has been investigated in the literature [14], [15]. Moghaddas-Tafreshi et al. [16] suggested a novel method for optimal operation of an energy hub integrated with RES, PHEVs, a fuel cell, and a thermal storage system. A new method was represented to forecast the power consumption uncertainty using information gap decision theory (IGDT) under risk-averse and risk-seeking strategies. An optimal load dispatch method was reported for a community energy hub considering both thermal and electrical demand response programs integrated with a robust optimization method to optimize the price uncertainty by Ma et al. [17]. In Ref. [18], a mixed-integer linear programming model was defined for energy hub considering operational limitation. The proposed method was implemented on the different working models such as dispatch, synchronization, de-synchronization, and soak. But in this study, DRP has not been implemented on and the effect of the using energy storage was not analyzed, too. Khodaei et al. [19] suggested a multi-objective method based on the proposed cost and pollution function. Using the trade-off method based on a fuzzy decision-making technique and using peak load management, operation cost and emission problems are minimized. Tian et al. [20] studied the stochastic procedure of energy hub in the presence of uncertainty parameters such as wind turbine speed, electricity market price, and load. Then downside risk limits were proposed to optimize the risk-in-cost considering the electricity market and thermal DRP. A nonlinear power flow was implemented in a simple energy hub include a heat exchanger, CHP, a furnace, and a storage unit to optimize the operation cost in [21]. As the CP unit and furnace are generated pollution emission so it is necessary to analyze this problem. He et al. [22] are applied the robustness and opportunity functions of IGDT to investigate the effect of load and market price uncertainty on the performance of simple hub energy. In [23], stochastic operation of energy hub was reported with considering the uncertain environment and also downside risk restraints to minimize the risk-in-cost where heating market and heat demand response program (DRP) was considered. Miao et al. [25] proposed CCHP based microgrid in the presence of renewable energies for hub energy, and then the environmental issue was considered as a bi-objective model. Compromise programming and a fuzzy approach were used to solve the proposed model. A real-time demand response program was applied to reduce operation costs and emissions. As renewable sources have been an uncertainty factor so it is a good option for considering the uncertainty factor and analyzing the risk level by CVaR method.
The other important option to improve the performance of hub energy is using a suitable method to solve energy hub functions. In [25], the performance of an energy hub integrated with the wind turbine, and a novel storage and demand response in the distribution network was optimized using two objective functions including operation costs, pollution emission considering the uncertainty of wind speed, price, and load. In this paper, the main purpose is reliability increasing in addition to reducing costs and pollution problem, so it would be better to use risk-aversion method to improve the reliability of the system against potential hazards. Optimal scheduling of smart residential energy hub (SREH) by using a compressed air energy storage under market prices, solar radiation, and DRP uncertain parameters was characterized by applying a risk-constrained two-stage stochastic scheduling model [26].
Ref. [27] reported an optimal bidding strategy model of energy hub by employing specific features such as multi-disciplinary and flexibility in the energy market. Besides, a stochastic model was represented considering wind energy. The model considers the uncertainty of prices both in day-ahead and real-time markets. In this work, wind turbine is the only renewable source while the using of PV source as well as fuel cell can be a good option. Pakdel et al. [28] studied the impact of heating market and heat demand response program on energy hub operation under price uncertainty. In [24], the researcher investigated the concept of the hydrogen-based smart micro-energy hub (SMEH) considering integrated demand response (IDR) and fuel cell-based hydrogen storage system (HSS). According to all of the mentioned researches, the energy hub plays an important role in supplying energy of the power system. It able to minimize all the mentioned problems such as economic cost, pollution emission problem and increasing reliability. In [29], the CVaR method was implemented to model the potential risk of SMEH scheduling cost as a constraint and a stochastic model was used to model the wind generation and load demand uncertainty. But other uncertain parameter and demand response program were not considered. In [30], a comprehensive multi-objective model was reported to minimize both the energy procurement cost and risk level in energy hub. For controlling the pernicious effects of the uncertainties, conditional value at risk (CVaR) method was used as risk management tool in energy hub performance. The proposed model was formulated as a mixed integer nonlinear programming (MINLP) problem and solved using GAMS.
Based on the literature, there are some drawbacks to the use of renewable sources, fuel cells, novel hybrid storage, and the implementation of various demand response methods for gas, water, electrical, and thermal loads. Important issues include the investigation of uncertainty factors with the most accurate methods and analysis of risk level and its impact on energy hub performance, and considering the impact of reserve market participation on system costs and pollution. Moreover, the EHs participate in the reserve market to cope with their mismatch between day-ahead offering and real energy production/consumption providing for the power market. The reserve market is designed to balance the energy generation and consumption incautiously. In fact, purchasing energy from the reserve market at a lower price increases the reliability of the system, reduces emission pollution and costs. The reserve market is considered as an active storage. Finding an accurate optimal operation method to achieve the best result is another major issue.

C. CONTRIBUTION
The proposed structure of the energy hub is illustrated in Fig. 1. In this work, the effect of E-fuel energy storage on the energy hub's optimal operation in the energy and reserve markets is examined. The optimum performance of the energy hub system is investigated from both economic and emission views in the presence of EDRP, TDRP, GDRP, and WDRP to reduce the operation cost and enhance environmental performance. Furthermore, the conditional risk at value (CVaR) method is adopted to describe the uncertainty of wind speed, solar irradiance, and load in objective functions combined with the robust optimization theory. Renewable and non-renewable sources and a fuel cell have also been utilized simultaneously. The reserve market cooperates with the energy market to deal with the unsought faults and events, which is unable to supply loads or purchase with a high price. The augmented ε-constraint and fuzzy methods are also adopted to solve the proposed multi-objective problem. The main contributions of this study are listed as follows: • A novel storage model as E-fuel energy storage inside TSS and HSS is used in the energy hub structure to optimize energy management; • The CVaR method is introduced to describe the risk level of the energy hub performance, in addition to considering the uncertainty of the wind speed, solar irradiance, and load; • The reserve market and energy market are considered for load demand satisfaction.
• The multi-objective problem is solved by a hybrid multi-step method based on the augmented ε-constraint and fuzzy methods;

II. PROBLEM DESCRIPTION
The EH is depicted in Fig. 1 consists of a wind turbine, upstream grid, water, and natural gas sources to supply the demand for electricity, cooling, heat, water, PV, and fuel cell as inputs sources. Moreover, the CHP and boiler are considered as a converter in interior devices, and a novel ESS system is equipped for each end consumer. The EH can be employed as both states of supplier and consumer, selling energy if there is an extra capacity from the EH, and purchasing energy when the EH system is not activated or when energy prices are low. The installation of a novel ESS system helps save the energy generated from the EH devices. Part of the gas purchased from the network is used to supply the CHP and boiler, which is then converted into electricity and heat. The purchased water is utilized to meet the water demand of residential areas. Thermal, electrical, and a hybrid fuel cell storage is considered as the storage system to provide energy in critical conditions. Moreover, TDRP, GDRP, WDRP, and energy DRP are used to flatten the load curve and minimize the pollution emission and total operation cost. The uncertainty in the day-ahead framework is taken into account to enhance the reliability and accuracy of the system. Moreover, the CVaR method is adopted to cope with unwanted fault and manage the risk level of the hub system. The main optimization objectives considered in this paper are economic cost and pollution emission cost. All of these issues are classified, as shown in Fig. 2.

A. MODEL FORMULATION
The following section describes the model of energy hub elements including solar, wind, fuel cell, and the proposed energy hub models.

1) PHOTOVOLTAIC SYSTEM
The electric power generated via the PV in each time step is formulated based on the hourly solar irradiance flux as follows [31]. where: where, Tnoc and TRef are 43 • C and 25 • C, respectively. Also, is 4.5 × 10−3 • C for mono-Si with η ref ≈ 0.12 and Tref ≈ 25•c [32].

2) WIND TURBINE
The energy generated by the fuel cell in each time step t∈ {1, 2, . . . , NT} is formulated based on the hydrogen consumption in peck time as follows [34]:

3) FUEL CELL
The energy generated by the fuel cell in each time step t∈ {1, 2, . . . , NT}, H FC 2t , is formulated based on the hydrogen consumption in peck time as follows [34]: Several constraints related to fuel cell performance are: Eq. (6) is defined to cope with the charging and discharging processes simultaneously.

4) E-FUEL ENERGY STORAGE
Herein, a new energy storage is introduced in the energy hub structure that operates based on an E-fuel charger and an E-fuel cell. The E-fuel energy storage performance is similar to the hydrogen storage system's performance and, in addition to the advantages of HSS, it has a highly safe and easy transfer and storage. Its schematic representation is plotted in Fig. 3. It consists of a fuel cell and an E-fuel charger. The E-fuel charger is charged electrically to save the extra electricity obtained from wind turbines and solar photovoltaics in liquid fuel state, while the E-fuel cell can produce electricity in the high-demand condition using the energy saved from the e-fuels. External reservoirs are used to save the charged and discharged E-fuels and by trucks, pipelines, or ships which are used in E-fuel storage, it can be transported. The easy storage and convenient transfer with E-fuel make the E-fuel system an exceptional choice for both off-grid and on-grid scales for power supplies for an energy hub. One of the main advantages of the proposed storage is that the charge and discharge mode can be operated simultaneously unlike all current electrical energy storage systems because the E-fuel charger and the E-fuel cell function independently [35].
The hydrogen adsorption energy of E-fuel and the Gibbs free energy for hydrogen adsorption are given by: The Eq. (8) can be approximated to be: The CHP unit is the other important element of the EH. It plays an important role in the hub because it can integrate all the electricity, heat energy, and natural gas carriers. The CHP unit has an applicable operation zone, as shown in Fig. 4. The action of the CHP operation zone income is as follow, with the produced heating, the CHP power production is restricted and the output heating of CHP can be limited by the generated electricity. This non-convex zone makes the entire problem of optimal operational planning a non-convex problem. Therefore, the overall optimization of these types of problems cannot be assured. In this paper, the convex operation zones for this unit are considered by a novel triple division. These operation zones are illustrated in Fig. 4 by Zone 1, Zone 2, and Zone 3. Convexity is performed by including three binary variables for the CHP. Notice that these zones are convex, meaning that each point on each line in any zone belongs to that zone [36].

B. OBJECTIVE FUNCTIONS
Economic and environmental aspects are two important factors for optimizing the performance of the energy hub operation, which are calculated below.

1) ECONOMIC OBJECTIVE FUNCTION
The first objective function is to minimize the cost of the EH and its cooperation with the energy and reserve markets, as given by: where the first and second terms represent the operation cost of the EH, and the cost of the reserve market, respectively. The operation cost of EH is calculated as follow: (11) where the Cost of each component is given by:  The cost of reserve market is calculated as follow:

2) EMISSION OBJECTIVE FUNCTION
Herein, three types of polluted gases are considered (CO 2 , SO 2 , and NO 2 ), which are produced by energy hub units. Thus, the pollution emission problem is another important factor considered as the second task in the operation of the energy hub, and expressed as: The total amount of the EH emission is represented by OF 2 . The emission of the imported power from the upstream grid, CHP unit, DG, boiler, and consumed residential gas is represented by EM grid , EM chp , EM DG , EM bo , and EM rl , which is formulated by Eqs. (13-1)-(13-6), respectively.

C. LIMITATION CONDITION
To optimize the operation of the EH, it is necessary to meet the limitations of power demand, gas, water, and thermal balance, storage operation, load variation, and system reserve.

1) ENERGY BALANCE CONSTRAINT
To optimal scheduling operation of an energy hub, the demand of energy in each time step must be equal to the total produced energy as follow: 2) THERMAL BALANCE CONSTRAINT The thermal balance constraint of EH must be equal to the total heat produced by the boiler, CHP, and the stored energy source as modeled by Eq. (15).

3) GAS BALANCE CONSTRAINT
The imported natural gas from the upstream grid is divided to use in three carriers including boiler, CHP, and load which is formulated as follow: It should be noted that a nominal range is defined for total purchased gas from the upstream grid.

4) WATER BALANCE CONSTRAINT
The amount of imported water from the upstream grid must be satisfied the demand load in each time steps: It should be noted that a nominal range is defined for total purchased gas from the upstream grid.

7) MODEL OF BOILER
Here, a boiler is considered a thermal resource, and a nominal capacity is defined for the boiler operating range.

8) MODEL OF CHP SYSTEM
The CHP operation is considered with 3 operation zones, and is nonlinear in some points of heat and electricity generation. Part of the extracted gas from the gas network is allocated for the consumption of CHP units to produce electric power and heat. Since the thermal production by the CHP system is often compatible with its electrical production, thermal operation restriction will be guaranteed if electrical operation restriction is ensured. To simplify the nonlinearity of the CHP unit performance, 3 operation zones are defined according to Fig. 4. constraints related to the CHP with three operational zones are defined as follows: The big-M method is implemented for mathematical formulation. The maximum and minimum sufficiently are represented byM and M which they are inactive at the anticipated time.

9) TSS MODEL
Herein, various types of energy storage are exploited to minimize energy losses, including TSS, ESS, and E-fuel energy storage. The stored heat level of TSS is modeled by Eq. (27).
Constraints of saved heat level, input and output heat amount are through TSS are formulated through: Heat waste of TSS is given by: Because the charging and discharging operation of the storage device cannot be done at the same time, thus:

10) ESS MODEL
Another critical issue in energy hub structure is to prevent waste of the produced energy through the grid. So, electrical storage is exploited in this part to save extra produced electricity in the grid. The mathematical model and constraint of ESS are expressed as [20].
Constraints of saved electric level, the amount of input and output electrical through ESS are formulated as:

11) E-FUEL STORAGE MODEL
The E-fuel energy storage is operated by saving the extra produced water of the grid and then converting water to energy, which increases hub energy reliability and minimizes the operation cost. The mathematical model and constraint of E-fuel are expressed as [21]: Restricting the saved heat level, the amount of input and output heat through E-fuel storage is given by:

D. DEMAND RESPONSE PROGRAMS
To satisfy the demand load in EH, it is necessary to implement a strong method for managing the load profile. Therefore, herein, a demand response program is adopted to satisfy the electrical load. An efficient method in demand response programing is the time-of-use (TOU) method that is applied here [38], [39]. In this method, the load consumption hours are analyzed, and then part of the load is transferred from the peak consumption hours to the low-consumption hours to minimize cost and pollution emission. Thus, EDRP based on TOU is formulated as follow: The constraint to prevent the process of simultaneous increasing and decreasing of the load level is modeled as: Also, to prevent the process of simultaneous increase and decrease of the thermal, water, and gas load level the following constraints should be satisfied: In the EH structure, the output power of the PV, WT, and load demand is stochastic, but EH scheduling is performed before obtaining the actual output of the WT and PV. Consequently, in day-ahead scheduling, a scenario simulation method is utilized to identify the day-ahead forecasting output result of the WPP and PV to overcome the uncertainty of WT, PV, and load demands. In this study, three factors of load demand, wind turbine speed, and PV irradiance have uncertainty.

1) LOAD UNCERTAINTY
The energy market includes three kinds of load: adjustable, interruptible, and sensitive. The sensitive load is hard to involve in the demand response and is also relatively fixed. On the other hand, the other two types of loads are flexible, causing strong uncertainty. In this study, load demand is divided into two parts, with and without uncertainty, as follows [36]: Note that load uncertainty and load are modeled by a normal distribution L u t ∼ [0, δ 2L,t ] function [36], and L t ∼ [L e t , δ 2L,t ] function, respectively. Thus, using the Monte Carlo method, 1000 scenarios with the normal distribution function are generated, and in the next step, using the k-means clustering method, the number of scenarios is reduced to 10.

2) WIND SPEED UNCERTAINTY
The power output of WT is a function of wind speed, and the natural wind speed has high uncertainty. Using the Monte Carlo method, 1000 scenarios are generated by the Weibull distribution function implemented to simulate the uncertainty of wind speed as follows [36]: In the next step, by the k-means clustering method in GAMS, the number of scenarios is reduced to 10 scenarios to the best 10 scenarios within the stochastic circumstances. Next, the 10 scenarios are implemented as the wind speed data in the EH structure to achieve the most accurate scheduling.

3) PV UNCERTAINTY
Like a wind turbine, the power output of the PV is a function of solar irradiance and due to climate changes, the irradiance of PV will introduce high uncertainty. First, 1000 scenarios are generated with the Monte Carlo method by the Beta distribution function to simulate the uncertainty of PV irradiance as follows [36]: where ω and ψ are calculated as follow: Then, by the k-means clustering method in GAMS, the number of scenarios is reduced to 10 scenarios to obtain better scenarios within the stochastic circumstances. Subsequently, these scenarios are implemented as the solar irradiance data in the EH structure to achieve the most accurate scheduling.

F. CONDITIONAL VALUE AT RISK MODEL OF ENERGY HUB
Based on the previous discussion, renewable energy sources and load have uncertainty. Thus, the power supply flexibility of the energy hub can be increased using the ESS, fuel cell, and DRP. Moreover, the negative effect of uncertainty on energy hub performance is minimized by analyzing the uncertainty factor. One of the best options for an optimal decision scheme is risk level measurement [38]. The CVaR method is a type of risk-based method, which means that the loss surpasses the conditional mean value of VaR, reflecting the average potential loss by exceeding the VaR value. The graphic diagram of CVaR is plotted in Fig. 5 [39]. Here, the Conditional Value at Risk (CVaR) method is used, which the portfolio vector of uncertainty factor and random factor is calculated as follow, respectively: [38]: Also, the joint probability density function of y R n is p(y) and g is portfolio sets, G ⊂ R n , f (g, y) is the loss function and equal to OF OP . the probability of f (G, Y ) does not surpass the threshold α. The distribution function is formulated as follow Eq. (67): where, y ∈ , which is a subset of n-dimension real-number space R n and represents the feasible set of the portfolio. If α γ is the loss of (G, Y) so, the VaR and CVaR analysis of EH operation is given by: φ γ (y) is the CVaR value that the loss is bigger than α γ (y). Since Eq. (69) is difficult to solve, it is simplified to Eq. (70) using the estimation algorithm. In this formula, F(G, α) is replaced with φ(G), and the value of CVaR is formulated as follow: where, (f (G, y) − α) + is chosen the maximum value of {f (G, y) − α, 0}. Also, Eq. (70) is simplified as Eq. (71):

III. MULTI-OBJECTIVE SOLUTION METHOD
In the multi-objective optimization problem, considering all consistent and inconsistent objectives subjected to equality and inequality limitations, the optimization process is performed simultaneously. The three main reasons for implementing optimal methods are: • Identifying the non-dominated points of the solution in the objective functions and the dominated points of the solution in the decision-making space; • Creating a Pareto front in the objective functions' space; • Making variation points of the Pareto front to determine the Pareto optimal solutions in the decision-making space and provide more mandate to the determiner. Various types of intelligent algorithms have been introduced in this field, but one of their shortcomings is their low convergence [2]. Different methods have also been introduced in the literature for multi-objective optimization, such as the weight sum method, ε-constraint method, and augmented ε-constraint method.
In this study, a hybrid method is adopted to optimize system scheduling, consisting of two steps. In the first step, the multi-objective problem is converted into a single-objective problem using the augmented ε-constraint [40] (The main objective is chosen between the objectives, and the other objectives are defined as the constraints [40]), and then, the best response is selected from the Pareto front using the fuzzy mechanism [41].
where, x and δ are the minor (75)value in the interval [10 −6 , 10 −3 ], respectively; also, q 2 , q 3 , and q p are the equal intervals in f 2 , f 3 , and f p . s 2 , s 3 , . . . , s p are the surplus variables. Also, r 1 and r p are the minimum and maximum ranges of p objective function (pay-off table), respectively.

B. FUZZY METHOD
Finally, the fuzzy method is used to find the best solution from the Pareto curve and is considered as a variable in the interval of [0, 1]. Afterward, the max-min fuzzy technique converts both conflicting objective functions into their normalized forms [41]: After the normalization process, a comprehensive comparison between the solutions of the objective functions is performed, and then their minimum value is chosen.
The best cooperation solution to supply a trade-off between two contradictory objective functions is the maximum value of the selected minima, as follows: Finally, the decision-maker tries to maximize minimum satisfaction between the total objective functions. The flowchart of the proposed hybrid method is depicted in Fig. 6.

IV. NUMERICAL RESULTS AND SIMULATION
The proposed EH framework, as given in Fig. 1, includes three types of storage (TSS, ESS, and E-fuel) as well as renewable and fuel cell resources for supplying electricity, heating, gas, and water demands. It is investigated from economical and emission aspects under system uncertainties while considering the risk level using the CVaR method.

A. DATA
The data required for the proposed energy hub management are divided into two parts (deterministic data and probabilistic data) for accurate scheduling. Initially, the probabilistic data in 500 scenarios are created using the Monte Carlo method, and then, using K-means clustering, these scenarios are reduced to 10. Fig. 7 displays the energy, gas, heat, and water demand of the energy hub (kW). In Fig. 8, the electricity price of the upstream grid is plotted. Table 1 lists the required parameters of ESS, TSS, and E-fuel storage. The information related to a different unit in EH (boiler, etc.) is given in Table 2. In Table 3 and Fig. 9, the parameters of the CHP unit are listed. The emission factor of NO 2 , SO 2 , and CO 2 in each unit of EH in kg/kWh is presented in Table 4. The operation cost of numerous grids and units are given in Table 5. Moreover, the parameters of the wind turbine, PV,  and fuel cell are presented in Table 6. Finally, Fig. 10 depicts the 10 scenarios reduced by the k-means clustering method. The proposed energy hub problem with stochastic limitations   is nonlinear. Thus, it is solved using the equivalent MINLP problem and is codified with CPLEX in GAMS [42].

B. RESULTS AND EVALUATION
To evaluate the proposed method, the following 4 cases are defined to determine the effect of the novel E-fuel storage reserve market, while considering the risk level: Case 1. Without considering risk level and E-fuel storage; Case 2. Considering risk level without E-fuel storage; Case 3. Considering risk level and reserve market without E-fuel storage; Case 4. Considering risk level, E-fuel storage, and reserve market.
In all cases, wind turbine, PV, and fuel cell are considered as sources; load market price, wind speed, and PV irradiation have uncertainty; and TDRP, EDRP, GDRP, and WDRP are implemented. The Pareto curve of the proposed method for four cases is illustrated in Fig. 11, and the best optimal response is shown with a red point in all 4 cases selected by augmented ε-constraint and fuzzy method.       Table 7 and shown by a red point in Fig. 11). Furthermore, the output power of  the wind turbine, PV, and fuel cell based on the minimum operation cost and pollution emission are obtained and illustrated in Fig. 12. Most of the energy needed for load demand is supplied by the wind turbine. Finally, the percentage of electricity generation participation by all the units in the EH operation is depicted in Fig. 13.
Case 2: The variance between this case and Case 1 is applied by the CVaR method for considering the risk level. The Pareto curve for this case is displayed in Fig. 11. The operation cost and pollution emission are $1873.468 and    11683.880 kg, respectively (listed in Table 7 and shown by a red point). The imported power from a different unit such as FC, PV, WT, CHP, ESS, and upstream grid is plotted in Fig. 14. Evidently, in the peak time of load demand, the power generated by the units is increased so that employing a good energy storage is necessary for energy management. The operation cost emission factor for   this case is reduced by about 7.34% and 0.54% compared to Case 1 because of using the CVaR method. Moreover, in Table 8, the changes in operating cost for different values of α and β are presented. Based on Fig. 15, when β increases, the cost and pollution emission irregularly change, and by increasing α, these values regularly decrease. According to the results of Table 6, β = 0.99 and α = 0.15 are the best solutions.
Case 3: In this case, the novel energy storage which is the cold E-fuel is exploited in the energy hub structure, and its performance is analyzed. The electric load and power imported from a different unit are plotted in Fig. 16. The operation cost and pollution emission of this case are equal to $1849.885 and 11704.599 kg, respectively (listed in Table 7 and shown by a red point in Fig. 11). It is clear that, by implementing this storage, the operation cost and emission index are significantly decreased by about 1.26% and 0.68% compared to Case 2, respectively. This storage helps manage the energy demand by injecting energy into  EH in peak demand time and discharging in low demand time.
Case 4: In this case, all the proposed cases are added to the energy hub framework, such as CVaR, novel storage, and the reserve market. In the day-ahead market, the price of energy and reserve markets is predicted, and the purchasing power by the upstream network and reserve market is shown in Fig. 17. Based on the proposed solution, the total operation cost and the emission index are $ 1752.676 and 10829.530 kg, respectively. The amount of power generated by other units is reduced compared to the other cases because a part of load demand is purchased in the reserve market. Total operation cost and emission index are reduced by about 5.25%, and 6.67% compared to Case 3, respectively. The energy purchased from the upstream network and reserve market, and the price of electricity imported from the upstream network and reserve market, are shown in Fig. 17. Gas imported from the gas network, and gas generated for CHP and boiler operation for all 4 cases are shown in Fig. 18.
Finally, the amount of charge and discharge power for Case 4 are plotted in Fig. 19. By analyzing Fig. 19, it can be stated that most energies are produced by the novel energy storage because of its low cost and structure. Furthermore, based on Table 7, Case 4 has better results (lower cost and pollution) than other cases because of using a novel storage and implementing the reserve market.
Due to cooperation with the reserve market, when the cost of energy purchase is low, the required amount of energy  is purchased from the reserve market and used during the peak hours of energy consumption, which leads to less energy production by units and a significant reduction in pollution and costs.
In order to show the robustness of the proposed method, the results of the case 4 are compared with the others methods [16], [27], [28] which are listed in Table 9. It can be seen that, the proposed method is better than the other similar methods due to using a hybrid strategy in the choice of the best result in optimization procedure.

V. CONCLUSION
In this paper, a bi-objective optimization model was implemented for the robust scheduling of an optimal operation energy hub in the day-ahead and reserve markets. Operation cost and pollution emission were considered as a bi-objective model considering TDRP, WDRP, GDRP, and EDRP using the TOU demand response program to achieve the best result in an uncertain environment. The uncertain parameters were assumed to be wind speed, solar irradiance, and loads. The augmented ε-constraint and fuzzy methods were decomposed to simplify the bi-objective optimization problem. One of the main contributions of this paper is implementing the CVaR method to minimize the EH risk of high costs in worst scenarios as an appropriate risk estimation. Furthermore, considering different uncertain parameters increases the operation cost and pollution emission, so implementing appropriate storage is a good option.
Three types of storage, i.e., TSS, ESS, and E-fuel storage, were applied in the energy hub structure. Wind turbine, PV, and fuel cell, in addition to the grid, gas, and water network were considered as renewable resources. To assess the efficiency of the proposed method, four case studies were evaluated. In Case 1, the base structure of the energy hub was considered. In Case 2, using the CVaR method, in addition to different DRP, the operation cost, and pollution emission index were decreased by 7.34% and 0.54% compared to Case 1, respectively. In the next step, using the novel E-fuel storage and fuel cell as a resource, the operation cost and pollution emission index were minimized by 1.26% and 0.68% compared to Case 2, respectively. Finally, the reserve market and energy market cooperated with the energy hub for bidding and sailing energy. In this mode, the operation cost and pollution emission index were impressively reduced by about 5.25%, and 6.67%, respectively, compared to Case 3. It can be concluded that the operation cost and pollution emission scenario are highly minimized by implementing the CVaR method in an uncertain environment and by using a novel E-fuel storage and reserve market participation. He has published more than 140 articles in international journals and conference proceedings. His research interests include the application of artificial intelligence and heuristic optimization algorithms to power system design, FACTS device, power system analysis and control, renewable energy, control and management of microgrids, and smart grid. He was selected as a Distinguished Researcher of the Azarbaijan Shahid Madani University, in 2016.

HASSAN HAES ALHELOU (Senior Member, IEEE) is currently a Faculty Member with
Tishreen University, Lattakia, Syria. He is also with University College Dublin, Ireland. He is included in the 2018 and 2019 Publons list of the top 1% best reviewer and researchers in the field of engineering. He has published more than 140 research articles in the high-quality peer-reviewed journals and international conferences. He has participated in more than 15 industrial projects. His current research interests include power systems, power system dynamics, power system operation and control, dynamic state estimation, frequency control, smart grids, micro grids, demand response, load shedding, and power system protection. He