Optimal Capacity Planning of Power to Hydrogen in Integrated Electricity–Hydrogen–Gas Energy Systems Considering Flexibility and Hydrogen Injection

With increasing penetration of renewable energy, it is important to source adequate system flexibility to maintain security of supply and minimize renewable generation curtailment. Power to hydrogen (P2H) plays an important role in the low-carbon renewable dominated energy systems. By blending green hydrogen produced from renewable power into the natural gas pipelines, it is possible to help integrate large-scale intermittent generation and smooth the variability of renewable power output through the interconnection of the natural gas network, hydrogen energy network, and electric network. A two-stage stochastic mixed-integer nonlinear planning framework for P2H sizing and siting is proposed in this paper, considering system flexibility requirements. The problem is then reduced to a mixed-integer second-order cone (MISOC) model through convex transformation techniques in order to reduce the computation burden. Then, a distributed algorithm based on Bender’s decomposition is applied to obtain the optimal solution. A modified hybrid IEEE 33-node and Gas 20-node system is then used for simulation tests. The results showed that investment of P2H can significantly reduce the total capital and operational costs with lower renewable generation curtailment and electricity demand shedding. Numerical tests demonstrated to demonstrate the validity of the proposed MISOC model.

1 Introduction

Over the recent years, various renewable energies, such as solar power and wind power, have seen rapid development. However, the connection of large-scale renewable energy to the electric network has dramatically changed the characteristics and increased the uncertainty of power flow, posing significant challenges to system operations in terms of power balancing and load following. Flexible resources, such as demand-side response, large-scale energy storage, power to hydrogen, etc., are the key to integrating large-scale renewable generation.

For environmental reasons, coal-fired generations are gradually being phased out, leading to a scarcity of flexible resources. As a result, energy storages, such as pumped storage, compressed air energy storage, chemical battery storage, etc., have developed by leaps and bounds in recent years. However, the development of large-scale pumped storage and compressed air storage is constrained by geographic and geological conditions. Chemical battery storage is not limited by geography, but their discharge time is 1 min–8 h (Mongird et al., 2019) and is suitable for medium- or short-term (Zablocki, 2019) use. As hydrogen has the appealing characteristics of being stored for the longer term and on a large scale, hydrogen energy storage has the prospects of being used as a long-term or seasonal energy storage. With the high penetration of renewable energy, the green hydrogen produced by P2H using surplus renewable energy can be stored on a large scale or used as an industrial raw material.

As a result, over the past few years, the development of green hydrogen has attracted the increasing attention of researchers and governments due to its characteristics of a high energy density and the ability to store it in large quantities. The total abandoned renewable energy has decreased from 312.523 to 66.797 MWH with the presence of power to hydrogen (Heris et al., 2020). To reduce CO₂ emissions, the cost-optimal design and operation of power to gas in Germany are studied in depth (Welder et al., 2018). In addition, the value of hydrogen storage in pipes and salt caverns has been analyzed (Welder, 2018).

In addition, Ref. Dolci et al. (2019) and Hu et al. (2020) present the power-to-hydrogen pathways, and the technical advances and barriers to green hydro development are also given. The green hydrogen-based power-to-gas facility is investigated in Germany and Texas in Glenk and Reichelstein (2019), which shows that green hydrogen production is already cost-competitive in niche applications. With ultra-high penetration of renewable energy, an optimal planning model (Pan, 2020) for electricity–hydrogen integrated energy systems considering power to hydrogen and seasonal storage is constructed in Pan et al., (2020). Li et al. (2019) have developed operation optimization modeling of power to hydrogen in active distribution networks coordinated with the district heating network (Li et al., 2019). Gils et al. (2019) have investigated the transition pathway toward a climate-neural energy system in Germany and concluded that hydrogen plays a critical role in the seasonal balance of power supply and demand but cannot completely replace other coupling options (Gils et al., 2021).

Hydrogen is an important industrial raw material, and it is only economically beneficial to convert hydrogen into methane when subsidies are provided (Robinius et al., 2017). Liu et al. (2020) indicate that both power to hydrogen (P2H) and power to methane (P2M) help to reduce operational costs and decrease CO₂ emissions. From the perspective of energy system operation, the impact of P2H on the system flexibility is analyzed (Cloete and Hirth, 1922) in Cloete et al. (2020) and Ge et al. (2020) . With the increasing penetration level of renewable energy, the traditional grid is evolving into a new energy system with widespread participation of variable generation and flexible demand (Gea-Bermúdez et al., 2021). The large-scale, long-duration energy storage devices will be one of the key enablers in the low-carbon renewable-dominated energy systems. To improve the management efficiency and deal with the stochastic nature of renewable energy, this paper aims to investigate the site selection and optimal capacity determination approach of large-scale energy storage technologies, also known as P2H, and their role in enhancing system flexibility.

The contributions and innovations of this article are as follows: first, a two-stage stochastic mixed-integer planning framework for P2H sizing and siting is proposed, considering the uncertainty of the renewable energy output and the investment cost of the P2H equipment. The model’s objective was to minimize the investment cost, operation cost, renewable output curtailment cost, and electricity demand shedding cost. Second, by blending permissible green hydrogen produced from surplus renewable energy into the natural gas pipelines, the effects of the injected hydrogen concentration on the gas pipeline transmission coefficient, that is, gas specific gravity, net calorific value, gas compressibility factor, etc., were investigated. Third, based on the system flexibility supply–demand balance mechanism, we quantified the enhanced effect of P2H on hourly system flexibility.

The rest of this article is organized as follows: the low-carbon-emission hybrid renewable energy system structure is proposed in Section 2, and the modeling of renewable energy and load demand uncertainties are also described. In Section 3, a two-stage stochastic planning model of P2H has been constructed, considering system flexibility and hydrogen injections. The nonlinear nonconvex terms are simplified to obtain a mixed-integer second-order cone (MISOC) programming model in Section 4. In Section 5, a distributed algorithm for the model is designed based on Bender’s decomposition and cut plane. Simulations are carried out to verify the effectiveness of the proposed MISOC model and solution algorithm in Section 6. The conclusion is illustrated in Section 7.

2 Integrated Electricity–Hydrogen-Gas System

2.1 System Description

The integrated electricity–hydrogen-gas energy system is illustrated in Figure 1. It consists of the power system, natural gas distribution system, and carbon-free hydrogen system. The interaction between the power system and the natural gas system is realized mainly through a combined cycle gas turbine (CCGT). The CCGT generation unit with a fast-responding, high efficiency uses natural gas from the distribution gas system to generate electricity. The connection of large-scale renewable energy to the distribution network makes it possible to produce hydrogen by electrolyzing water with renewable energy. It is worth noting that hydrogen is more accessible to store for long periods than electricity. Therefore, hydrogen storage technology can be used as a medium- or long-term energy storage option for excess renewable energy. In this article, green hydrogen produced by surplus renewable energy sources plays a vital role in the integrated electricity distribution system, natural gas system, and hydrogen energy system. Power-to-hydrogen devices make full use of surplus renewable energy, and it is then blended and injected into the existing gas distribution pipeline with a maximal volumetric limit proportion of 15% hydrogen.

FIGURE 1

FIGURE 1. Basic structure diagram of a hybrid renewable energy system.

In future high-renewable energy penetration scenarios, the uncertainty and intermittence of renewable energy output will significantly impact supply power reliability. To deal with these problems, the system’s flexibility needs to be enhanced urgently. Based on the characteristics of flexible resources in the distribution network, we have developed a system flexibility supply and demand balance model. In addition, considering the direct injection of hydrogen produced by P2H into the natural gas pipelines, the hydrogen concentration on the gas specific gravity, the calorific value, and the compressibility factor have been investigated. Finally, a two-stage stochastic planning model is constructed with the objective of minimizing the investment and operation costs during the planning horizon. In the proposed hybrid energy system, we consider several renewable energy generation units, such as wind turbines and photovoltaics, equipped in the system. To promote renewable energy consumption and enhance the system’s flexibility, the investment cost of P2H and the renewable energy curtailment cost are included as one of the optimization objectives. Scenario generation and reduction methods have been utilized to address renewable energy sources and load uncertainty.

