Probabilistic energy flow and risk assessment of Electricity–Gas systems considering the thermodynamic process

Energy xxx (xxxx) xxx

Contents lists available at ScienceDirect

Energy journal homepage: www.elsevier.com/locate/energy

Probabilistic energy flow and risk assessment of ElectricityeGas systems considering the thermodynamic process Shiyuan Bao a, Zhifang Yang a, *, Juan Yu a, **, Wei Dai a, Lin Guo a, Hongxin Yu b a

State Key Laboratory of Power Transmission Equipment & System Security and New Technology, College of Electrical Engineering, Chongqing University, Chongqing, 400030, China b State Grid Chongqing Electric Power Research Institute, Chongqing, 401123, China

a r t i c l e i n f o

a b s t r a c t

Article history: Received 29 March 2019 Received in revised form 25 July 2019 Accepted 1 October 2019 Available online xxx

The existing analyses of the electricityegas systems generally do not regard the nodal gas temperature as a state variable. However, the temperature variations caused by the thermodynamic process are important because of the following two reasons: 1) gas temperature influences the pressure and gas flow; 2) gas temperature is an important indicator for hydrate formation in gas networks, which greatly jeopardizes the gas transmission. We herein propose a probabilistic energy flow method considering the thermodynamic process in gas networks. To build the energy flow model, the algebraic pipeline flow model under the non-isothermal condition is derived according to a set of partial differential equations reflecting the thermodynamic process. The models of pressure-regulating stations and compression stations considering the thermodynamic processes are presented in an electricityegas analysis for the first time. Risk indices are proposed to quantify the risks of hydrate formation and state variables that exceed the limits. The simulation results of two test cases demonstrate the effectiveness of the proposed method. © 2019 Elsevier Ltd. All rights reserved.

Keywords: Electricityegas systems Probabilistic energy flow Thermodynamic process Pressure-regulating stations Hydrate formation

1. Introduction 1.1. Research motivations In recent years, the coupling relationship between electrical and natural gas systems has grown significantly with the increasing utilization of facilities such as gas-fired units [1], CHP units [2] and P2G devices [3]. Thus, the unified analysis of electricityegas systems is imperative. Owing to massive uncertainties in electricityegas systems [4], probabilistic energy flow analysis becomes an effective tool to reflect the violation risks of state variables such as voltage and gas pressure. The components of gas networks undergo thermodynamic processes when transfers of heat, work, and energies occur within them. In the thermodynamic process, the gas temperature and pressure will change [5]. However, the existing studies regard gas temperature as a constant by ignoring the thermodynamic process. The thermodynamic process should be properly considered for

* Corresponding author. ** Corresponding author. E-mail addresses: [email protected] (Z. Yang), [email protected] (J. Yu).

three main reasons: 1) The gas temperature significantly differs from one place to another within a system. For example, the temperature at the beginning of the Central Asia-China Gas pipeline is 45
C, which is determined by the gas source; however, it is close to the ambient temperature at the end of the pipeline, which is lower than 5
C in the winter. Hence, proper consideration of thermodynamic process helps to accurately describe the distribution of temperature within a gas network. 2) The gas pressure and temperature are tightly coupled with each other. Hence, properly considering the thermodynamic process in the energy flow model also plays an important role in avoiding considerable errors when calculating gas pressures. 3) The formation of hydrates results in several adverse effects to an electricity-gas network, and a simple way to verify hydrate formation requires the information of gas temperature and pressure simultaneously [6]. The adverse effects are explained as follows. First of all, the hydrates formed in pipelines interrupt the gas transmission [7]. Considering the interdependencies between electricity and natural gas systems, the reliability of power supply will be reduced when the interruption happens at the branch supplying gas to power systems. Secondly, severe damage to components can be caused by hydrates [6]. On the one hand, the hydrate blockage will cause a local high pressure

https://doi.org/10.1016/j.energy.2019.116263 0360-5442/© 2019 Elsevier Ltd. All rights reserved.

Please cite this article as: Bao S et al., Probabilistic energy flow and risk assessment of ElectricityeGas systems considering the thermodynamic process, Energy, https://doi.org/10.1016/j.energy.2019.116263

2

S. Bao et al. / Energy xxx (xxxx) xxx

Nomenclature

Indices i, j m, n p, r, c

Index of electricity system nodes Index of natural gas system nodes Index of pipelines, pressure-regulating stations and compression stations in a natural gas system

Parameters and Variables Amn Cross-section area of pipeline mn Bij Imaginary elements at the ith row and jth column of bus admittance matrix Constant related to work efficiency of compression Bmn k station connecting node m to node n Energy consumption of compression station BHPmn connecting node m to node n c Scale parameter of a Weibull distribution CP Gas heat capacity at constant pressure Cmn Discharge coefficient of pressure-regulating station connecting node m to node n Dmn Inner diameter of pipeline mn FD,n Normal gas load at node n FG,n Gas source injection at node n FGAS,n Gas load of gas-fired generators at node n Gas flow of compression station connecting node m Fcmn to node n Gas flow of pipeline mn Fpmn Frmn Gij GHV k Kvmn Lmn Mpmn N Nmax Ne Nm p pk,n pl pn pu p0 PD,i PG,i PGAS,i PW,i PHFn PNHPVn PNHTVn PNLPVn QC,i

Gas flow of pressure-regulating station connecting node m to node n Real elements at the ith row and jth column of bus admittance matrix Gross heating value of gas-fired generators Shape parameter of a Weibull distribution Flow coefficient of pressure-regulating station connecting node m to node n Length of pipeline mn Mass flow of pipeline mn Number of samples already calculated in a Monte Carlo simulation process Maximum number of samples in a Monte Carlo simulation process Number of nodes in the power system Number of nodes in the natural gas system Natural gas pressure The kth element in set Pn Lower limit of nodal gas pressure Gas pressure at node n Upper limit of nodal gas pressure Pressure of the standard condition (101 kPa) Active power load at node i Active power output of coal-fired generator at node i Active power output of gas-fired generator at node i Active power output of wind farm at node i Probability of Hydrate Formation at node n Probability of Node High Pressure Violation at node n Probability of Node High Temperature Violation at node n Probability of Node Low Pressure Violation at node n Parallel reactive compensator output at node i

QD,i QG,i