2.2 Uncertainty Modeling and Scenario Generation

A scenario-based approach is adopted to deal with the uncertainty and stochastic characteristics of renewable output and electrical load demand (Dagoumas and Koltsaklis, 2019). For the sake of simplicity, renewable generation output and load can be modeled as follows:

$\begin{array}{c}{P}_{load,t}=P_{load,t}^{\text{fore}}+\text{Δ}{P}_{load,t}t=1,\dots ,24\\ P_{wind,t}=P_{wind,t}^{\text{fore}}+\text{Δ}{P}_{wind,t}t=1,\dots ,24\end{array}(1)$

The forecast errors, $\Delta P_{wind,t},and\Delta P_{load,t}$, are assumed to be Gaussian distribution. The probability density functions of normally distributed random variables are discretized in Figure 2​

$\begin{array}{l}\mathrm{\Delta }{P}_{i,t}^{load}={\alpha }_{i,t}^{load}\text{i}\in \text{interval}1\text{,}\cdots \text{interval}7\\ \mathrm{\Delta }{P}_{j,t}^{wind}={\beta }_{j,t}^{wind}\text{j}\in \text{interval}1\text{,}\cdots \text{interval}7\end{array}(2)$

where $\Delta P_i^{load},and\Delta P_j^{wind}$ are the probabilities of intervals 1–7 at time $t$ related to prediction errors, and their specific values are ${\alpha }_{i,t}^{load},and{\beta }_{j,t}^{wind}$.

FIGURE 2

FIGURE 2. Wind out power/load forecast error.

Then, the corresponding cumulative distribution functions (CDFs) are calculated separately using the probability density functions for the discrete intervals described above. Based on Monte Carlo simulations, two random variables, $\eta \in \left[\mathrm{0,1}\right]\text{ and }\lambda \in \left[\mathrm{0,1}\right]$, are generated. The randomly generated variables $\eta \text{ and }\lambda$ are compared with the cumulative probability distribution functions values above, respectively, to determine in which interval $\eta \text{ and }\lambda$ lie in. If $CDF\left(\mathrm{int}erval1\right)\le \eta <CDF\left(\mathrm{int}erval2\right)$, then, the binary variable $B_{2,t,s}^{Load}$ equals to 1, and the other binary variables $B_{(1,t,s)}^{Load},\dots ,B_{(7,t,s)}^{Load}$ related to load prediction errors equal to 0; the same applies to $\lambda$

$\begin{array}{l}\text{ scenario }\left(\text{s}\right)={\left[B_{(1,t,s)}^{Load},\dots ,B_{(7,t,s)}^{Load},B_{(1,t,s)}^{Wind},\dots ,B_{(7,t,s)}^{Wind}\right]}_{t=1,\dots ,24}\\ {\displaystyle \sum _{\mathrm{int}erval=1}^7{B}_{(\mathrm{int}erval,t,s)}^{Load}=1}\forall \text{t,}\forall \text{s}\\ {\displaystyle \sum _{\mathrm{int}erval=1}^7{B}_{(\mathrm{int}erval,t,s)}^{Wind}=1}\forall \text{t,}\forall \text{s}\end{array}(3)$

The probabilities of the generated scenarios are listed as follows, and the sum of the probabilities of all scenarios equals to 1

${\pi }_s=\frac{{\displaystyle \prod _{24}^{t=1}\left({\displaystyle \sum _7^{i=1}\left(B_{i,t,s}^{Load}\times {\alpha }_{i,t}^{load}\right)\times {\displaystyle \sum _7^{j=1}{B}_{j,t,s}^{Wind}\times {\beta }_{j,t}^{wind}}}\right)}}{\begin{array}{l}{\displaystyle \sum _{N_s}^{s=1}\left({\displaystyle \prod _7^{t=1}\left(B_{i,t,s}^{Load}\times {\alpha }_{i,t}^{load}\right)\times {\displaystyle \sum _7^{j=1}\left(B_{j,t,s}^{Wind}\times {\beta }_{j,t}^{wind}\right)}}\right)}\\ {\displaystyle \sum _{s=1}^{N_s}{\pi }_s}=1\end{array}}(4)$

Due to a large number of scenarios, it will lead to the burden of computation. It is generally essential to balance the solving tractability and the modeling accuracy. Therefore, the scenario reduction techniques using the fast forward algorithm (Hajiabbas and Mohammadi-Ivatloo, 2020) are listed as follows:

$\begin{array}{l}{P}_s={\displaystyle \sum _{t=1}^{24}{P}_{t,s}^{load}}+{\displaystyle \sum _{t=1}^{24}{P}_{t,s}^{wind}}\forall \text{s}\\ v\left(s,s^{,}\right)=\left|P_s-P_{s^{,}}\right|\forall \text{s}\end{array}(5)$

where $P_{t,s}^{load}\text{ and }{P}_{t,s}^{wind}$ are the load demand and renewable energy output at time $t$ under the typical scenario $s$, respectively. $v\left(s,s^{,}\right)$ is a matrix of $N_s$ multiplied by $N_s$ , which is the distance between $P_s$

$\begin{array}{c}{s}_1=\arg \left\{\begin{array}{c}\min \\ s^{\prime }\in \Omega \end{array}{\displaystyle \sum _{s\in \Omega }{\pi }_{s}v\left(s,s^{"}\right)}\right\}\\ {\text{Ω}}_S={\text{Ω}}_S\cup \left\{s_1\right\}\\ {\text{Ω}}_J=\text{Ω}\mathrm{\setminus }\left\{s_1\right\}\end{array}(6)$

After obtaining the first selected scenarios through the above steps, the distance between the chosen and nonselected scenarios is calculated, and the set of scenes chosen is updated

$\begin{array}{c}{s}_i=\arg \left\{\begin{array}{c}\min \\ S^{\prime }\in \Omega \end{array}\begin{array}{c}\mathrm{\Sigma }\\ s\in {\text{Ω}}_J\mathrm{\setminus }\left\{s^{"}\right\}\end{array}{\pi }_s\begin{array}{c}\min \\ s{s}^{,,}\in {\text{Ω}}_S\cup \left\{s\right\}\end{array}v\left(s,s^{,,}\right)\right\}\\ {\text{Ω}}_S={\text{Ω}}_S\cup \left\{s_i\right\}\\ {\text{Ω}}_J={\text{Ω}}_J\mathrm{\setminus }\left\{s_i\right\}\end{array}(7)$

Finally, the probabilities of the reduced scenarios are calculated. The probability of a nonselected scenario is added to the probability of a selected scenario that is close to it. The specific calculation is as follows:

${\pi }_s^{\ast }={\pi }_s+\begin{array}{c}\mathrm{\Sigma }\\ s^{\mathrm{,}}\in J\left(s\right)\end{array}{\pi }_{s^{,}}S\in {\text{Ω}}_s{S}^{\prime }\in {\text{Ω}}_J(8)$

3 Optimal P2H Planning Considering Flexibility Requirements in an Integrated Multi-Vector Energy System

In this section, a two-stage stochastic dynamic mixed-integer nonlinear programming model is constructed for optimal P2H planning to minimize the investment and annual operation costs (Bramstoft et al., 2020) of integrated multi-vector energy systems.

3.1 The Objective Function

The optimal P2H planning problem considering the system flexibility requirement can be formulated as the following stochastic model:

$\min f=C_{inv}+C_{ope}(9)$

$C_{inv}={\displaystyle \sum _{i\in {\Omega }_{p2h}}\frac{r{(1+r)}^{{\tau }_{i,\xi }}}{{(1+r)}^{{\tau }_{i,\xi }}-1}{I}_i{P}_i^{\max }{c}_i}(10)$

$\begin{array}{l}{C}_{ope}=D{\displaystyle \sum _{s=1}^{N_s}\pi \left(S\right)}{\displaystyle \sum _{t=1}^T[{\displaystyle \sum _{i\in {\Omega }_{ele}}{e}_i^{ele}{P}_{i,t,s}^{ele}}+}\\ {\displaystyle \sum _{n\in {\Omega }_{saurce}}{g}_n^{saurce}{P}_{n,t,s}^{saurce}}+{\displaystyle \sum _{j\in {\Omega }_{renew}}{c}_j^{renew}\text{Δ}{P}_{j,t,s}^{renew}}\\ +{\displaystyle \sum _{k\in {\Omega }_{load}}{c}_k^{load}\text{Δ}{P}_{k,t,s}^{load}}+{\displaystyle \sum _{m\in {\Omega }_{gas}}{c}_m^{gas}\text{Δ}{P}_{m,t,s}^{gas}}]\end{array}(11)$

subject to

$0\le P_i^{\max }\le {\overline{P}}_i^{\max }\forall \text{i}(12)$

$I_i\in \left\{\mathrm{0,1}\right\}\forall \text{i}\in {\text{Ω}}_{p2h}(13)$

$0\le {\displaystyle \sum _{i\in {\Omega }_{p2h}}{I}_i}\le N^{\max }(14)$

${\displaystyle \sum _{i\in {\Omega }_{p2h}}{I}_i{P}_i^{\max }}\le P^{\max }(15)$

$P_{i,t,s}^{ele},P_{n,t,s}^{saurce},\text{Δ}{P}_{j,t,s}^{renew},\text{Δ}{P}_{k,t,s}^{load},\text{Δ}{P}_{m,t,s}^{gas}\in {\text{Ω}}_o\forall \text{t,}\forall \text{s }(16)$

The objective function (9) comprises two terms: P2H investment cost $C_{inv}$ and operation cost $C_{ope}$. $C_{ope}$ includes the costs of purchasing electricity and gas, the renewable energy curtailment cost, and the electricity and gas load shedding cost. It is worth noting that the generation cost of CCGT is not included in the operation cost as it is already calculated in the output costs of the gas sources. Equation 12 sets the maximum capacity investment of each P2H. Equations 14, 15 limit the number of installed P2H devices and the total investment cost. The binary variable $I_i$ is a decision variable on whether a prospective P2H device is constructed at node $i$ or not. $D,e_i^{ele},c_j^{renew},c_k^{load}$ in Eq. 12 are the number of days in the year, the time-of-use electricity price, the cost of abandoned renewable energy, and the load (Liu et al., 2020), respectively.

3.2 Operation Constraints

3.2.1 Power Distribution Network

The power flow equations of the radial distribution network can be presented by the DistFlow model, which can be seen as follows (Li et al., 2020):

$P_{ij}+P_j^g-r_{ij}{i}_{ij}={\displaystyle \sum _{k\in \pi \left(j\right)}{P}_{jk}+P_j^d\forall t,\forall s}(17)$

$Q_{ij}+Q_j^g-x_{ij}{i}_{ij}={\displaystyle \sum _{k\in \pi \left(j\right)}{Q}_{jk}+Q_j^d\forall t,\forall s}(18)$

$U_j=U_i-2\left(r_{ij}{P}_{ij}+x_{ij}{Q}_{ij}\right)+{\left(z_{ij}\right)}^2{i}_{ij}\forall t,\forall s(19)$

$i_{ij}{U}_i=P_{ij}^2+Q_{ij}^2\forall t,\forall s(20)$

Eq. 17 above can be modified with high penetration renewable energy connection as follows:

$P_{ij}+P_j^g-r_{ij}{i}_{ij}+P_j^{renew}-\text{Δ}{P}_j^{renew}+P_j^{CCGT}={\displaystyle \sum _{k\in \pi \left(j\right)}{P}_{jk}+P_j^d-\mathrm{\Delta }{P}_j^{load}+P_j^{p2h}}(21)$

where $P_{ij},Q_{ij}$ are the active and reactive powers through line $l\left(ij\right)$, respectively. Definitions of the several optimization variables $U_i,r_{ij},i_{ij},x_{ij},z_{ij},P_j^d,Q_j^g,Q_j^d$ are in Ref. Heris et al., 2020.

We can depict the boundary constraints in the power distribution network as follows:

$0\le P_{j,t,s}^{CCGT}\le P_j^{CCGT,MAX}\forall j,\forall t,\forall s(22)$

$0\le P_{j,t,s}^{renew}\le P_{j,t,s}^{renew,\max }\forall j,\forall t,\forall s(23)$

$0\le P_{j,t,s}^{ele}{\le }_{j,t,s}^{ele,\max }\forall \text{j,}\forall \text{t,}\forall \text{s}(24)$

$0\le \text{Δ}{P}_{j,t,s}^{load}\le P_j^d\forall \text{j,}\forall \text{t,}\forall \text{s}(25)$

$0\le P_{j,t,s}^{p2h}\le P_j^{\max }\forall j,\forall t,\forall s(26)$

$P_{ij}^{\min }\le P_{ij,t,s}\le {\overline{P}}_{ij}^{\max }\forall t,\forall s(27)$

$Q_{ij}^{\min }\le Q_{ij,t,s}\le {\overline{Q}}_{ij}^{\max }\forall t,\forall s(28)$

$U_j^{\min }\le U_{j,t,s}\le U_j^{\max }{\text{ U}}_{ref}=V_{ref}^2(29)$

Equation 21 constrains the power balance at node $j$, considering the curtailment of renewable energy and shedding load. Equations 22, 23 enforce limits on the output of the CCGT generation unit and each renewable energy source. Equation 24 depicts the boundary of the amount of electricity purchased. Equations 25–29 give the lower and upper bounds on the variables $\text{Δ}{P}_{j,t,s}^{load},P_{j,t,s}^{p2h},P_{ij,t,s},Q_{ij,t,s},U_{j,t,s}$.

3.2.2 Natural Gas Distribution Network

The steady-state flow equation Eq. 20 related to the pressure of gas and the flow rate in the distributed gas network is depicted as follows:

$Q_{ij}=1.1494\times {10}^{-3}\left(\frac{T_b}{P_b}\right){\left[\frac{\left(\left|P_i^2-P_j^2\right|\right)}{G{T}_{f}LZf}\right]}^{0.5}{D}^{2.5}(30)$

where $G,T_f,L,Z,f$ indicate the gas specific gravity, the average flowing temperature of the gas, the pipe length (km), the gas compressibility factor, and the friction factor, respectively. Terms $T_b,P_b$ are the base temperature and pressure, respectively.

It is important to note that the gas thermal–physical properties change with its composition (Deymi-Dashtebayaz et al., 2019). When green hydrogen produced from renewable energy sources is injected into a natural gas pipeline, the gas thermo-physical properties vary with the hydrogen concentration, which can be described below:

$Z_i=a\cdot v_i+b\text{ i}\in {\Omega }_{saurce}(31)$

$\text{L}\text{H}{\text{V}}_i^{gas}=c\cdot v_i+d\text{i}\in {\Omega }_{saurce}(32)$

$G_i=e\cdot v_i+f\text{ i}\in {\Omega }_{saurce}(33)$

$0\le v_i\le 15\%\text{ i}\in {\Omega }_{saurce}(34)$

where $LH{V}_i^{gas}$ is the net calorific value of gas mixtures at node $i$ and $v_i$ is the hydrogen concentration by volume injected into the pipeline. In addition, $a,b,c,d,e,f$ are constants. From Eqs 31–33 above, we can also see that $Z_i,LH{V}_i^{gas},G_i$ are each approximately linearly related to the hydrogen concentration $v_i$ at gas node i (Dolci, 2019). Equation 34 enforces the upper and lower bounds on the concentration of hydrogen in a natural gas pipeline. No additional investment in gas distribution pipelines and end-use appliances is needed within this concentration range.

In addition to the above constraints, certain boundary constraints need to be considered

$P_{i,t,s}^{\min }\le P_{i,t,s}\le P_{i,t,s}^{\max }\forall \text{i,}\forall \text{t,}\forall \text{s}(35)$

$Q_{l,t,s}^{\min }\le Q_{l,t,s}\le Q_{l,t,s}^{\max }\forall l\left(ij\right)\text{,}\forall \text{t,}\forall \text{s}(36)$

$0\le P_{n,t,s}^{saurce}\le P^{saurce,\max }\text{ n}\in {\Omega }_{saurce},\forall t,\forall s(37)$

$0\le \text{Δ}{P}_{m,t,s}^{gas}\le P_{m,t,s}^{gas,load}\text{ m}\in {\Omega }_{gas},\forall t,\forall s(38)$

$P_{j,t,s}={\beta }_{ij}\cdot P_{i,t,s}\forall \text{t,}\forall \text{s}(39)$

Equations 35, 36 enforce limits on the pressure of the gas node and the flow rate in each pipeline. Equations 37, 38 indicate the upper and lower bounds for natural gas source output and shedding gas load.

Pressure drops occur during the transmission of natural gas. Therefore, a compressor is placed in the natural gas pipeline to maintain the nodal pressure at the required level. Equation 39 indicates the gas pressure relationship at nodes i and j. ${\beta }_{ij}$ is the ratio of the compressor inlet pressure to the outlet pressure.

Considering the effect of hydrogen injection, the energy balance equation for each natural gas distribution network node is formulated as follows:

$\begin{array}{l}{P}_{j,t,s}^{saurce}+(1/\alpha )P_{j,t,s}^{h2\_saurce}+{\displaystyle \sum _{i\in \Omega s}{Q}_{ij,t,s}}={\displaystyle \sum _{k\in \Omega r}{Q}_{jk,t,s}}+\\ {\displaystyle \sum _{GT\in {\Omega }_{cg}}{P}_{GT,t,s}^{gas}}+P_{j,t,s}^{gas,load}-\text{ Δ}{P}_{j,t,s}^{gas},\forall t,\forall s\end{array}(40)$

$\alpha =\frac{LH{V}^{C{H}_4}}{LH{V}^{H_2}}(41)$

$P_{i,t,s}^{h2\_saurce}=C_{ij}\cdot P_{j,t,s}^{h2}(42)$

where $LH{V}^{C{H}_4},LH{V}^{H_2}$ are the net calorific values of natural gas and hydrogen, respectively, while $\alpha$ is their ratio. Equation 40 enforces the nodal gas balance for the gas distribution network. In this equation, the shedding gas load and hydrogen injection are considered. Due to the decrease in the Wobbe index of the gas mixture, it is worth considering the differential properties of the calorific value of natural gas and hydrogen. Equation 42 indicates that the green hydrogen produced by P2H is coupled to the gas source and then blended with natural gas at a specific concentration.

3.2.3 Flexibility Supply and Demand Balance

When renewable energy’s share is relatively low, the volatility and uncertainty have less impact on the energy system. At this stage of the planning process, the flexibility supply and demand balance in the system is usually not considered for the time being. However, with a high proportion of renewable energy, insufficient flexibility supply can lead to a large number of renewable energy curtailments, which prevents the system from operating correctly (Koltsaklis and Dagoumas, 2018). Therefore, the system’s flexibility is quantified in the process of optimizing the size of P2H.

Flexibility in the integrated energy system is characterized by multiple time scales, locations, and directions (Agbonaye et al., 2021). From the prospect of time scale, system flexibility can be categorized into very short-term ramp ($\le 15\min$), short-term ramp ($15\min -4h$), and peak shaving ($24h$) (Lu et al., 2018). This article focuses on the supply–demand of system flexibility at an hourly scale.

Flexibility resources such as P2H, CCGT, etc., are introduced to supply flexibility. Of these, P2H, CCGT, and system power purchase can provide both upward and downward flexibility, while wind abandonment can only provide downward flexibility, and shedding load can only supply upward flexibility. The detailed formulae are shown below:

$F_{p2h,t,s}^{s\_up}={\displaystyle {\sum }_{p2h\in {\Omega }_{p2h}}\min (R_{p2h}^{down}\tau ,P_{p2h,t,s}^{p2h}-P_{p2h,t,s}^{\min })}\forall t,\forall s(43)$

$F_{p2h,t,s}^{s\_down}={\displaystyle {\sum }_{p2h\in {\Omega }_{p2h}}\min (R_{p2h}^{up}\tau ,P_{p2h}^{\max }-P_{p2h,t,s}^{p2h})}\forall t,\forall s(44)$

$F_{GT,t,s}^{s\_up}={\displaystyle {\sum }_{GT\in \Omega gt}\min (R_{GT}^{up}\tau ,P_{GT}^{CCGT\_MAX}-P_{GT,t,s}^{CCGT})}\forall t,\forall s(45)$

$F_{GT,t,s}^{s\_down}={\displaystyle \sum _{GT\in \Omega gt}\min (R_{GT}^{down}\tau ,P_{GT,t,s}^{CCGT}-P_{GT}^{CCGT\_MIN})}\forall t,\forall s(46)$

$F_{buy,t,s}^{s\_up}={\displaystyle \sum _{buy\in {\Omega }_{ele}}\min (R_{buy}^{up}\tau ,P_{buy,t,s}^{ele,\max }-P_{buy,t,s}^{ele})}\forall \text{t,}\forall \text{s}(47)$

$F_{buy,t,s}^{s\_down}={\displaystyle \sum _{buy\in {\Omega }_{ele}}\min (R_{buy}^{down}\tau ,P_{buy,t,s}^{ele}-P_{buy,t,s}^{ele,\min })}\forall \text{t,}\forall \text{s}(48)$

$F_{k,t,s}^{s\_up}={\displaystyle \sum _{k\in {\Omega }_{load}}\mathrm{\Delta }{P}_{k,t,s}^{load}}\forall \text{t,}\forall \text{s}(49)$

$F_{j,t,s}^{s\_down}={\displaystyle \sum _{j\in {\Omega }_{renwe}}\text{Δ}{P}_{j,t,s}^{renew}}\forall \text{t,}\forall \text{s}(50)$

$F_{sum,t,s}^{s\_up}=F_{p2h,t,s}^{s\_up}+F_{GT,t,s}^{s\_up}+F_{buy,t,s}^{s\_up}+F_{k,t,s}^{s\_up}\forall t,\forall s(51)$

$F_{sum,t,s}^{s\_down}=F_{p2h,t,s}^{s\_down}+F_{GT,t,s}^{s\_down}+F_{buy,t,s}^{s\_down}+F_{j,t,s}^{s\_down}\forall \text{t,}\forall \text{s}(52)$

where $R_{p2h}^{up},R_{GT}^{up},R_{buy}^{up}$ represent the upward ramp rates of P2H, CCGT, and the power purchase, respectively, and $R_{p2h}^{down},R_{GT}^{down},R_{buy}^{down}$ are the downward ramp rates of P2H, CCGT, and the power purchase, respectively. $F_{sum,t,s}^{s\_up},F_{sum,t,s}^{s\_down}$ denote the total upward and downward flexibility capacities provided by the various flexibility resources, respectively.

Equations 43, 45, 47, 49 define the upward flexibility capacity, and Eqs 44, 46, 48, 50 explain the downward flexibility capacity. It is also observed that Eqs 51, 52 give the total flexibility supply capacity.

The flexibility demand arises mainly from fluctuations in the net load. Therefore, the proposed hybrid system’s upward and downward flexibility demand can be described by the following equations:

$P_{t,s}^{net}={\displaystyle \sum _{k\in {\Omega }_{load}}{P}_{k,t,s}^d-}{\displaystyle \sum _{j\in {\Omega }_{renew}}{P}_{j,t,s}^{renew}}\forall \text{t,}\forall \text{s}(53)$

$F_{t,s}^{D\_up}=F_{t+1,s}^{net}-F_{t,s}^{net}(F_{t+1,s}^{net}\ge F_{t,s}^{net})(54)$

$F_{t,s}^{D\_down}=F_{t,s}^{net}-F_{t+1,s}^{net}(F_{t+1,s}^{net}\le F_{t,s}^{net})(55)$

The system flexibility supply–demand equilibrium equation is as follows:

$F_{t,s}^{up}=F_{sum,t,s}^{s\_up}-F_{t,s}^{D\_up}\ge 0\forall \text{t,}\forall \text{s}(56)$

$F_{t,s}^{down}=F_{sum,t,s}^{s\_down}-F_{t,s}^{D\_down}\ge 0\forall \text{t,}\forall \text{s}(57)$

Equations 56, 57 emphasize that the system’s ability to supply flexibility is greater than the demand for flexibility.

3.2.4 Modeling of Energy Coupling Devices

A schematic of the coupling relation between electricity–gas-hydrogen energy networks is illustrated in Figure 3. P2H plays a crucial role in coupling multi-energy networks, and CCGT links the gas networks and power grid effectively.

FIGURE 3

FIGURE 3. Schematic of the coupling relation between electricity–gas–hydrogen energy networks.

The simplified P2H and CCGT models are given in the following equations:

$P_{j,t,s}^{CCGT}={\eta }_{GT}\cdot P_{GT,t,s}^{gas}\text{ GT}\in {\text{Ω}}_{cg},j\in {\text{Ω}}_{gt},\forall t,\forall s(58)$

$P_{j,t,s}^{h2}={\gamma }_j\cdot P_{j,t,s}^{p2h}j\in {\text{Ω}}_{p2h},\forall t,\forall s(59)$

where ${\eta }_{GT},{\gamma }_j$ are the energy conversion efficiencies of CCGT and P2H, respectively.

4 Model Simplications

The model constructed above is a mixed-integer nonlinear nonconvex programming model, which is difficult to solve quickly and efficiently. It is therefore essential to simplify the model to reduce the computational burden.

The nonlinear terms $I_i{P}_i^{\max }$ , $i_{ij}{U}_i$ and terms of form $\min (X,Y)$ are dealt with in Appendix A.

There are quadratic terms in Eq. 30, so Eq. 30 is nonlinear and nonconvex. The results tend to fall into the local optimum, and it is not easy to obtain the global optimum exactly. Here, Eq. 30 can be reformulated as follows:

$\rho \cdot Q_{ij}^2=\xi \cdot \left(P_i^2-P_j^2\right)(60)$

$\xi =\{\begin{array}{l}1{\text{ P}}_i\ge {\text{P}}_j\\ -1{\text{ P}}_i<{\text{P}}_j\end{array}(61)$

$\rho =\frac{{\left(G{T}_{f}LZf{P}_b\right)}^2}{T_{\mathit{b}}^2{D}^5{\left(1.1494\times {10}^{-3}\right)}^2}(62)$

Let ${\pi }_i=P_i^2$, ${\pi }_j=P_j^2$. Then, Eq. 62 can be rewritten as follows:

$\rho \cdot Q_{ij}^2=\xi \cdot \left({\pi }_i-{\pi }_j\right)(63)$

Let $x_l$ be the gas flow direction. When $x_l=1$, it means $Q_{ij}\ge 0$. When $x_l=0$, it denotes $Q_{ij}<0$. Therefore, the transmission characteristic constraint and capacity constraint denoted by Eq. 65 can be equivalently expressed as follows:

$(2x_l-1){\pi }_i+(1-2x_l){\pi }_j=\rho \cdot Q_{ij}^2(64)$

$-Q_{ij}^{\max }(1-x_l)\le Q_{ij}\le x_l(65)$

Equation 66 can be further relaxed to

$(2x_l-1){\pi }_i+(1-2x_l){\pi }_j\ge \rho \cdot Q_{ij}^2(66)$

For the line $l=(i,j)$, new variables $z_{li}=x_l{\pi }_i$ and $z_{lj}=x_l{\pi }_j$ are introduced due to the product of variables $x_l$ and $\pi$. Constraint Eq. 68 is further replaced with equations in Singh and Kekatos, 2020. The nonlinear, nonconvex terms are converted into mixed-integer second-order cone constraints with the above transformation

$x_l{\underset{¯}{\pi }}_i\le z_{li}\le x_l{\overline{\pi }}_i(67)$

${\pi }_i+(x_l-1){\overline{\pi }}_i\le z_{li}\le {\pi }_i+(x_l-1){\underset{¯}{\pi }}_i(68)$

$x_l{\underset{¯}{\pi }}_j\le z_{lj}\le x_l{\overline{\pi }}_j(69)$

${\pi }_j+(x_l-1){\overline{\pi }}_j\le z_{lj}\le {\pi }_j+(x_l-1){\underset{¯}{\pi }}_j(70)$

$2z_{li}-2z_{lj}+{\pi }_j-{\pi }_i\ge \rho \cdot Q_{ij}^2(71)$

$-Q_{ij}^{\max }(1-x_l)\le Q_{ij}\le x_l{Q}_{ij}^{\max }(72)$

5 Distributed Algorithms Based on Bender’s Decomposition and Cut Plane

Since large-scale mixed-integer programming models are difficult to solve in an efficient and quick manner, this section develops a distributed algorithm based on Bender’s decomposition (Gharaei et al., 2019). The basic idea of the distributed algorithm is to group the constraints and variables of a complex problem into more minor and easily tractable master problems and subproblems, which is shown in Figure 4.

FIGURE 4

FIGURE 4. Structure of the distributed algorithm.

The section below details the master problem and the subproblem based on Bender’s decomposition and cut plane. The master problem with respect to the problem Eq. 9 is formulated as follows (Conejo et al., 2006):

$\min {\displaystyle \sum _{i\in {\Omega }_{p2h}}\frac{r{(1+r)}^{{\tau }_{i,\xi }}}{{(1+r)}^{{\tau }_{i,\xi }}-1}{Z}_i{c}_i}+\eta (73)$

$C_{ope}+{\displaystyle \sum {\lambda }_i^{\left(k\right)}(P_i^{\max }-P_i^{\max \left(k\right)}})\le \eta (74)$

This is subject to constraints Eq. 12–Eq. 15, Eq. A3, and Eq. 74.

The subproblem related to Eq. 9 is depicted below:

$\min C_{ope}(75)$

$P_i^{\max }=P_i^{\max \left(k\right)}:{\lambda }_i(76)$

The dual variables${\lambda }_i$with respect to equation constraint (79) are given following a colon. The subproblem is subject to constraints Eqs 18, 19, Eqs 21–29, Eqs 31–42, Eqs 43–59, Eqs 67–72, Eq. 76, and Eqs A1, A2.

The detailed interactions between the master problem and the subproblem are presented below.

Step 1. Initialization. Set $LB=-\infty ,UB=+\infty ,k=0$.

Step 2. Solve the aforementioned master problem Eq. 76 to obtain the optimal values $Z_i^{\ast },{\eta }^{\ast }$ and update lower bound $LB={\displaystyle \sum _{i\in {\Omega }_{p2h}}\frac{r{(1+r)}^{{\tau }_{i,\xi }}}{{(1+r)}^{{\tau }_{i,\xi }}-1}{\left(Z_i^{\ast }\right)}^k{c}_i}+{\left({\eta }^{\ast }\right)}^k$ .

Step 3. Call the commercial solver Gurobi to solve the subproblem in (78) and update the upper bound $\min \left(UB,{\displaystyle \sum _{i\in {\Omega }_{p2h}}\frac{r{(1+r)}^{{\tau }_{i,\xi }}}{{(1+r)}^{{\tau }_{i,\xi }}-1}{\left(Z_i^{\ast }\right)}^k{c}_i}+{\left(C_{ope}^{\ast }\right)}^k\right)$ .

Step 4. If $UB-LB\le \epsilon$, the optimal values are obtained and the algorithm terminates. Otherwise, update the iteration counter $k=k+1$, and go to Step 2.It is noted that if $k=0$, then the constraint Eq. 75 is not incorporated in the master problem. The flowchart of the distributed optimization algorithm to solve the two-stage stochastic dynamic MISOC program is shown in Figure 5.

FIGURE 5

FIGURE 5. Flowchart of the distributed algorithm.

6 Case Study and Discussion

6.1 Description of Test Systems

In this section, simulations are carried out to verify the effectiveness of the proposed MISOC model and solution algorithm. The configuration of the simulation test systems is presented in Figure 6, which is composed of a modified IEEE 33-bus distributed network and a Belgian 20-node gas network comprising two gas compressors, four gas sources, and 1 CCGT. The wind turbine capacity installed at bus 15 is 3 MW, and the other three 1 MW wind units are located at buses e18, e22, and e26. The installed capacity of the CCGT unit at bus e33 is 1 MW, which obtains its gas supply from node g8. In addition, two gas compressors are employed in the distribution gas network to manage the gas pressures in the gas network. The CCGT unit is the coupling device between the electricity and gas system. It is worth noting that the variables $\rho$ of pipelines are considered constants, calculated in terms of a permissible hydrogen concentration of 10 vol%.

FIGURE 6

FIGURE 6. Configuration of the hybrid renewable energy system.

The time-of-use electricity price is illustrated in Figure 7. The predicated electrical load, the wind output power, and the gas load of the hybrid renewable energy system are shown in Figure 8. Based on Monte Carlo simulation and scenario reduction techniques, these data were used to generate the 10 operation scenarios for the study. The system net load is depicted in Figure 9. With high penetration renewable energy connection, it can be seen that the peak-valley difference of the demand becomes greater. At the same time, the net load profile varies more steeply. Therefore, it is urgent to enhance the system’s flexibility and promote the consumption of renewable energy.

FIGURE 7

FIGURE 7. Time-of-use electricity price.

FIGURE 8

FIGURE 8. Predicted values of electrical load, wind output, and the gas load.

FIGURE 9

FIGURE 9. Power demand load, wind power, and the net load of the system.

The other parameters used in the simulation are depicted in Table 1, and some of them can be found in Ref. Menon, 2005.

TABLE 1

TABLE 1. Other parameters used in the simulation.

6.2 Simulation Results

In order to verify the performance and effectiveness of the two-stage MISOC model and the distributed solution algorithm proposed in this article, two cases are designed for comparative analysis. The simulations are carried out on a laptop with Intel i7-9750H CPU and 16 GB memory.

Case 1. Basic scenario, optimal operation of the hybrid renewable energy system without P2H. The case mainly analyzes the total operation costs without P2H.

Case 2. Comprehensive scenario, coordinated operation, and planning of the hybrid renewable energy system with P2H, the flexibility, and hydrogen injections are also taken into account. The case primarily discusses the total investment and operation costs of the system, and it also focuses on the impact of P2H on the system flexibility.

6.2.1 Planning and Operation Result Analysis

The annual total costs of the two cases are provided in Table 2, including investment costs, purchased electricity costs, purchased gas costs, renewable energy curtailment costs, the electricity load-shedding cost, and the gas load-shedding cost. A comparison of the data in the sixth and seventh columns of Table 2 shows that after the construction of P2H, the renewable output curtailment cost drops from M$ (million dollars) 0.6524 to M$0.027705, and the electricity demand shedding cost decreases from M$0.21681 to M$0. It can also be seen that the annual total cost of case 2 is reduced by M$0.794 compared to case 1 by reducing the renewable energy curtailment costs and the load-shedding costs. Although the capital costs of P2H exist in Case 2, the annual total cost is still lower than that of case 1. In the meantime, we conclude that the abandoned wind power reaches 0.1814 and 0.1981 MW in Case 2, which occurs at wind turbines 1 and 2 at 3 a.m. 4 a.m., respectively. This phenomenon occurs largely because the wind speed is higher, and wind turbines produce more power when the electrical load is lower.

TABLE 2

TABLE 2. Comparison of planning and operation costs.

The planning results of P2H in Case 2 are shown in Table 3. It is shown that the installed capacity of P2H at buses 15, 18, and 26 is 469.7, 803, and 189 kW, respectively. The total investment cost of P2H is M$1.8998. The life cycle of P2H is assumed to be 10 years. Therefore, the annual investment cost is M$0.28313, considering the interest rate. It is worth noting that the capacity of the wind generator equipped at bus 15 is 3 MW, which is larger than the other three. Meanwhile, the electricity demand in this area is high, and the wind power output is mainly supplied to the electrical consumers. Thus, the capacity of P2H installed at bus 15 is lower than that of bus 18.

TABLE 3

TABLE 3. Case 2 results.

Take scenario 5 as an example. The optimal generator outputs of two cases in the hybrid renewable energy system are shown in Figures 10, 11. In this scenario, the renewable energy penetration level is 40%. As demonstrated from Figures 4, 10, wind generators, the distribution system operator, and the CCGT unit supply the electrical load demand to the end users, and the electrical load shedding mainly occurs at 11 p.m. P2H participates in the operation of the system as a flexible unit in Case 2. As shown in Figure 11, in Case 2, four wind generators, one CCGT unit, and the distributed system operator serve as the source to supply electricity to the consumers and P2H. The higher wind speeds from 1 to 7 a.m. make the wind turbine generators output more power. However, the power demand of end users is low during this time. P2H produces hydrogen using surplus wind power to reduce renewable curtailment. At 8 in the morning, the electrical load demand increases sharply, and the output of the wind turbines decreases, at which point the CCGT unit starts to supply the electrical load. The electrical load peaks at 10 a.m. and 11 a.m. When the CCGT unit output reaches its maximum power, P2H does not consume power. The electricity price is high at this stage, so the electricity purchase is less. With the gradual reduction of power load demand, the output power of the CCGT unit decreases. From 10 to 12 p.m., the electricity load gradually decreases. On the contrary, the wind generators’ output increases. The P2H units consume surplus renewable energy to produce hydrogen during this period, and the CCGT generator does not supply electricity to consumers.

FIGURE 10

FIGURE 10. Optimal generator outputs of case 1.

FIGURE 11

FIGURE 11. Optimal generator outputs of case 2.

As can be seen from Figure 11, case 2 yields smaller load-shedding and renewable energy curtailment than case 1. P2H plays a vital role in promoting the consumption of renewable energy.

The power output of the CCGT unit, the power purchased from the grid, and the four wind generators of case 2 are shown in Figure 12.

FIGURE 12

FIGURE 12. Power output of each unit of case 2.

6.2.2 System Flexibility Analysis

This part investigates the hybrid renewable energy system’s flexible supply and demand balance. The effect and role of P2H in system flexibility supply are evaluated by comparing the results of the two cases.

The total upward/downward flexibility supply and demand of case 1 are comparatively studied in Figures 13, 14. As shown in Figure 13, the blue line represents the total upward flexibility supply, and the orange line indicates the upward flexibility demand. We can also learn from Figure 13 that the demand for upward flexibility is greater at 6, 7, and 8 a.m. and reaches its maximum at 7 a.m. At 7 in the morning, the upward flexibility adequacy is less than 0, which indicates that upward flexibility supply is less than upward flexibility demand. Figure 14 illustrates the total downward flexibility supply and demand of case 1. We can learn from the following figure that there is a greater need for downward flexibility from 8 to 11 p.m. During this period, the supply of downward flexibility has not been able to meet the demand for downward flexibility, and there has been a shortage of downward flexibility.

FIGURE 13

FIGURE 13. Total upward flexibility supply–demand of case 1.

FIGURE 14

FIGURE 14. Total downward flexibility supply–demand of case 1.

The total upward/downward flexibility supply and demand of case 2 are shown in Figures 15, 16, respectively. It is shown in Figure 15 that there is a greater need for upward flexibility during the period from 6 to 8 a.m. However, after the construction of P2H, the lack of upward flexibility no longer exists with the coordinated optimization of flexibility resources. As can be seen in Figure 15, the same phenomenon is present in downward flexibility. By comparing Figure 13 with Figure 15 and Figure 14 with Figure 16, we can conclude that the system flexibility can be effectively enhanced after the construction of P2H, which helps promote the consumption of renewable energy.

FIGURE 15

FIGURE 15. Total upward flexibility supply–demand of case 2.

FIGURE 16

FIGURE 16. Total downward flexibility supply–demand of case 2.

7 Conclusion

This article proposes a low-carbon-emission hybrid renewable energy system to utilize power, hydrogen, and natural gas in an integrated mode. The P2H element serves as a key multi-energy coupling device connecting the distributed power network and gas and hydrogen networks.

In order to cope with the uncertainty of power load and renewable energy output, multiple operating scenarios are generated based on Monte Carlo simulations and the fast forward algorithm. Then, a two-stage stochastic planning model is constructed to consider the hybrid energy system flexibility and hydrogen injections. The nonlinear, nonconvex terms in the model are simplified, and a MISOC programming model is constructed to reduce the computation burden and further improve the solving efficiency. The results obtained by solving the above MISOC model using the distributed algorithm validate the effect of the construction of P2H.

Simulation results show that P2H can effectively reduce renewable energy curtailment, thus promoting the renewable energy consumption. On one hand, the upward/downward flexibility of the system has been enhanced by coordinating P2H and other flexibility resources. On the other hand, by injecting the hydrogen produced by P2H into the gas-distributed pipeline, the cost of hydrogen is effectively decreased, while reducing the system’s CO₂ emissions.

Due to fluctuations in renewable energy output, the operational efficiency of electrolysis equipment is challenging to achieve in practice, and the actual available capacity of the equipment is reduced.

Data Availability Statement

The original contributions presented in the study are included in the article/Supplementary Material, further inquiries can be directed to the corresponding author.

Author Contributions

JW, writing the original draft, conceptualization, formal analysis, methodology and writing programs; PZ, writing, review and editing, conceptualization and investigation; YL, funding acquisition, project administration, supervision; JL, conceptualization, writing, review and editing.

Funding

This work was supported by funds provided via the National Key Research and Development Program of China (Grant No. 2018YFE0208400). The authors are also grateful for the support from the Science and Technology Project of State Grid Corporation of China (Key Technologies of Novel Integrated Energy System Considering Cross-border Interconnection).

Conflict of Interest

The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.

Publisher’s Note

All claims expressed in this article are solely those of the authors and do not necessarily represent those of their affiliated organizations or those of the publisher, the editors, and the reviewers. Any product that may be evaluated in this article or claim that may be made by its manufacturer is not guaranteed or endorsed by the publisher.

References

Agbonaye, O., Keatley, P., Huang, Y., Ademulegun, O. O., and Hewitt, N. (2021). Mapping Demand Flexibility: A Spatio-Temporal Assessment of Flexibility Needs, Opportunities and Response Potential. Appl. Energ. 295, 117015. doi:10.1016/j.apenergy.2021.117015

CrossRef Full Text | Google Scholar

Bobo, L., Venzke, A., and Chatzivasileiadis, S. (2021). Second-order Cone Relaxations of the Optimal Power Flow for Active Distribution Grids: Comparison of Methods. Int. J. Electr. Power Energ. Syst. 127, 106625. doi:10.1016/j.ijepes.2020.106625

CrossRef Full Text | Google Scholar

Bramstoft, R., Pizarro-Alonso, A., Jensen, I. G., Ravn, H., and Münster, M. (2020). Modelling of Renewable Gas and Renewable Liquid Fuels in Future Integrated Energy Systems. Appl. Energ. 268, 114869. doi:10.1016/j.apenergy.2020.114869

CrossRef Full Text | Google Scholar

Cloete, S., and Hirth, L. (1922). Flexible Power and Hydrogen Production: Finding Synergy between CCS and Variable Renewables. Energy 192, 116671. doi:10.1016/j.energy.2019.116671

CrossRef Full Text | Google Scholar

Conejo, A. J., Baringo Morales, L., Kazempour, S. J., and Siddiqui, A. S. (2016). Investment in Electricity Generation and Transmission. Springer Nat. doi:10.1007/978-3-319-29501-5

CrossRef Full Text | Google Scholar

Conejo, A. J., Castillo, E., Minguez, R., and Garcia-Bertrand, R. (2006). Decomposition Techniques in Mathematical Programming:Engineering and Science Applications. Springer, 558.

Google Scholar

Dagoumas, A. S., and Koltsaklis, N. E. (2019). Review of Models for Integrating Renewable Energy in the Generation Expansion Planning. Appl. Energ. 242, 1573–1587. doi:10.1016/j.apenergy.2019.03.194

CrossRef Full Text | Google Scholar

Deymi-Dashtebayaz, M., Ebrahimi-Moghadam, A., Pishbin, S. I., and Pourramezan, M. (2019). Investigating the Effect of Hydrogen Injection on Natural Gas Thermo-Physical Properties with Various Compositions. Energy 167, 235–245. doi:10.1016/j.energy.2018.10.186

CrossRef Full Text | Google Scholar

Dolci, F., Thomas, D., Hilliard, S., Guerra, C. F., Hancke, R., Ito, H., et al. (2019). Incentives and Legal Barriers for Power-To-Hydrogen Pathways: An International Snapshot. Int. J. Hydrogen Energ. 44, 11394–11401. doi:10.1016/j.ijhydene.2019.03.045

CrossRef Full Text | Google Scholar

Ge, P., Hu, Q., Wu, Q., Dou, X., Wu, Z., and Ding, Y. (2020). Increasing Operational Flexibility of Integrated Energy Systems by Introducing Power to Hydrogen. IET Renew. Power Generation 14, 372–380. doi:10.1049/iet-rpg.2019.0663

CrossRef Full Text | Google Scholar

Gea-Bermúdez, J., Jensen, I. G., Münster, M., Koivisto, M., Kirkerud, J. G., Chen, Y.-k., et al. (2021). The Role of Sector Coupling in the green Transition: A Least-Cost Energy System Development in Northern-central Europe towards 2050. Appl. Energ. 289, 116685. doi:10.1016/j.apenergy.2021.116685

CrossRef Full Text | Google Scholar

Gharaei, A., Karimi, M., and Hoseini Shekarabi, S. A. (2019). Joint Economic Lot-Sizing in Multi-Product Multi-Level Integrated Supply Chains: Generalized Benders Decomposition. Int. J. Syst. Sci. Operations Logistics 7, 309–325. doi:10.1080/23302674.2019.1585595

CrossRef Full Text | Google Scholar

Gils, H. C., Gardian, H., and Schmugge, J. (2021). Interaction of Hydrogen Infrastructures with Other Sector Coupling Options towards a Zero-Emission Energy System in Germany. Renew. Energ. 180, 140–156. doi:10.1016/j.renene.2021.08.016

CrossRef Full Text | Google Scholar

Glenk, G., and Reichelstein, S. (2019). Economics of Converting Renewable Power to Hydrogen. Nat. Energ. 4, 216–222. doi:10.1038/s41560-019-0326-1

CrossRef Full Text | Google Scholar

Hajiabbas, M. P., and Mohammadi-Ivatloo, B. (2020). Optimization of Power System Problems _Methods, Algorithms and MATLAB Codes. Springer 262. doi:10.1007/978-3-030-34050-6

CrossRef Full Text | Google Scholar

Heris, M.-N., Mirzaei, M. A., Asadi, S., Mohammadi-Ivatloo, B., Zare, K., Jebelli, H., et al. (2020). Evaluation of Hydrogen Storage Technology in Risk-Constrained Stochastic Scheduling of Multi-Carrier Energy Systems Considering Power, Gas and Heating Network Constraints. Int. J. Hydrogen Energ. 45, 30129–30141. doi:10.1016/j.ijhydene.2020.08.090

CrossRef Full Text | Google Scholar

Hu, G., Chen, C., Lu, H. T., Wu, Y., Liu, C., Tao, L., et al. (2020). A Review of Technical Advances, Barriers, and Solutions in the Power to Hydrogen (P2H) Roadmap. Engineering 6, 1364–1380. doi:10.1016/j.eng.2020.04.016

CrossRef Full Text | Google Scholar

Koltsaklis, N. E., and Dagoumas, A. S. (2018). State-of-the-art Generation Expansion Planning: A Review. Appl. Energ. 230, 563–589. doi:10.1016/j.apenergy.2018.08.087

CrossRef Full Text | Google Scholar

Li, J., Lin, J., Song, Y., Xing, X., and Fu, C. (2019). Operation Optimization of Power to Hydrogen and Heat (P2HH) in ADN Coordinated with the District Heating Network. IEEE Trans. Sustain. Energ. 10, 1672–1683. doi:10.1109/tste.2018.2868827

CrossRef Full Text | Google Scholar

Li, R., Chen, L., Yuan, T., and Li, C. (2020). Optimal Dispatch of Zero-Carbon-Emission Micro Energy Internet Integrated with Non-supplementary Fired Compressed Air Energy Storage System. J. Mod. Power Syst. Clean Energ. 4, 566–580. doi:10.1007/s40565-016-0241-4(2016).20Menon

CrossRef Full Text | Google Scholar

Liu, J., Zeng, P. P., Xing, H., Li, Y., and Wu, Q. (2020). Hierarchical Duality-Based Planning of Transmission Networks Coordinating Active Distribution Network Operation. Energy 213, 118488. doi:10.1016/j.energy.2020.118488

CrossRef Full Text | Google Scholar

Lu, Z., Li, H., and Qiao, Y. (2018). Probabilistic Flexibility Evaluation for Power System Planning Considering its Association with Renewable Power Curtailment. IEEE Trans. Power Syst. 33, 3285–3295. doi:10.1109/tpwrs.2018.2810091

CrossRef Full Text | Google Scholar

Menon, E. S. (2005), Gas Pipeline Hydraulics R1. Continuing Education and Development, Inc., 21.

Google Scholar

Mongird, K., Fotedar, V., Viswanathan, V., and Koritarov, V. (2019). Energy Storage Technology and Cost Characterization Report.

Google Scholar

Pan, G., Gu, W., Lu, Y., Qiu, H., Lu, S., and Yao, S. (2020). Optimal Planning for Electricity-Hydrogen Integrated Energy System Considering Power to Hydrogen and Heat and Seasonal Storage. IEEE Trans. Sustain. Energ. 11, 2662–2676. doi:10.1109/tste.2020.2970078

CrossRef Full Text | Google Scholar

Robinius, M., Welder, L., and Ryberg, D. S. (2017). “Power-to-hydrogen and Hydrogen-To-X:Which markets?Which Economic potential?Answers from the Literature,” in 14th International Conference on the European Energy Market (EEM).

Google Scholar

Singh, M. K., and Kekatos, V. (2020). Natural Gas Flow Solvers Using Convex Relaxation. IEEE Trans. Control. Netw. Syst. 7, 1283–1295. doi:10.1109/tcns.2020.2972593

CrossRef Full Text | Google Scholar

Welder, L., Ryberg, D. S., Kotzur, L., Grube, T., Robinius, M., and Stolten, D. (2018). Spatio-temporal Optimization of a Future Energy System for Power-To-Hydrogen Applications in Germany. Energy 158, 1130–1149. doi:10.1016/j.energy.2018.05.059

CrossRef Full Text | Google Scholar

Zablocki, A. (2019). Fact Sheet Energy Storage. Washington, DC: Environmental and Energy Study Institute.

Google Scholar

Appendix AModel Simplifications

Term $i_{ij}{U}_i$ in Eq. 20 can be relaxed to the following second-order cone programming constraints (Bobo et al., 2021):

$\Vert \begin{array}{l}2P_{ij,t,s}\\ 2Q_{ij,t,s}\\ i_{ij,t,s}-U_{i,t,s}\end{array}\Vert \le i_{ij,t,s}+U_{i,t,s}\forall t,\forall s(\mathrm{A1})$

Equations 43–48 are of the form $Z=\min (X,Y)$ , which can be dealt with the big M method

$\begin{array}{l}Y-X\le M\omega \\ X-Y\le M(1-\omega )\\ Z\ge Y-M\omega \\ Z\le X+M\omega \\ Z\ge X-M(1-\omega )\\ Z\le Y+M(1-\omega )\end{array}(\mathrm{A2})$

where $\omega$ is the binary variable, and $\omega =1$ indicates that $X\le Y$ .

The nonlinear term $I_i{P}_i^{\max }$ in the objective function (10) is the product of variables $I_i$ and $P_i^{\max }$, which can be equivalently replaced with the following constraints (Conejo et al., 2016):

$\begin{array}{l}{Z}_i\le I_i\cdot {\overline{P}}_i^{\max }\\ Z_i\ge I_i\cdot {\underset{¯}{P}}_i^{\max }\\ Z_i\le P_i^{\max }+{\underset{¯}{P}}_i^{\max }\cdot (1-I_i)\\ Z_i\ge P_i^{\max }-{\overline{P}}_i^{\max }\cdot (1-I_i)\end{array}(\mathrm{A3})$

where $Z_i$ is an auxiliary variable used to substitute for $I_i{P}_i^{\max }$. ${\overline{P}}_i^{\max },$ and ${\underset{¯}{P}}_i^{\max }$ are the upper and lower bound values of $P_i^{\max }$, respectively. It can be observed that when the binary variable $I_i$ equals to 0, $Z_i=0$, and when $I_i=1$, $Z_i=P_i^{\max }$ .

Nomenclature

Binary Variable

I The P2H device is built or not (1:‘built’, 0:‘not’)

x Gas flow direction, (1:‘the gas flow rate is positive’ and 0: ‘otherwise’)

Sets

$\Omega$ The initial set of all generated scenarios

${\Omega }_{cg}$ A set of gas nodes connecting CCGT

${\Omega }_s$ A set of selected scenarios

${\Omega }_J$ A set of remaining scenarios after the selected scenarios are removed

${\Omega }_{p2h}$ A set of candidate locations for P2H devices

${\Omega }_{ele}$ A set of nodes that purchase electricity

${\Omega }_{renew}$ A set of nodes equipped with renewable energy sources

${\Omega }_{source}$ A set of nodes equipped with natural gas supply units

${\Omega }_o$ A set of feasible solutions to the operation optimization problem

${\Omega }_{load}$ A set of electric demand nodes

${\Omega }_{gas}$ A set of gas demand nodes

${\Omega }_{gt}$ A set of electric nodes equipped with gas turbines

Continuous variables

$P_{load,t}$ $P_{wind,t}$ Actual electric demand and renewable energy output at time t

$P_{wind,t}^{fore}$ $P_{load,t}^{fore}$ The predicted output of renewable energy/demand at time t

$\Delta P_{wind,t}$ $\Delta P_{load,t}$ The prediction errors about renewable energy output and electric demand at time t

$P_{t,s}^{net},F_{t,s}^{D\_up}$ $F_{t,s}^{D\_down}$ The net load and the total upward/downward flexibility demand capacity at time t for the scenario s

$P_i^{\max }$ Capacity of candidate P2H(MW)

$P_j^{p2h}$ The consumed electricity by P2H

$\text{Δ}{P}_{k,t,s}^{load}$ The shedding electric load (MW)

$\text{Δ}{P}_{m,t,s}^{gas}$ The shedding gas load (m³)

$P_{i,t,s}^{ele}$ The power purchased from TSO

$P_{n,t,s}^{saurce}$ The output of the gas source n

$P_{j,t,s}^{gas,load}$ The gas load at node $j$ , at time $t$, for scenario $s$

$\text{Δ}{P}_{j,t,s}^{renew}$ The curtailment of renewable energy

$P_j^{renew},P_j^{CCGT},P_j^g$ The renewable energy output, the CCGT output, and the power purchase at node $j$

$P_{GT,t,s}^{gas}$ The gas consumed by CCGT

$Q_{ij}$ Gas flow rate through the pipeline $l_{ij}$

$P_{i,t,s}$ The pressure of gas node $i$ at time t under the scenario $s$

$P_{j,t,s}^{h2}$ , $P_{i,t,s}^{h2\_saurce}$ The output of hydrogen produced by P2H and the hydrogen coupled to the gas source

$\begin{array}{l}{F}_{p2h,t,s}^{s\_up},F_{GT,t,s}^{s\_up}\\ F_{buy,t,s}^{s\_up},F_{k,t,s}^{s\_up}\end{array}$ The upward flexibility capacity provided by P2H, CCGT, power purchase, and the shedding load

$\begin{array}{l}{F}_{p2h,t,s}^{s\_down},F_{GT,t,s}^{s\_down}\\ F_{buy,t,s}^{s\_down},F_{j,t,s}^{s\_down}\end{array}$ The downward flexibility capacity provided by P2H, CCGT, power purchase, and the renewable curtailment

$F_{t,s}^{up},F_{t,s}^{down}$ The system upward and downward flexibility adequacy

Parameters

$r,{\tau }_{i,\epsilon }$ The discount rate/life cycle of the device $\epsilon$ at node $i$

$c_i$ unit investment cost of candidate device $\epsilon$ ($/MWh)

${\pi }_s$ The probability of each scenario

$N_s$ The total number of all scenes

$C_{ij}$ The incident matrix related to the sites of P2H and gas sources