Reactive power load at node i Reactive power output of coal-fired generator at node i QGAS,i Reactive power output of gas-fired generator at node i R Gas constant of natural gas T Natural gas temperature Tav Average gas temperature TG,n Temperature of gas source at node n Tk,n The kth element in set Tn Tn Final temperature of a mixture of several gas injections at node n Ts Temperature of the ambient environment Tu Upper limit of nodal gas temperature T0 Temperature of the standard condition (273.15 K) Gas temperature at node n of branch (pipeline, Tnmn compression station, or pressure-regulating station) mn T(x) Gas temperature at a distance of x from the inlet end of a pipeline Umn Heat transfer coefficient of pipeline mn Vi Voltage magnitude at node i x Distance from the inlet end of a pipeline Z Compressibility factor natural gas ai,n, bi,n, gi,n Energy conversion parameters of gas-fired generators connecting electrical node i to gas node m amn, bmn,gmn Energy conversion parameters of compression station connecting node m to node n d Standard derivation of a probability distribution hJT JouleeThompson coefficient qi Voltage angle at node i qij Voltage angle difference between node i and node j l Friction coefficient m Expectation of a probability distribution r Natural gas density r0 Natural gas density under the standard condition (T0 ¼ 273.15 K, p0 ¼ 101 kPa) ttmn Gas consumption of compression station connecting c node m to node n y Natural gas velocity 4mn Ratio of the regulator’s orifice diameter and the pipeline’s internal diameter for pressure-regulating station connecting node m to node n cmn Polytropic index of compression station connecting node m to node n Functions and Sets Fflow($) Expression of the flow model of a gas component Fthermal($) Expression of the thermal model of a gas component pcritical(Tk,n) Critical gas pressure of hydrate formation when the gas temperature is Tk,n, Pn Set of all possible values for pressure at node n P(a) Probability of event a P(a|b) Conditional probability of event a given that event b has occurred Tcritical(pk,n) Critical gas temperature of hydrate formation when the gas pressure is pk,n Tn Set of all possible values for temperature at node n

Please cite this article as: Bao S et al., Probabilistic energy flow and risk assessment of ElectricityeGas systems considering the thermodynamic process, Energy, https://doi.org/10.1016/j.energy.2019.116263

S. Bao et al. / Energy xxx (xxxx) xxx

at the upstream of the block point, which leads to a pipeline burst. On the other hand, the movement of hydrates in gas networks can cause serious deformation to downstream elements such as valves, nozzles and so on. Thirdly, the pipeline burst and gas leak is a potentially fatal peril to field workers [6]. Last but not the least, it is costly to solve the problem of hydrates. Field workers usually handle the hydrate blockage by releasing the gas, which causes a huge waste of gas sources. Generally speaking, the cost used to prevent hydrate formation accounts for 5%e8% of the total production cost [8]. In addition, the temperature and pressure decrease simultaneously in pressure-regulating stations, thus aggravating the risk of hydrate formation [9]. In conclusion, the gas temperature variation caused by the thermodynamic process is crucial in precisely analyzing the states of natural gas systems. To address the problems above, we propose a probabilistic energy flow analysis method properly considering the thermodynamic process. The proposed method achieves this goal by: 1) taking gas temperature as a variable during the whole modelling process; 2) expanding the application range of the model. The influence of the gas temperature on the probability characteristics of energy flow is represented, and the risks such as hydrate formation are quantified. 1.2. Literature review As a common way to deal with the massive uncertainties in power systems, the probabilistic power flow analysis is studied both in modelling [10] and in fast calculation method [11]. Probabilistic energy flow analysis in electricityegas systems has been investigated in several studies. The analysis can be divided into two steps: the establishment of energy flow model and the calculation of risk indices. For the probabilistic energy flow model: Ref. [12] considers the gas-fired units as the coupling components of electricity and gas networks, and a piecewise linearization method is applied to achieve a better fitting precision of probability distribution; The electricity and gas networks in Ref. [13] are coupled by gas-fired units and moto-compressors, and the “removing compressor technique” is proposed to simplify the topological structure of the gas networks; Three aspects of couplings are considered in Ref. [14], i.e. the gas-fired units, the electric-driven compressors, and the energy hubs integrated with P2G devices. However, they focus on the nodal flow balance with gas pressure and flow as state variables. In these models, the gas temperature is always equal to the ambient temperature while the thermodynamic process in gas networks is ignored. Ref. [15] considers the gas temperature variation in pipelines for the first time, adding nodal thermal balance in the energy flow model. Based on it, Ref. [16] further studies the energy flow calculation considering the gas temperature variation, in which distributed slack node technique is applied in electricity networks. However, a contradiction occurs in that the flow model of the pipelines is derived with gas temperature as a constant. The accuracy of the analysis is thereby limited. For studies focusing on natural gas systems, the thermodynamic process in gas pipelines is observed by obtaining the gas temperature and pressure simultaneously. A set of partial differential equations is solved by numerical calculation methods such as the finite difference approximations for the partial derivatives [17] or the RungeeKutta method [18]. However, for the energy flow model of electricityegas systems, the partial differential equations are no longer feasible. The important components in natural gas systems also include pressure-regulating stations and compression stations. Pressureregulating stations fix the outlet pressure to satisfy the requirements of terminal gas users. In Ref. [19], the max and min gas expansion ratios of the pressure-regulating stations are mentioned

3

as constraints of the proposed optimization model. Ref. [20] classifies the pressure-regulating stations and compression stations as active branches that can enforce the given pressure on a specified side. However, the detailed model in electricityegas studies are non-existent. The influence of pressure-regulating stations on the nodal temperature and pressure of the system cannot be considered effectively. Meanwhile, the temperature increase [21] in a compression station caused by the thermodynamic process is also ignored in studies of electricityegas systems. Ignoring the thermodynamic process in pressure-regulating stations and compression stations can yield significant errors to the results of probabilistic energy flow analysis. It is noteworthy that the gas model proposed in the unified analysis of electricity-gas networks is much simpler than the model based on fluid dynamics. The reasons are as follows: 1) The unified analysis of electricity-gas systems is expanded from the power system analysis and aims at improving the quality of power system management. For natural gas systems, only some basic macroscopic properties need to be considered, such as the quantity balance of gas supplied, the gas pressure and temperature. Hence, the simplified model is sufficient for the application scenario. 2) As the power system is often modeled as a set of static nonlinear equations, the dynamic behaviors of both electricity and gas networks (which are normally depicted by differential equations) are hard to be considered. In conclusion, considering the necessity and difficulty of modelling based on fluid dynamics, the energy flow model proposed in this paper is simplified to algebraic. For risk indices, in Ref. [22], the probabilities of voltage and load exceeding the limits in electricity networks can be obtained; similarly, the probability of pressure exceeding the limits in gas networks can be calculated by the results in Ref. [14]. However, indices quantifying the risk of hydrate formation (which requires a comprehensive analysis of gas temperature and pressure) are lacking. There are many methods proposed to verify the hydrate formation in a deterministic natural gas system. Ref. [6] summarizes the gas-gravity plot and the empirical correlation method. In these methods, the composition of natural gas is assumed to be unchanged during the transmission process. The composition of natural gas exploited from different sources varies, and different composition corresponds to different critical condition curve depicted by gas pressure and temperature. In some other methods, the change in gas composition during the process of hydrate formation is considered. In Ref. [23], the fugacity of the gas phase (which associates the chemical potential of each gas component with other variables) is obtained by different equations of state to determine the hydration pressure and temperature. The chemical potential of each component in the natural gas at equilibrium stays the same in each phase. Ref. [24] points out that the hydrate formation is a complex multiphase flow problem and uses the CFDsoftware Fluent to analyze it. To model the effects of hydrate condensation and accumulation, the 3-D CFD method is used. Ref. [25] establishes a new model consisting of mass, momentum and energy balance equations incorporated the hydrate formation process. Thus, the temperature, pressure and diameter in the pipeline will change over time. The hydrate formation region and the hydrate blockage degree at different transport time are determined through this method. Comparing all the methods mentioned above, the gas-gravity plot and the empirical correlation method are relatively rough. However, they are suitable for the unified analysis of electricity-gas systems. The reason is that the gas model is described by a set of nonlinear equations, which are convenient for the analysis of electricity-gas systems providing gas pressures and temperatures; once the gas composition is determined, the hydrate formation can be verified only by that information.

Please cite this article as: Bao S et al., Probabilistic energy flow and risk assessment of ElectricityeGas systems considering the thermodynamic process, Energy, https://doi.org/10.1016/j.energy.2019.116263

4

S. Bao et al. / Energy xxx (xxxx) xxx

1.3. Contributions

2.1. Derivation of the pipeline flow model

This study proposes a probabilistic energy flow analysis method of electricity-gas systems. The primary contributions of this paper can be summarized as follows:

The derivation starts from (1)e(2), which are transformed from a set of partial differential equations representing the primary conservation principles of gas flow (mass conservation, momentum conservation, and energy conservation) and the equation of state [26]. For a horizontal pipeline mn, equation (1) is the flow model in differential form considering the thermodynamic process, which describes the relationship among the temperature, pressure, and flow. It is obtained from the equation depicting the momentum conservation by the following steps: 1) regarding the pipeline elevation as 0; 2) substituting the gas density r and gas velocity y by y ¼ Mpmn /

1) A probabilistic energy flow model of electricityegas networks is proposed. In this model, a pipeline flow model under the nonisothermal condition is derived to reflect the thermodynamic process. Additionally, the thermal model of pressure-regulating stations and compression stations are incorporated in electricityegas systems for the first time. The gas temperature is regarded as a state variable, and the other state variables influenced by the coupling relationship among them can be represented. 2) A series of indices are proposed to assess the operation condition of the electricityegas systems. Particularly, we present the Probability of Hydrate Formation based on the probabilistic energy flow results to reflect the influence of temperature on the operational security. According to these indices, the risk level of hydrate formation and the quality of the nodal gas temperatures and pressures in an electricityegas system can be quantified. The remainder of this paper is organized as follows. Section II presents the models of gas components considering the thermodynamic processes in them. The probabilistic energy flow model of electricityegas systems and risk indices are proposed in Section III. Simulation results are presented in Section IV, and Section V concludes the discussion. 2. Component models considering thermodynamic processes in gas networks Component models are the key in probabilistic energy flow calculations. Unlike the existing models, this study considers the thermodynamic processes in different components to ensure the accuracy of the calculation results. For pipelines, the gas temperature at the inlet end (Tm) is significantly different from that at the outlet end (Tnmn ). This phenomenon is caused by the thermodynamic process in the pipelines. Meanwhile, the states of variables (Tm, Tnmn , pm, pn, and Fpmn ) are interactive with each other owing to the thermodynamic process, which is depicted by the pipeline model, which is shown in Fig. 1. Similar to pipelines, the thermodynamic processes in pressure-regulating stations and compression stations must be considered to represent the interaction among the state variables. In this section, the detailed models of the components considering the thermodynamic process are established as follows. The model of each component consists of a flow model and a thermal model.

(Amnr) ¼ Mpmn ZRT/(Amnp) (which is transformed from the equa-

tion of state and the equation depicting mass conservation). Equation (2) originating from the energy conservation represents the gas temperature at a distance of x from the inlet end of a pipeline [27]. The purpose of this derivation is to obtain an algebraic flow model of pipelines while preserving the thermodynamic characteristic.

dp dx  ZRT ¼ l mn p 2D

M mn p ZRT

!2 þ

Amn p

TðxÞ ¼ Ts þ ðTm  Ts Þea

mn

x

 hJT

M mn 1 p ZRT d 2 Amn p

!2 (1)

 mn ðpm  pn Þ  1  ea x mn mn a L

(2)

where

amn ¼

pDmn U mn F mn p r0 CP

The difference between this study and the existing models is that apart from gas pressure p, the temperature T in (1) is also regarded as a variable and retained in the differential sign d($), as shown in (3). The reasons for this operation are as follows: 1) the gas temperature changes significantly along the pipeline, 2) the gas pressure interacts with the temperature, and 3) the risk of hydrate formation depends on gas temperature and pressure simultaneously.



dp dx ZRT ¼ l mn p 2D  2 M mn p ZR ¼  2 2 Amn p

M mn p ZRT Amn p

lT 2

!2 þ

mn 1 M p ZR 2 Amn

2

!2 d !

T2

!

p2 (3)

dx 2T dp þ 2TdT  p Dmn

In equation (3), the third term on the right side (2T2dp/p) are ignored because it is much smaller than the first term and is generally ignored in engineering practice. Then, by substituting (2) into the T in front of dx, (3) can be transformed into (4):



2 Mmn ZR  l  p mn Ts þ ðTm  Ts Þea x pdp ¼  2 Dmn 2 Amn  mn ðpm  pn Þ   hJT mn mn 1  ea x dx þ 2dT a L

(4)

In (4), Mpmn ¼ Fpmn r0, R ¼ p0/(r0T0) under the standard condition, which is generally used in engineering practice. Fig. 1. Diagram of the thermodynamic process in a pipeline.

Please cite this article as: Bao S et al., Probabilistic energy flow and risk assessment of ElectricityeGas systems considering the thermodynamic process, Energy, https://doi.org/10.1016/j.energy.2019.116263

S. Bao et al. / Energy xxx (xxxx) xxx

5

Finally, integrating between the inlet and outlet ends:

p2m



p2n

¼

F mn p

!2

r0 Zp0

Amn



   h  mn  amn Lmn ðTm  Ts Þ þ hJT ðpm  pn Þ  JT mn amn Lmn 1  e  T L þ ðp  p Þ þ 2 T  T s m n m n Dmn amn ðamn Þ2 Lmn

T0

l

Thus far, the algebraic pipeline flow model (5) under the nonisothermal condition is derived. Compared with the existing models, the gas temperature is regarded as a variable in the whole derivation process; thus, the maximum preservation of the thermodynamic characteristic in pipelines is achieved. Meanwhile, the thermal model (6) can be obtained from (2) when x equals Lmn, which comprises the pipeline model together with the flow model: a T mn n ¼ Ts þ ðTm  Ts Þe

mn mn

L

 hJT

 ðpm  pn Þ  mn mn 1  ea L mn mn a L

(6)

The models of pressure-regulating stations and compression stations are first applied to the probabilistic energy flow calculations. The detailed models are introduced in the following: For the pressure-regulating station model, the thermal model (7) [9] represents the relationship between pressure loss and temperature decrease caused by the JouleeThomson effect (the thermodynamic process in pressure-regulating stations). The flow model (8) [28] represents the relationship between gas flow and state variables at nodes, i.e., nodal gas temperature and pressure (when pn/pm is larger than 0.5).

rffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi  ffi 1  ð4mn Þ4 1  ðCmnÞ2  C mn ð4mn Þ2 Tm  Tnmn ¼ hJT rffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi ðpm  pn Þ  ffi 1  ð4mn Þ4 1  ðCmnÞ2 þ C mn ð4mn Þ2 (7) sffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi ðpm  pn Þpn r0 Tm

(8)

For the compression station model, the thermal model (9) [21] and flow model (10) and (11) [29] are shown as below. The compression is between two ideal situations (isothermal compression and adiabatic compression).

FG;n TG;n  "  FG;n 

X

¼ Tm

cmnmn1 c

(9)

 mn BHP mn ¼ Bmn k F c Tm

 3 2  Z cmnmn1 c 7 6 pn  15 4 cmn  1 pm

cmn



mn tmn þ bmn BHP mn þ gmn BHP mn c ¼a

2

(10)

(11)

3. Probabilistic energy flow method of electricity-gas systems Based on the presented component models considering the thermodynamic process in gas networks, a probabilistic energy flow model and a series of risk indices are proposed in this section.

3.1. Probabilistic energy flow model The modeling of probabilistic energy flow includes the modeling of nodal balance equations and random factors. The nodal balance equations (12) and (16) include the corresponding equations of natural gas systems, electrical systems, and coupling components in an electricityegas system. For natural gas systems, the nodal thermal balance equation (12) implies that the temperature of gas flowing out from a node is the weighted mean of temperatures injected to the node; the nodal flow balance equation (13) implies that the flow into a node is equal to that out from the same node. For electrical systems, the active and reactive power flow balances at each node are represented in (14) and (15), respectively. For the coupling components, i.e., the gas-fired generators herein, the conversion relationship between gas flow and electrical power is represented in (16).

  X X     mn mn mn sgn2 m; n F mn sgn2 m; n F mn sgn2 m; n F mn p Tn þ c Tn  r Tn

m2n

X

pn pm

The models of the components (including pipelines, pressure regulating-stations, and compression stations) in natural gas systems are established. The thermodynamic process in gas networks are considered comprehensively and accurately.

2.2. Models of pressure-regulating stations and compression stations

mn F mn r ¼ 514K v

 T mn n

(5)

m2n

m2n

  X X    

sgn2 m; n F mn sgn2 m; n F mn sgn2 m; n F mn Tn ¼ 0; p þ c  r

m2n

m2n

" FG;n  FGAS;n  FD;n 

X

F mn p þ

m2n

X m2n

(12) cn ¼ 1; 2; …; Nm

m2n

sgn1 ðm; nÞF mn c þ

X

sgn2 ðm; nÞtmn c 

m2n

X

# sgn1 ðm; nÞF mn ¼ 0; r

cn ¼ 1; 2; …; Nm

(13)

m2n

Please cite this article as: Bao S et al., Probabilistic energy flow and risk assessment of ElectricityeGas systems considering the thermodynamic process, Energy, https://doi.org/10.1016/j.energy.2019.116263

6

S. Bao et al. / Energy xxx (xxxx) xxx

PG;i þ PW;i þ PGAS;i  PD;i  Vi

Ne X

  Vj Gij cosqij þ Bij sinqij

j¼1

¼ 0;

c i ¼ 1; 2; …; Ne

QG;i þ QGAS;i þ QC;i  QD;i  Vi

(14) Ne X

  Vj Gij sinqij  Bij cosqij

j¼1

¼ 0;

FGAS;n ¼

ci ¼ 1; 2; …; Ne

ai;n þbi;n PGAS;i þ gi;n P 2GAS;i GHV

(15)

(16)

where

 sgn1 ðm; nÞ ¼  sgn2 ðm; nÞ ¼

þ1; pm < pn 1; pm > pn þ1; pm < pn 0; pm > pn

m2n denotes that node m is adjacent to node n. For the pipelines and pressure-regulating stations, the pressure at the inlet node m is higher than that at the outlet node n, while the situation is opposite for the compression stations. For the modeling of random factors, the uncertainties of wind power generation PW,i as well as electrical load PD,i, QD,i and gas load FD,n are considered. The uncertainties of wind speeds are modeled as Weibull distributions. A Weibull distribution can be depicted by the shape parameter k and scale parameter c [30], which can be estimated by many ways including the graphic method, the maximum likelihood method, the moment method and so on [31]. Without the loss of generality, the variations in loads are modeled as normal distributions. Other distributions of loads, such as those considering the correlation between electrical loads and wind intensity [32], can be handled as well. So far, the probabilistic energy flow model are established (5) and (16), and the state variables are X ¼ [Tn, pn, qi, Vi]T. Through this model, the uncertainties in the electricity-gas systems can be analyzed.

Fig. 2. Critical condition curve of hydrate formation.

p

¼

2Pn

k;n X

      P Tn ¼ Tk;n ,P pn > pcritical Tk;n Tn ¼ Tk;n

Tk;n 2Tn

(17) Other risk indices: State variables exceeding the limits can yield problems (such as the damage to the performance of anti-corrosive materials in system facilities, the burst of pipelines, and the decrease in transmission capacity) to the safe operation of electricityegas systems. Hence, PNHTVn, PNHPVn, and PNLPVn are established to quantify the risks of gas temperature exceeding high limits, gas pressure exceeding high limits, and gas pressure exceeding low limits at node n, respectively. These risk indices are expressed as:

PNHTVn ¼

X

  P Tn ¼ Tk;n

(18)

  P pn ¼ pk;n

(19)

  P pn ¼ pk;n

(20)

T k;n 2Tn Tk;n > Tu

PNHPVn ¼

X pk;n 2Pn pk;n > pu

3.2. Risk indices After calculating the results of the probabilistic energy flow model, a series of risk indices quantifying the risks of hydrate formation and state variables exceeding the limits are further proposed as follows: Probability of Hydrate Formation at node n (PHFn): This index is proposed to reflect the risk of hydrate formation at each node of an electricityegas system. Assuming that the composition of gas is unchanged during the hydrate formation process, whether hydrates are formed or not only depends on the temperature and pressure. Using the empirical correlation in Ref. [33] as the criterion, once the composition of gas is determined, the corresponding critical condition curve is selected, which is illustrated as Fig. 2. It can also be seen that when states are near to the critical state, the formation of hydrates is sensitive to both pressure and temperature. Therefore, even small error in state variables can lead to wrong assessment results. When the states of variables are located at the left upper side of the curve (with high pressure and low pressure), hydrates are formed; meanwhile, the states at the right lower side are safe. Based on the facts above, the PHFn can be expressed as:

       P pn ¼ pk;n ,P Tn < Tcritical pk;n pn ¼ pk;n

X

PHFn ¼

PNLPVn ¼

X pk;n 2Pn pk;n < pl

3.3. Process of the probabilistic energy flow analysis method The calculation process of the probabilistic energy flow method can be summarized in the flow chart shown in Fig. 3. The detailed solution steps are as follows: Step 1: Inputting initial system data. In this step, the unchanged topological structure of the electricity-gas network as well as the parameters of wind power generation and electricity/gas load probability distributions are obtained. N is initially set as N ¼ 1. Step 2: Sampling state. In this step, the random variables (the electrical loads, gas loads, and wind speeds) are sampled by the Monte Carlo method, then a

Please cite this article as: Bao S et al., Probabilistic energy flow and risk assessment of ElectricityeGas systems considering the thermodynamic process, Energy, https://doi.org/10.1016/j.energy.2019.116263

S. Bao et al. / Energy xxx (xxxx) xxx

In this paper, the maximum number of samples Nmax is taken as the convergence criterion of the Monte Carlo method. If the number of samples already calculated N is less than Nmax, the solution process returns to step 2 and the value of N updates by N]Nþ1; otherwise, the process stops and the values of indices are printed as the final results.

Start Input initial system data, N=1 Do sampling for random variables and obtain a system state

4. Case studies 4.1. System description and models compared

Calculate the unified steady-state energy flow model

No

Is the power flow converged?

N=N+1

Yes Calculate the proposed risk indices

Is the convergence criterion met?

7

No

Yes Print the results End Fig. 3. Probabilistic energy flow analysis process for electricityegas systems using the proposed method.

deterministic system state is obtained. Step 3: Calculating energy flow. A unified steady-state system energy flow model (12) and (16) is established on the basis of component models (5) and (11). The whole model which describes the energy flow of a deterministic system state is a set of nonlinear equations and is solved by the Newton-Raphson method [34]. Then the results of the sampled system state are acquired. Step 4: Verifying the convergency of energy flow calculation. Considering the shortcomings of the model and the NewtonRaphson method, it is possible that the energy flow calculation cannot be converged. Therefore, the purpose of this step is to verify the convergence of the energy flow calculation conducted in Step 3. If it is converged, effective results are obtained and the proposed risk indices can be calculated in the next step; otherwise, the process goes straight to Step 6. Step 5: Calculating the proposed risk indices. Based on the effective results of energy flow calculated in Step 3 after each sampling, the proposed risk indices, i.e., PHFn, PNHTVn, PNHPVn, and PNLPVn (17) and (20) can be calculated. As the sample-calculation-sample-calculation process continues, the values of indices are updated. Step 6: Verifying the convergence of the Monte Carlo method.

The proposed probabilistic energy flow method is tested using electricityegas systems including IEEE14-NGS13 (the IEEE 14-bus system and a 13-node gas network) and IEEE118-NGS31 (the IEEE 118-bus system and a 31-node gas network). The detailed parameters are reported in Ref. [35]. The IEEE14-NGS13 is represented as Case 1, and the 13-node gas network (NGS13) is shown in Fig. 4. In this case: 1) A pressureregulating station is installed at the end of the pipeline connecting gas node 12 to gas node 13; 2) Electrical node 1 is the slack node of the electrical system and is connected to gas node 8; 3) Gas node 1 is the slack node of the natural gas system. The IEEE118-NGS31 is represented as Case 2, and the 31-node gas network (NGS31) is shown in Fig. 5. In this case: 1) Two pressure-regulating stations are installed at the end of two pipelines; one of them connects gas node 9 to gas node 10, and another connects gas node 23 to gas node 24; 2) Electrical node 69 is the slack node of the electrical system and is connected to gas node 31; 3) Gas node 1 is the slack node of the natural gas system. To demonstrate the effectiveness of the proposed probabilistic energy flow method, three models are implemented in the two cases, and the maximum iteration number of the Monte Carlo simulation reaching 10000 is used as the convergence criterion. The details of the compared models are presented in Table 1. 4.2. Verification of the proposed model The value of the proposed model is verified from two aspects: the effectiveness of the thermodynamic processes in the pipelines and compression stations considered herein, and the influence of the pressure-regulating stations. For the first aspect, the calculation results of M1 and M2 are compared. As shown in Table 2, when calculating the energy flow of 10000 sampled states using M1 and M2, we find an interesting phenomenon: only 96.79% and 85.68% of the samples can be calculated successfully by M1 in Case 1 and Case 2, respectively, while the proportions are 100% by M2 in both cases. This phenomenon

7

S2 4

5

3

8

to E3 10

9 2

6

to E1 11

12

1 S1

S

13

Gas Source Compression Station Pressure-regulating Station

Fig. 4. Single line diagram of the test natural gas system (Case 1).

Please cite this article as: Bao S et al., Probabilistic energy flow and risk assessment of ElectricityeGas systems considering the thermodynamic process, Energy, https://doi.org/10.1016/j.energy.2019.116263

8

S. Bao et al. / Energy xxx (xxxx) xxx

23 21 S5 31 to E69

30

28

27

26

20

25

19

29

24 to E49 16

22 18

15

14

13

17 S4

12 9

10

7

6

5

11 4

3

2

1 S1

8 to E12 S3

S2

S

Gas Source Compression Station Pressure-regulating Station

Fig. 5. Single line diagram of the test natural gas system (Case 2).

Table 1 Details of compared models. Thermodynamic Process of Components a

Pipelines Compression Stationsb Pressure-regulating Stationsc

M1 (Proposed in Ref. [15])

M2 (Proposed in this paper)

M3 (Proposed in this paper)

✓  

✓ ✓ 

✓ ✓ ✓

a

In M1, the thermodynamic process in pipelines is considered not properly enough, the detailed explanation is shown in Section I. In M1, the existence of compression stations is considered while their thermodynamic process is ignored. c The existence of pressure-regulating stations is ignored in M1 and M2. Thus, the branches 12e13 in Case 1, 9e10 and 23e24 in Case 2 are pipelines without pressureregulating stations. b

indicates the limitation of application range of M1. The temperature in the pipeline flow model in M1 is represented as Tav ¼ Tsþ(Tm-Tn)/ ln((Tm-Ts)/(Tn-Ts)) [15], in which there is a logarithmic function ln($) whose dependent variable must be larger than zero. Hence, M1 is not applicable when the case (Tm-Ts)/(Tn-Ts)<0 happens. However, the pipeline flow model in M2 (5) is derived by regarding the temperature as a variable in the whole process, thus avoiding the existence of logarithmic function ln($). The reasons for the case (TmTs)/(Tn-Ts)<0 are: 1) The temperature at the inlet end of a pipeline Tm is usually higher than ambient temperature Ts, and because of the heat exchange with surrounding environment, the gas temperature is getting close to Ts as gas flows along a pipeline, which is indicated by the first two terms of (6); 2) after the gas temperature reaches Ts, it will continue to decrease due to the friction during the transmission, which is indicated by the last term of (6). Table 3 shows the gas temperature at each node of Case 1 calculated by M1 and M2, respectively. It is observed that the gas temperatures at different nodes in M1 and M2 vary from 309 K to 285.67 K and from 318.53 K to 287.04 K, respectively. It indicates that with any certain model, the gas temperature at different nodes of a system differs from each other significantly, and the difference can be larger than 30 K. Additionally, the difference between gas temperatures calculated by M1 and

Table 2 Percentage of samples calculated successfully by M1 and M2. Model

Case 1

Case 2

M1 M2

96.79% 100%

85.68% 100%

M2 at node 2 is 9.53 K, which is large. The reason for the gas temperature difference of a node calculated by M1 and M2 is the different consideration of thermodynamic process. In a word, proper consideration of thermodynamic process helps to simulate the distribution of temperature in gas networks authentically. The average differences in the nodal gas temperatures and pressures between M1 and M2 in Case 2 are calculated, as shown in Fig. 6. It is observed that the average differences in gas temperatures are large at nodes adjacent to the compression stations (above 9 K at gas nodes 2, 4, 11, 12, 18, 19, 25, and 26). The reason is that the thermodynamic process in the compression stations is ignored in M1, where the temperatures at inlet and outlet nodes of these components are equal to temperatures of gas sources. This reemphasizes the importance of properly considering the

Table 3 Gas temperature calculated by M1 and M2 in Case 1 (K). Node

M1

M2

1 2 3 4 5 6 7 8 9 10 11 12 13

309.00 309.00 297.58 304.00 304.00 296.90 289.93 287.27 309.00 309.00 297.80 290.39 285.67

309.00 318.53 305.46 303.25 310.57 299.98 291.81 288.10 303.24 314.73 301.41 292.62 287.04

Please cite this article as: Bao S et al., Probabilistic energy flow and risk assessment of ElectricityeGas systems considering the thermodynamic process, Energy, https://doi.org/10.1016/j.energy.2019.116263

S. Bao et al. / Energy xxx (xxxx) xxx

Fig. 6. Average differences in nodal gas pressures and temperatures (Case2).

thermodynamic process to gas temperature calculation. As for gas pressures, the average differences are large at the nodes downstream (the average difference is 14 psi at node 31). This is because of the coupling relationship among gas temperature, pressure and flow in the energy flow model considering the thermodynamic process. Hence, the calculation results of the gas pressures vary with the gas temperatures. The risk indices of some typical nodes calculated by M1 and M2 are shown in Table 4. In Case 1 and Case 2, the PNHTVn of every node is 0. The maximum differences in PHFn are 18% and 27.2% in Case 1 and Case 2, respectively. The maximum differences in PNHPVn are 12.6% and 4.4%, respectively. The maximum differences in PNLPVn are 9.6% and 6.7%, respectively. The values of these indices are obtained based on the calculation results of gas temperatures and pressures. It can be concluded that the chosen model greatly influences the assessed risk level of systems. Besides, when observing node 29 and 31 of Case 2, we find that even the average temperature differences are smaller than 1 K according to Fig. 6, the differences in PHFn reach to 13.7% and 19.6%, respectively. It indicates the high sensibility of hydrate formation to gas temperature and pressure. In conclusion, the difference of results calculated by M1 and M2 are important for two reasons: 1) the computational robustness is guaranteed; 2) risk indices such as PHFn are sensitive to both pressure and temperature. Therefore, the energy flow model should properly consider the thermodynamic process to ensure the accuracy of the calculation results, which is meaningful to the operation management and decision-making of systems. For the second aspect, the influence of the pressure-regulating stations is primarily reflected in the probability of hydrate

Table 4 Risk indices calculated by M1 and M2 (%). Test Cases

Case 1

Case 2

Node

5 8 13 5 22 29 30 31

M1

M2

PHFn

PNHPVn

PNLPVn

PHFn

PNHPVn

PNLPVn

0.0 27.4 21.0 0.0 4.9 33.7 2.8 39.2

46.1 1.1 0.0 13.2 0.0 0.3 0.1 0.0

0.0 0.0 22.1 0.0 0.0 0.1 2.0 9.8

0.0 13.8 3.0 0.0 8.1 47.4 30.0 58.8

33.5 0.6 0.0 8.8 0.0 0.5 0.1 0.0

0.0 0.2 31.7 0.0 0.7 3.3 7.4 16.5

9

formation by M2 and M3. The probabilities of hydrate formation at the outlet nodes of the pressure-regulating stations are shown in Fig. 7, with the fixed outlet pressures varying from 590 psi to 650 psi in Case 1, and from 630 psi to 690 psi in Case 2. As shown, the PHFn of M3 (with outlet pressure fixed to 590 psi) is 28.49% higher than that of M2 at node 13 in Case 1; the PHFn of M3 (with outlet pressure fixed to 690 psi) is 12.47% higher than that of M2 at node 24 in Case 2. These results indicate that when installing a pressure-regulating station at the end of a pipeline, the PHFn at the outlet node becomes higher. This is because that the outlet pressure and temperature decrease simultaneously owing to the thermodynamic process in a pressureregulating station. Therefore, for the safe operation of the systems, we should emphasize the hydrate formation at the pressureregulating stations. Additionally, when the fixed values of the outlet pressures vary, the PHFn at different nodes in the same system could change in different trend, which can be observed from (b) and (c) of Fig. 7. When the outlet pressures are set from 630 psi to 690 psi, the PHFn of M3 decreases from 16.73% to 15.42% at node 10 in Case 2; however, it increases from 12.07% to 17.5% at node 24 in Case 2. It is because that the outlet temperature changes in the same trend with the outlet pressure according to (7), and the formation of hydrates is sensitive to both pressure and temperature. When the pressures are fixed at the same value, the slight difference in temperatures can lead to different results. 4.3. Influence of electrical systems on natural gas systems The systems in this study is operated according to the “electricity-orientated mode”. The influence of electrical systems with different wind penetrations and wind characteristics (depicted by the expectation value m and standard derivation d of wind speeds, which are calculated by the shape parameter k and scale parameter c of a Weibull distribution) on natural gas systems is tested by calculating the proposed risk indices. The values of PHFn, PNHPVn, and PNLPVn are calculated when the wind penetration rate varies from 0% to 35% with different m and d of wind speeds, as shown in Figs. 8 and 9. Obviously, as the wind penetration increases, the PHFn and PNHPVn increase while the PNLPVn decreases. The reason is that a higher wind penetration rate generates more power at the wind farms, thus decreasing the gas flows supplied to gas-fired generators. Consequently, the nodal gas pressures become higher owing to the coupling relationship among the state variables. Meanwhile, when the wind penetration rate is 35% and the d of the wind speeds is 1.73, the PHFn of node 31 in Case 2 increases from 37.58% to 46.61% as the m varies from 3.92 to 5.92. For the PHFn of the same node with the same wind penetration rate, it increases from 37.58% to 40.98% when the m is 3.92, and the d varies from 1.73 to 3.73. This indicates that the probabilistic characteristics of wind speeds also influence the values of the proposed risk indices. Consequently, the risks of natural gas systems are influenced by the power generated at wind farms, thereby necessitating the integrated analysis of electrical and natural gas systems. 5. Conclusions Herein, a probability energy flow analysis method of electricityegas systems considering the thermodynamic process in gas networks is proposed. The effectiveness of this method is verified. The following conclusions are drawn: 1) For the energy flow model, a pipeline flow model under the non-isothermal condition is derived, and the thermodynamic

Please cite this article as: Bao S et al., Probabilistic energy flow and risk assessment of ElectricityeGas systems considering the thermodynamic process, Energy, https://doi.org/10.1016/j.energy.2019.116263

10

S. Bao et al. / Energy xxx (xxxx) xxx

(a)

(b)

(c)

Fig. 7. PHFn calculated by M2 and M3. (a) PHFn of node 13 in Case 1. (b) PHFn of node 10 in Case 2. (c) PHFn of node 24 in Case 2.

(a)

(b)

Fig. 8. PHFn calculated by M3. (a) PHFn of node 8 in Case 1. (b) PHFn of node 31 in Case 2.

(a)

(b)

Fig. 9. PNHPVn/PNLPVn calculated by M3. (a) PNHPVn of node 5 in Case 1. (b) PNLPVn of node 28 in Case 2.

processes in the pressure-regulating stations and compression stations is considered in the electricityegas systems for the first time. Compared with existing gas flow models used in electricity-gas system analysis, two advantages of the proposed model are demonstrated by the case studies. First, the proposed model performs better in computational robustness as it can successfully obtain the calculation results of all samples while existing models cannot be applied in some cases. Additionally, it

considers the thermodynamic processes in different gas components more properly than existing models. 2) A series of risk indices including the probabilities of hydrate formation and state variables exceeding the limits are proposed to quantify the risks in the natural gas systems. Compared with previous work, the uncertainties of the unified system are considered, the risk of hydrate formation is illustrated as probabilities while existing methods can only verify the hydrate

Please cite this article as: Bao S et al., Probabilistic energy flow and risk assessment of ElectricityeGas systems considering the thermodynamic process, Energy, https://doi.org/10.1016/j.energy.2019.116263

S. Bao et al. / Energy xxx (xxxx) xxx

formation of each node at a deterministic state. These indices are helpful to illustrate the effectiveness of the proposed models and can provide guidance for the operation and decisionmaking of electricityegas systems. 3) The cases studied demonstrate that the wind penetration rate and probabilistic characteristics of wind speeds affect the risks of natural gas systems significantly. The tight coupling relationship between electrical and natural gas systems is therefore indicated. Acknowledgments This work was supported in part by the Science and Technology Project of State Grid Corporation of China under grant 5100201999333A-0-0-00, part by the Fundamental Research Funds for the Central Universities under grant No. 2018CDQYDQ0006, and part by the Academician-mentoring Scientific and Technological Innovation Project of Chongqing under grant No. cstc2018jcyjyszxX0001.

[14]

[15]

[16]

[17]

[18]

[19]

[20]

[21]

References

[22]

[1] Dokic SB, Rajakovic NL. Security modelling of integrated gas and electrical power systems by analyzing critical situations and potentials for performance optimization. Energy 2018. available online, https://doi.org/10.1016/j.energy. 2018.04.165. [2] Sheikhi A, Bahrami S, Ranjbar AM. An autonomous demand response program for electricity and natural gas networks in smart energy hubs. Energy 2015;89:490e9. [3] Trifonov TO. Coordination of battery energy storage and power-to-gas in distribution systems. Prot Control Mod Power Syst 2017;2(2):421e8. [4] Badal FR, Das P, Sarker SK, Das SK. A survey on control issues in renewable energy integration and microgrid. Prot Control Mod Power Syst 2019;4(4): 87e113. [5] Thermodynamic Process. [Online] Available: https://en.wikipedia.org/wiki/ Thermodynamic_process#Flow_process. [6] Ghavipour M, Ghavipour M, Chitsazan M, Najibi SH, Ghidary SS. Experimental study of natural gas hydrates and a novel use of neural network to predict hydrate formation conditions. Chem Eng Res Des 2013;91(2):264e73. [7] Nasrifar K, Moshfeghian M, Maddox RN. Prediction of equilibrium conditions for gas hydrate formation in the mixtures of both electrolytes and alcohol. Fluid Phase Equilib 1998;146(1e2):1e13. [8] Sloan ED. CA Koh. Clathrate hydrates of natural gases. Boca Raton, FL, USA: CRC Press; 2007. Mercel Dekker Inc N Y, 1990. [9] Ashouri E, Veysi F, Shojaeizadeh E, Asadi M. The minimum gas temperature at the inlet of regulators in natural gas pressure reduction stations (CGS) for energy saving in water bath heaters. J Nat Gas Sci Eng 2014;21:230e40. [10] Yu J, Lin W, Kamel S, Li WY, Zhang L. An equivalent model considering frequency characteristics and renewable uncertainties for probabilistic power flow. IET Gener, Transm Distrib 2018;12(22):5939e48. [11] Lin W, Yang ZF, Yu J, Bao SY, Dai W. Toward fast calculation of probabilistic optimal power flow. IEEE Trans Power Syst 2019;34(4):3286e8. [12] Hu Y, Lian H, Bie Z, Zhou B. Unified probabilistic gas and power flow. J. Mod. Power Syst. Clean Energy 2017;5(3):400e11. [13] Qi S, Wang X, Wang Y, Tian S, Wang R. Integrated probabilistic energy flow analysis in natural gas and electricity coupled systems considering the

[23] [24] [25]

[26] [27] [28] [29]

[30]

[31] [32]

[33] [34]

[35]

11

randomness of wind power. In: IEEE conference on energy internet and energy system integration; 2017. p. 1e6. Nov. 26-28 2017, Beijing, China. Chen S, Wei Z, Sun G, Cheung KW, Sun Y. Multi-linear probabilistic energy flow analysis of integrated electrical and natural-gas systems. IEEE Trans Power Syst 2017;32(3):1970e9. Martinez-Mares A, Fuerte-Esquivel CR. Integrated energy flow analysis in natural gas and electricity coupled systems. In: IEEE north American power symposium; 2011. p. 1e7. Aug. 4-6 2011, Boston, USA. Martinez-Mares A, Fuerte-Esquivel CR. A unified gas and power flow analysis in natural gas and electricity coupled networks. IEEE Trans Power Syst 2012;27(4):2156e66. Farzaneh-Gord M, Rahbari HR, Bajelan M, Pilehvari L. Investigation of hydrate formation in natural gas flow through underground transmission pipeline. J Nat Gas Sci Eng 2013;15:27e37. pez-Benito A, Elorza-Tenreiro FJ, Gutie rrez-Pe rez LC. Steady-state nonLo isothermal flow model for natural gas transmission in pipes. Appl Math Model 2016;40(23e24):10020e37. Biskas PN, Kanelakis NG, Papamatthaiou A, Alexandridis I. Coupled optimization of electricity and natural gas systems using augmented Lagrangian and an alternating minimization method. Int J Electr Power Energy Syst 2016;80: 202e18. Pambour KA, Erdener BC, Bolado-Lavin R, Dijkema GPJ. SAInteA novel quasidynamic model for assessing security of supply in coupled gas and electricity transmission networks. Appl Energy 2017;203:829e57. Mokhatab S, Poe WA, Speight JG. Handbook of natural gas transmission and processing. Burlington, VT, USA: Gulf Professional Publishing; 2006. Fan M, Vittal V, Heydt GT, Ayyanar R. Probabilistic power flow analysis with generation dispatch including photovoltaic resources. IEEE Trans Power Syst 2013;28(2):1797e805. Bazvand M, Farhangian H, Najaflou P. Comparison of several methods for obtaining hydrate formation pressure. Pet Sci Technol 2018;36(11):801e6. Naseer M, Atter WB. hydrate formation in natural gas pipelines. WIT Trans Eng Sci 2011;70:261e70. Liu WY, Hu JQ, Li XF, Sun Z, Sun FR, Chu HY. Assessment of hydrate blockage risk in long-distance natural gas transmission pipelines. J Nat Gas Sci Eng 2018;60:256e70. Chaczykowski M, Zarodkiewiczb P. Simulation of natural gas quality distribution for pipeline systems. Energy 2017;134:681e98. Coutler TL, Bardon MF. Revised equation improves flowing gas temperature prediction. Oil Gas J 1979;26:107e8. Calculation of Pressure Regulators for Gas. [Online] Available: https://www. mankenberg.de/en/173/calculation-of-pressure-regulators-for-gas.html. Li Q, An S, Gedra TW. Solving natural gas load flow problems using electric load flow techniques. In: 2003 proc of the north American power symposium. USA: Rolla; 2003. p. 1e7. Seguro JV, Lambert TW. Modern estimation of the parameters of the Weibull wind speed distribution for wind energy analysis. J Wind Eng Ind Aerodyn 2000;85(1):75e84.  SA, Dinler A. A new method to estimate Weibull parameters for wind Akdag energy applications. Energy Convers Manag 2009;50(7):1761e6. Zhang J, Xiong G, Meng K, Yu P, Yao G, Dong Z. An improved probabilistic load flow simulation method considering correlated stochastic variable. Int J Electr Power Energy Syst 2019;111:260e8. SS Byk VI. Fomina gas hydrates. Russ Chem Rev 1968;37(6):469e91. Yu J, Guo L, Ma MN, Kamel S, Li WY, Song XM. Risk assessment of integrated electrical, natural gas and district heating systems considering solar thermal CHP plants and electric boilers. Int J Electr Power Energy Syst 2018;103: 277e87. Data of the IEEE 14-NGS 13 System and the IEEE 118- NGS 31 System. [Online] Available: https://figshare.com/articles/Gas_Model_Considering_ Thermodynamic_Process_Theory_and_Application_in_Probabilistic_Energy_ Flow_Analysis/7478537.

Please cite this article as: Bao S et al., Probabilistic energy flow and risk assessment of ElectricityeGas systems considering the thermodynamic process, Energy, https://doi.org/10.1016/j.energy.2019.116263