Nodal-pressure-based heating flow model for analyzing heating networks in integrated energy systems

Nodal-pressure-based heating flow model for analyzing heating networks in integrated energy systems

Energy Conversion and Management 206 (2020) 112491 Contents lists available at ScienceDirect Energy Conversion and Management journal homepage: www...

2MB Sizes 0 Downloads 27 Views

Energy Conversion and Management 206 (2020) 112491

Contents lists available at ScienceDirect

Energy Conversion and Management journal homepage: www.elsevier.com/locate/enconman

Nodal-pressure-based heating flow model for analyzing heating networks in integrated energy systems

T



Dongwen Chena, Xiao Hub, Yong Lia, , Ruzhu Wanga, Zulkarnain Abbasa, Shunqi Zengc, Li Wangc a

Institute of Refrigeration and Cryogenics, Shanghai Jiao Tong University, Shanghai 200240, China School of Electronic Information and Electrical Engineering, Shanghai Jiao Tong University, Shanghai 200240, China c System Operation Department, Guangzhou Power Supply Group, Guangzhou 510060, China b

A R T I C LE I N FO

A B S T R A C T

Keywords: Nodal-pressure-based model District heating networks analysis and planning Heating flow Unified nodal pressure model Energy transmission loss Integrated energy system

District heating technology has prompted the development of integrated energy systems, where energy flow models have been widely used. However, as heating networks become complicated, the traditional mass flow model fails to provide efficient solutions. This study proposes a novel heating flow model based on nodal pressure in heating networks. First, transmission models based on nodal pressure are formulated. Subsequently, nodal conservation models of mass and heating flow are developed. Further, the transmission and nodal conservation models are combined, and thereby nodal pressure and user load are directly correlated. As the number of variables depends only on the nodes, the resulting model is more efficient at adapting to new network connections when the nodes do not change. Moreover, as loop pressure conservation equations are omitted in this model, fewer iterations are required. Owing to the small relative fluctuations of nodal pressure, initial iteration values can be set more effectively in the solution of the proposed model. A case study was conducted, and it was demonstrated that compared with the traditional mass flow model, the proposed model is more effective in adapting to varying network topologies, and the computational burden is significantly reduced. Moreover, the proposed model can be efficiently applied for analyzing energy transmission loss as well as system planning in integrated power–gas–heat energy systems.

1. Introduction Integrated energy systems have been extensively studied in recent years, as they can significantly improve the stability and efficiency of energy supply and consumption [1]. However, compared with traditional electricity networks, integrated energy systems contain a variety of energy networks, such as electricity, natural gas, and heat, so that the energy transmission analysis is more challenging to carry out. Accordingly, the development of efficient integrated energy network modeling schemes has become a critical issue [2,3]. With the evolution of integrated energy systems, the complexity of energy networks has increased significantly. Therefore it is urgent to develop efficient energy flow models that can be applied not only for network parameter analysis but also for network planning. The concept of energy flow, which describes energy transmission in networks, has been proposed to analyze and model energy networks [4]. Two critical related aspects are energy flow formulation and calculation. For the formulation of energy flow in gas networks, nodal pressure models have been developed, including flow resistance and



nodal mass flow conservation models [5]. The flow resistance equations in each pipeline and the nodal mass flow conservation equations are both described by nodal pressure differences. Generally, in nodal pressure models, the temperature is considered to yield more accurate results [6]. However, calculations have demonstrated that the influence of temperature can be neglected in gas distribution networks. Hence, traditional nodal pressure gas flow models have neglected the influence of temperature, and they are formulated based only on nodal pressure [7]. Based on this model, the influence of different components of injected gas [8] and the influences of the uncertainty of renewable energy/user [9] have been studied. Furthermore, the nodal pressure model has been used to schedule different integration levels of nodal electricity generator and gas compressor [10], as well as optimize scheduling with a new transactive approach [11]. Although nodal pressure models have been widely applied for gas network analysis, gas flow models consider the only mass transmission. The heating flow is the simultaneous transmission of mass and thermal energy; hence, nodal pressure gas flow models cannot be directly applied to district heating networks. In previous studies on heating

Corresponding author. E-mail address: [email protected] (Y. Li).

https://doi.org/10.1016/j.enconman.2020.112491 Received 6 November 2019; Received in revised form 7 January 2020; Accepted 8 January 2020 0196-8904/ © 2020 Elsevier Ltd. All rights reserved.

Energy Conversion and Management 206 (2020) 112491

D. Chen, et al.

Nomenclature

fi g Gij Ii Iij K lij m mij mh,ij

A. Greek letters ηcomp efficiency of compressor or pump ηsteam,H GB energy conversion efficiency from steam to heat by gas boiler ηgas,E GT energy conversion efficiency from natural gas to electricity by gas turbine ηsteam,H HB energy conversion efficiency from steam to heat by heat recovery boiler ηE,CH HP energy conversion efficiency from electricity to cooling by heat pump ηsun,E PV energy conversion efficiency from solar to electricity by photovoltaic panel λ coefficient of on-way friction γ calorific value Φg,i user’s power load of natural gas at node i ρ density ρE,ij resistivity

ph,j pi p0 Δpij Δpf,ij Δpg,ij Δph,ij ΔPE,ij ΔPP,ij ΔPH,ij Qi Re Si sij Tam T0 Ti

B. Parameters c Cp d Dg e ei

flow velocity specific heat capacity diameter diameter of pipeline base of natural logarithms real part of voltage at node i

networks, mass flow models, which are based on the mass flow in pipelines, have been developed [12–17]. These mass flow based heating flow models include models of flowing resistance and temperature transmission resistance, nodal conservation models of mass flow and heating flow (in which heating flow is described by mass flow), and models of loop pressure drop conservation. In certain simplified heating flow models based on mass flow, temperature is neglected [12,13].These models also neglect temperature transmission resistance and nodal heating flow conservation, and they are suitable for shortdistance modeling of the heating flow. A temperature transmission resistance model for heating networks was developed [14]. Then a heating flow model, considering temperature drop in energy networks, was proposed for radial networks [15]. By considering loop pressure equation, the heating flow model was improved for analyzing electricity-heating networks [16] and electricity-gas-heating networks [17]. The model can be applied for both loop and radial networks. In these mass flow-based heating flow model, the temperature transmission resistance, and the nodal conservation of the mass and heating flow were both described by the mass flow. In heating flow models based on mass flow, the number of independent variables describing the conservation of mass flow and loop pressure drop is equal to the number of branches. If the number of network loops increases but the nodes do not vary, the complexity of the formulation will increase accordingly. Therefore, it is desired to develop a heating model based on nodal parameters, where model complexity does not depend on the network loops or topology when the number of nodes is constant. The solution of an energy flow model is generally based on the Newton–Raphson algorithm [18]. The convergence rate of this algorithm, which refers to the number of iterations, primarily depends on the Jacobian matrix [19,20], whose size is determined by the number of independent variables [21,22]. As there are a large number of independent variables in traditional mass flow models, the convergence rate is reduced significantly. To resolve this, a variable reduction is required.

imaginary part of voltage at node i gravity acceleration real part of element in admittance matrix current at node i current in line between node i and node j surface roughness of pipeline length between node i and node j mass flow mass flow in pipeline between node i and node j normalized hot water mass flow in pipeline between node i and node j normalized pressure at node i pressure at node i standard atmospheric pressure total pressure drop pressure drop of on-way resistance normalized pressure drop in natural gas pipeline normalized pressure drop in heating pipeline active power loss pressure loss thermal loss reactive power at node i Reynolds number apparent power at node i element of correlation coefficient ambient temperature temperature after use temperature at node i

Moreover, the relative difference between the initial iteration values and the resulting data usually has a critical influence on the numerical solution [23]. In mass flow models, the initial iteration values in each pipeline are set to the maximum possible mass flow. However, the mass flow rates in various pipelines are quite different in district heating networks [16,17]. This will lead to a significant reduction of the convergence rate, and in some cases, convergence cannot even be achieved. Hence, mass flow models have apparent weaknesses in terms of their numerical solution. Energy flow models have been widely adopted in studying integrated energy system optimal dispatching [24–27]. Chen et al. [25] suggested a method with limited numbers of control actions to simplify the optimization of energy flow into a mixed-integer linear programming problem. Based on energy flow models, an optimally coordinated dispatch of heat source/fluid [26] and a coordinated sequential optimal method for energy storage equipment and conversion devices were proposed [27]. Energy flow models, based on nodal voltage and pressure, have also been effectively applied in network planning with different focuses [28–32]. By using the energy flow model, the combined electricity-gas network expansion planning has been conducted [28], and the minimized cost of meeting gas and electricity demand has been achieved [29]. Ma et al. [30] conducted generic optimal planning to obtain both the optimal structure configuration and energy management strategy based on energy flow. By applying the energy flow model and optimization in the process of maximizing the investment benefits and minimize line losses [31], a two-stage method to optimize the capacity of distributed generation and energy storage was proposed. Wang et al. [32] proposed an expansion planning model of the multi-energy system with the integration of an active electricity distribution network, natural gas network, and energy hub. However, to the author’s best knowledge, the combined electricity-gas-heating network planning method is still not available due to the lacking heating flow model based on nodal parameters. 2

3

Network panning

√ × × × √

Optimal dispatching

√ √ √ √ √ N-N B N-N-B N-N-B N-N-N

Independent variables*

Electricity-gas networks Heating networks Electricity-gas-heating networks

In this study, models of five widely used energy conversion devices, namely, turbine CHP units, gas boilers (GBs), electrical boilers (EBs), water pumps for circulating (WP), and gas compressors (GCs), are developed. In gas turbine CHP units, natural gas is burned to generate power, and then the wasted high-temperature tail gas is captured by a heat recovery boiler to generate heat. In GBs, natural gas is combusted for heating. The GB is typically accompanied by CHP units to improve the flexibility of dispatching. In EBs, electricity is converted into heat. The EB is typically used as the supplement of GB. WPs and GCs are used for raising the pressure in heating and gas networks, respectively. Both

[6–11] [12,13] [15] [16,17] Our model

2. Multi-energy conversion device modeling

Networks

The paper has been organized as follows. Section 2 presents models of typical energy conversion equipment used in integrated electricity–gas–heating networks. Section 3 presents the detailed formulation of nodal pressure models for district heating networks. Section 4 presents energy transmission loss models, which are essential in integrated electricity–gas–heating networks. Section 5 presents a case study in which the proposed model is compared with the traditional mass flow model. Then, the energy transmission loss in electricity–gas–heating networks is analyzed. Section 6 concludes the paper.

Refs

Table 1 A review of previous energy flow models for integrated energy system.

Applications

Contributions

(1) A new heating flow model, based on nodal pressure, is proposed for district heating networks to improve the performance of the traditional mass flow based heating flow model. In the model, the mass and thermal energy flows are described by the nodal pressure, and the loop pressure conservation is integrated into the nodal mass flow conservation. Thereby, nodal pressure in heating and gas networks can be analogized to nodal voltage in electricity networks. (2) Compared with traditional mass flow models, the proposed model allows a more efficient formulation with the same correctness and numerical results within the allowable error range. In the proposed model, the nodal pressure has been used to define the independent variables for the simplified formulation, and improved adaptability to network topology changes for network optimization or planning. (3) The proposed model is efficient in terms of solution convergence and iteration initialization. It has fewer variables and equations, as well as a smaller Jacobian matrix. Therefore, a higher convergence rate can be achieved. The proposed model can obtain the numerical results within the allowable error range when compared with the existing model. Moreover, the relative fluctuation of the nodal pressure at different nodes is smaller than the relative fluctuation of the mass flow in different pipelines. Thus, the initial iteration values can be set more effectively. (4) To analyze the energy transmission loss, the pressure drop energy loss in district heating and gas networks has been calculated along with the active power and heating energy loss in integrated electricity–gas–heat energy systems.

Nodal voltage-pressure model Mass flow-based heating network analysis model Mass flow-based heating flow model for radial networks considering both temperature and pressure drop Mass flow-based heating flow model for generic networks including loop and radial networks Nodal pressure-based heating flow model for generic networks including loop and radial networks

A review of previous works about the energy flow modeling in integrated energy systems is presented in Table 1. At the end of the table, the proposed model has been compared with the previous models. To perform network optimization or planning, the values of the objective function of different branch connection schemes should be compared when the nodes are constant. Hence, to choose a planning model, the condition of constant nodal parameters should be fully met. Therefore, energy flow models based on nodal parameters are suitable for network planning [28,33]. However, new variables or pressure drop conservation equations are used in energy flow models based on branch parameters. Traditional mass flow models are therefore ineffective, complex, and challenging to apply in the network optimization or planning. To overcome these issues, this study proposes an efficient heating flow model based on nodal pressure. As shown in Table 1, compared with the previous work, the main contributions of this study are the following:

* N refers to the model formulated based on nodal voltage/pressure; B refers to the model formulated based on branch mass flow.

Energy Conversion and Management 206 (2020) 112491

D. Chen, et al.

Energy Conversion and Management 206 (2020) 112491

D. Chen, et al.

3.1.1. Transmission equations of mass flow and temperature drop

WPs and GCs can convert electrical power to pressure potential energy. The two devices are used to increase the pressure in the pipeline. In CHP units, the energy conversion efficiency from natural gas to electricity is defined as follows: G, E ηCHP =

PE , v·γ

(1) Mass flow described by nodal pressure In a water pipeline, flow resistance includes on-way resistance due to friction and local resistance due to geometric changes in the pipeline. The on-way resistance is calculated by [17,34]:

(1)

where PE is the electrical outlet power, v is the volume flow, and γ is the calorific value of natural gas. Then, the energy conversion efficiency from natural gas to heat is defined as follows: G, H G, H ηCHP = (1 − ηCHP )·ηH , R ,

Δpf , ij =

Re =

P = H, v·γ

(3)

PH . PE

(4)

In WPs, the energy conversion efficiency from electricity to pressure is defined as follows:

ηWP =

(pout − pin )·v·ρ , PE

(5)

⎡ kG ·p v · kG − 1 in in ⎢ ηGC =

( ) ⎣ pout pin

PE

⎤ − 1⎥ ⎦,

(7)

4m , πdρμ

Δpf , ij =

where pout and pin are the WP outlet and inlet pressure, respectively, and ρ is the density of water. In GCs, the energy efficiency from electricity to pressure is defined as follows: kG − 1 kG

λ,

(8)

where μ is the kinematic viscosity of water, which is a function of temperature. Five flow states are summarized in Table 2 [34]: where K is the surface roughness of the pipeline. The surface roughness values for different pipeline materials are listed in Appendix A. Considering that steel is a widely used pipeline material, the roughness is 0.5 mm. For heating networks, considering that the maximum flow velocity would not exceed 2.5 m/s, the Reynolds number is within the range of the fourth flow state listed in Table 2. Thus, the onway resistance is calculated as follows:

where PH is the heating outlet power. In EBs, the energy conversion efficiency from electricity to heat is defined as follows: E,H ηEB =

π 2d5ρ

where m is the mass flow, lij is the length between nodes i and j, d is the diameter of the pipeline, ρ is the density of water, and λ is the coefficient of on-way friction (which is a function of mass flow). To determine the relationship between λ and m, the Reynolds number is introduced as follows:

(2)

where ηH,R is the heat recovery efficiency of the CHP unit. In GBs, the energy conversion efficiency from natural gas to electricity is defined as follows: G, H ηGB

8lij m2

17πρμd ⎞0.25 0.88lm2 ⎛ K 0.88lm2 K 0.25 ⎛ ⎞ . + ≈ π 2d5ρ ⎝ d m ⎠ π 2d5ρ ⎝ d ⎠

(9)

The mass flow for heating networks is expressed as:

mij =

1/8 1.066πd5/2ρ1/2 Δpf1/2 , ij ⎛d⎞ , 1/2 l ⎝K ⎠

(10)

where mij is the mass flow of the heating pipeline between nodes i and j. The relationship between the mass flow and the pressure drop between these two nodes is determined. In this model, the mass flow equation is described by nodal pressure by Eq. (10).

(6)

where kG is the adiabatic index of natural gas, which can be considered constant (1.30), and vin is the inlet volume flow.

(2) Temperature drop 3. Energy flow modeling for electricity-gas-heating distribution networks

In radial networks, when hot water flows from node i to node j, the relationship between the temperatures at these nodes is determined as follows [16,35]:

3.1. Heating flow model based on nodal pressure

Tj = Tam + (Ti − Tam) e

In this section, a nodal pressure model for district heating networks is developed in detail. The model is summarized into five steps, as shown in Fig. 1. The model for heating fluid networks is developed according to Fig. 1. The process comprises five steps: 1) Formulating the equations of flowing and temperature transmission resistance, in which mass flow and temperature drop in pipelines are both described by nodal pressure differences between the inlet and outlet. 2) Formulating the nodal conservation equations of mass and heating flow. Heating flow is described by mass flow and nodal temperature. 3) For convenience, the flow transmission equations are normalized. In the model, only coefficients and the nodal pressure are present. The pressure at node 1, which is used as the slack node, is a normalized parameter that equals 1.0. 4) By combining the normalized flow transmission and the nodal conservation models, the normalized pressure model for district heating networks is constructed. 5) The Newton algorithm is applied for numerically solving the model. The iteration schemes are applied in this part.



λ 0 lij Cp mij ,

Fig. 1. Nodal pressure model for district heating networks. 4

(11)

Energy Conversion and Management 206 (2020) 112491

D. Chen, et al.

Table 2 On-way resistance for five flow states.

Gj λ 0 lij 1/8 1.066πd5/2ρ1/2 Δpf1/2 − , ij ⎛ d ⎞ ·Cp·⎡Tam + (Ti − Tam) e Cp mij ⎢ 1/2 l ⎝K ⎠ ⎣

n

Re

λ

Δpf,ij

Re < 2300

16π d2ρ m

̂ 128 I ¼ml π d3

2300 < Re < 4000

0.025Re1/3

=

∑ i = 1, ≠ j

⎤ − T0⎥ − Φj = 0, ⎦

0.03175lm7/3 ̂ 1/3 π7/3d16/3ρ4/3 I ¼

4000 ≤ Re < 22.2(d/K)

8/7

0.03164 Re1/4

22.2(d/K)8/7 ≤ Re < 597 (d/K)9/8

0.11( +

Re ≥ 597(d/K)9/8

0.11( )

K d

68 1/4 ) Re

1/4 ̂ 0.88lm2 K 17πρ I ¼d ( + ) m π 2d5ρ d

0.88lm2 K 1/4 ( ) π 2d5ρ d

K 1/4 d

3.1.3. Normalization for pressure-related models Since the heating energy flow conservation equation is formulated based mass flow conservation, the normalization is only carried out for mass flow conservation equation. Assume that the following variables are defined as

where Tam is the ambient temperature, which is assumed constant throughout the year; λ0 is the thermal conductivity per unit length of the pipeline, which depends mainly on the thickness and type of thermal insulation material. It is worth noting that if hot water flows from j to node i, then mij < 0, which is consistent with Eq. (11). By substituting Eq. (10) into Eq. (11), the temperature drop is described by nodal pressure. Therefore, the transmission equations for both mass flow and temperature drop are described by nodal pressures.

(1) Mass flow equation at node

Cp (Tj − T0)

ph, i = pi / p1 kh, ij = sij

∑ j ∈ i, ≠ i

⎛d⎞ ⎝K ⎠



Φj Cp (Tj − T0 )

(15)

ph, i − ph, j |ph, i − ph, j |

, (16)

3.1.4. Nodal pressure model for district heating networks The nodal pressure model for node i is shown as: 1



1.066π

D 5ρ D 8 ⎛ ⎞ · l i, j ⎝ K ⎠

pi − pj

where Tj is the temperature at node j. T0 is the hot water temperature after use. Cp is the specific heat of water at constant pressure. Φj is the heating load of thermal energy. By combining the mass flow from other nodes in Eq. (10) and mass flow converted from user heating load in Eq. (12), the nodal mass flow conservation is formulated as follows:

Fi =

,

0 i, j ⎧ ⎫ ⎪ ⎪ ∑ mi, j [Tam + (Ti − Tam)·e mi, j Cp ] / ∑ mj, i − T0 ·m1· ph, i ⎨ j ∈ i, ≠ i ⎬ j ∈ i, ≠ i ⎪ ⎪ ⎩ ⎭

mh, ij = kh, ij

(12)

1/8

d 1/8 K

()

−λ l

Cp

|pi − pj |

1.066πd5/2ρ1/2 Δpf1/2 , ij l1/2

lij1/2 m1 Φj

i ∈ j, ≠ j

,

1.066πd5/2ρ1/2

Then Eq. (10) is simplified as follows:

Similar to Kirchhoff’s current law, the nodal mass flow conservation states that the algebraic sum of mass flow at any node equals zero for the steady flowing network. For heating networks, the mass flow includes two parts: mass flows from other nodes and user load. The mass flow from other nodes in each pipeline can be calculated by Eq. (10) based on nodal pressures. The user load of mass flow is calculated by converting the heating load into mass flow load. The converted mass flow is calculated by the following:

Φj

mh, ij = mij / m1

⎧ ⎪ ⎪ ⎪ ⎪ ⎨ ⎪ kh, ii = ⎪ ⎪ ⎪ ⎩

3.1.2. Mass flow and heating flow conservation at node

mj,0 =

(14)

where Φj is the heating requirement at node j. The heating energy supplied to the users originates from the temperature difference (i.e. Tj − T0). In this model, the nodal heating flow conservation model is described by nodal pressures.

̂ 1/4lm7/4 1.7898 I ¼ π7/4 d19/4ρ3/4

Φj

+

−λ 0 li . j

Cp ⎧ ∑ mi, j [Tam + (Ti − Tam)·e mi, j Cp ]/ ∑ mi, j − T0 ⎫ ⎬ ⎨ i ∈ j, ≠ j i ∈ j, ≠ j ⎭ ⎩ (17)

= 0,

The normalized nodal pressure model for Eq. (13) is shown as following: n

Fi =

∑ j = 1, ≠ i

= 0, (13)

kh, ij·(ph, i − ph, j ) |ph, i − ph, j |

+ kh, ii Δph, ii = 0, (18)

From node 1 to node n, the nodal pressure model for entire network is developed:

where ΔPf,ji is the pressure difference between nodes j and i. mjj is injected mass flow. In this model, the nodal mass flow conservation model is described by nodal pressures.

k

(p

−p

)

k

(p

−p

)

k

(p

−p

)

⎧ kh,11 Δph,11 + h,12 h,1 h,2 + h,13 h,1 h,3 + ...+ h,1n h,1 h, n |ph,1 − ph,2 | |ph,1 − ph,3 | |ph,1 − ph, n | ⎪ ⎪ kh,21 (ph,2 − ph,1) kh,23 (ph,2 − ph,3 ) kh,2n (ph,2 − ph, n ) ⎪ |p − p | + kh,22 Δph,22 + |p − p | + ...+ |p − p | h,2 h,1 h,2 h,3 h,2 h, n ⎪ ... kh, in (ph, i − ph, n ) ⎨ kh, i1 (ph, i − ph,1) + kh, i2 (ph, i − ph,2 ) + ...+k h, ii Δph, ii + ... + ⎪ |ph, i − ph,1 | |ph, i − ph,2 | |ph, i − ph, n | ⎪ ... ⎪ kh, n1 (ph, n − ph,1) kh, n2 (ph, n − ph,2 ) kh, n3 (ph, n − ph,3 ) + + ...+kh, nn Δph, nn ⎪ |p − p | + |ph, n − ph,2 | |ph, n − ph,3 | h, n h,1 ⎩

(2) Heating energy flow equation at node At each node, if the flow work for a node can be neglected, the heating flow conservation can be obtained. In the heating network, thermal energy is the primary user load, which is hundreds of times greater than flow work. Therefore, the nodal heating flow conservation is feasible. The inlet thermal power minus outlet thermal power equals zero. The heating energy flow equation at node j is formulated as follows:

=0 =0 , =0 =0 (19)

3.1.5. Application of Newton algorithm In the energy flow analysis, the temperature and pressure at every node, except node 1 (the station site), are unknown variables. Typically, 5

Energy Conversion and Management 206 (2020) 112491

D. Chen, et al.

the Newton algorithm is applied for solving energy flow equations [16,17]. By combining Eq. (14) and Eq. (18), the entire nodal pressure based heating flow model at any node i is expressed as following:

method for natural gas distribution networks under low pressure (no more than 75 kPa) are developed. The pressures at nodes i and j are related as follows [17]:

Fi = 0 , node i: ⎧ ⎨ ⎩Gi = 0

pi − pj = 1.22 × 10−2

(20)

The nodal pressure based heating flow model for entire network is expressed as following:

⎧ F1 = 0 ⎪ F2 = 0 ⎪ ⋯ ⎪ Fn = 0 ⎡ F = 0⎤ = , ⎣ G = 0 ⎦ ⎨ G1 = 0 ⎪G2 = 0 ⎪ ⋯ ⎪ ⎩Gn = 0

vg, ij − Li = 0 (28)

j = 1, ≠ i

where Li is the outflow natural gas. For practical applications, the volumetric flow must be converted from actual to standard conditions (i.e. pressure = 101.325 kPa and temperature = 0 °C) using this formula [8,36]:

(21)

Li =

Φgi

p0 273 + Ti γ pi + p0 273



∂G ⎟ ⎟ ∂T ⎠

(29)

where, γ is the calorific value of combustion, Φgi is the user energy load, Ti is the temperature at node i, Zi and Z0 are the gas compressibility factors under the actual conditions at node i and the standard conditions, respectively. The nodal pressure model for node i is formulated as:

(22)

When applying Newton algorithm, the Jacobian matrix is determined as follows: ∂F ∂T

,



(23)

j ∈ i, ≠ i

p1 Dg5

pi − pj

1.22 × 10−2lij v12

|pi − pj |

4. Energy transmission losses

⎛ Δp ⎞ = J −1 F , G ⎝ ΔT ⎠

4.1. Active power loss

()

(24)

The iteration equation is defined as ⎜

⎝T

k+1





=

k ⎛p ⎞ ⎜

k



⎝T ⎠



k ⎛ Δp ⎞ ⎜

k



⎝ ΔT ⎠

,

+ Lgi

p1 p0 273 + Ti =0 v1 pi + p0 273 (30)

If the station is defined as node 1, then node 2 to node n are used to solve the equations iteratively, and the errors from node 2 to node n are calculated as

k+1 ⎛p ⎞

(27)

n



T ⎧ F = [F2 F3 ⋯ Fn ] ⎪G = [G2 G3 ⋯ Gn ]T , ⎨ p = [p2 p3 ⋯ pn ]T ⎪ T ⎩ T = [T2 T3 ⋯ Tn ]

∂F

Dg5

where pi and pj denote the pressures at nodes i and j (Pa), respectively, vNg,i is the volumetric flow, Dg is the diameter of pipeline (m). The equations of mass flow conservation at node is formulated as following:

The bold F, G, p, T are defined as follows:

⎛ ∂p J = ⎜ ∂G ⎜ ⎝ ∂p

2 lij vNg , ij

Next, the energy transmission losses are calculated. For electricity networks, the line loss (called active power loss) is calculated as follows [31]:

(25)

ΔPE , ij = Iij2 lij ρE , ij / Sij

where the superscript, k, refers to the iteration number. The iteration stops when the absolute maximum error is less than the set error.

(31)

where ρE,ij is the resistivity and Sij is the cross-sectional area of the conducting wire.

3.2. Electrical power flow model based on nodal voltage 4.2. Pressure drop energy loss The analysis of electricity distribution networks is based on the classical power flow analysis method. The model has been developed for decades. An electricity flow model includes two parts: the active power flow model and the reactive power flow model. These models are defined as follows [6]: n

For natural gas and heating networks, the pressure loss is calculated as follows [36]:

ΔPf , ij = mij Δpf , ij / ρij

where ΔPf,ij is the flow work consumption due to pressure loss, mij is the mass flow, ρij is the density of the fluid flowing in the pipeline, and ηcomp is the compressor efficiency. For the heating network, the flow work consumption is supplied by the circulating pump. Considering the mechanical efficiency, the pressure drop energy loss is formulated as following,

n

⎧ Pi = ei ∑ (Gij ej − Bij f ) + f ∑ (Gij f + Bij ej ) j i j ⎪ j=1 j=1 n n ⎨ ⎪Qi = fi ∑ (Gij ej − Bij f j ) − ei ∑ (Gij f j + Bij ej ) j=1 j=1 ⎩

(32)

, (26)

For each node, the real voltage, ei, and the imaginary voltage, fi, are unknown. Gij and Bij are the real and imaginary part of admittance. Pi and Qi are the active and reactive electricity. By defining the voltage and apparent power of the slack node as node 1, all unknown variables at node can be calculated. The Newton-Raphson algorithm is used to solve this set of equations.

ΔPP, ij =

ΔPf , ij ηpump

=

mij ·Δpf , ij ηpump ·ρij

(33)

where ηpump is the mechanical efficiency of the circulating pump. 4.3. Heat loss

3.3. Gas flow model based on nodal pressure Besides the pressure loss in heating networks, the heating loss is also calculated as follows [36]:

Similar to the electricity distribution networks, the energy flow 6

Energy Conversion and Management 206 (2020) 112491

D. Chen, et al.

λ l

0 ij − ⎧ ⎡ ⎤⎫ ΔPH , ij = mij ·Cp· Ti − ⎢Tam + (Ti − Tam) e Cp mij ⎥ , pi > pj ⎨ ⎬ ⎣ ⎦⎭ ⎩

(34)

The heating loss comes from heat exchange between the pipeline and the ambient environment. Therefore, the temperature drops along with the flow. 5. Case study The case study is based on the integrated district power–gas–heat energy system shown in Fig. 2. The electrical power network uses the IEEE-14 node system, the gas network uses a 23-node system, and the district heating network uses a 14-node system. In Fig. 2, node 1 is set as a slack node to balance the power–gas–heating network. At node one, CHP units are used for generating power and thermal energy by consuming natural gas. Therefore, WPs are used in the heating network and GCs in the gas network. At nodes 2 and 3, EBs are used for supplying thermal energy. At node 5, a GB is used for supplying thermal energy.

Fig. 3. District heating network diagram.

variables to be determined. Therefore, two independent variables and equations at each node are defined. The nodal conservation equations for mass and heating flow are derived, in which the transmission equations of mass flow and temperature drop are described by the pressure at each node. In the mass flow conservation equations in (13), the sum of inlet mass flow minus the sum of outlet mass flow is zero. The parameter at node one is used as a known reference. Therefore, for the other 13 nodes, 13 independent mass flow conservation equations are derived. For the heating flow conservation equations in (14), 13 independent nodal heating flow conservation equations are derived. Therefore, 26 equations are formulated in the nodal pressure model. The diagram of the traditional mass flow model is shown in Fig. 5. The regulation of this model is explained in [16,17]. The model includes three types of conservation equations: nodal mass flow, nodal heating flow, and loop pressure drop. As in the nodal pressure model, the conservation equations for nodal mass flow and nodal heating flow are derived. The mass flow is used directly in this model, whereas it is described by nodal pressure in the proposed model. Moreover, the loop pressure drop conservation equations are derived, in which the pressure drop is described by the mass flow in each pipeline. In the formulation of the mass flow model, the nodal temperature and mass flow in the pipelines are unknown variables. As there are twenty branches, twenty variables for mass flow in each pipeline are introduced. Compared with the nodal pressure model, the mass flow model uses twenty variables to describe the pressure at thirteen nodes, and the twenty mass flows in the pipelines. Therefore, thirteen independent mass flow conservation equations at each node and seven

5.1. Formulation complexity comparison between the proposed model and mass flow model for district heating networks This section is concerned with formulation complexity, adaptability to network topology changes, convergence speed, and effectiveness in setting initial iteration values. The convergence speed is measured by the number of iterations required for meeting a stopping condition to avoid computing time errors caused by performance differences of different software platforms and computers. 5.1.1. Formulation complexity comparison The modeling-related parameters, namely, pipeline diameter, length, and unit and thermal conductivity in heating networks, are provided in Appendixes A and B. User loads are provided in Appendixes C. In the analysis of heating networks, a typical 14-node district heating network is developed from Fig. 2 and is shown in Fig. 3. Node 1 is the slack node, and nodes 2, 3, and 5 are heating stations, which supply heat energy to other users. The formulation diagram of the nodal pressure model is shown in Fig. 4. In this model, nodal pressure and temperature are two unknown

Fig. 2. Integrated power-gas-heat energy system. 7

Energy Conversion and Management 206 (2020) 112491

D. Chen, et al.

Fig. 4. Formulation diagram of nodal pressure model.

are always twice as many equations as unknown nodes because only two unknown group variables (nodal pressure and temperature) can describe the conservation of nodal mass flow, nodal heating flow, and loop pressure drop. Comparing with the traditional mass flow model, the number of new equations equals to that of the loops in the network. When there is no loop, the two models have the same number of equations; however, more loop pressure drop conservation equations are required in the traditional mass flow model in a network with loops. In summary, the nodal pressure model allows a more concise formulation.

loop pressure drop conservation equations are derived. Combining these with the thirteen nodal heating flow conservation equations yields thirty-three independent equations for the mass flow model. By comparing the number of equations, it is clear that the nodal pressure model has the apparent advantage of fewer equations. In the nodal pressure model, the nodal pressure is directly used to define the variables, and therefore the loop pressure drop conservation equations are redundant. The conservation equations for nodal mass flow and loop pressure are integrated into the nodal pressure model, and thus the number of independent equations can be reduced. The number of equations, in this case, is reduced from 33 to 26. A more concise energy flow model can be obtained by the proposed method. Furthermore, it can be concluded that in the proposed model, there

5.1.2. Adaptation to network topology changes An analysis of network topology changes is always conducted in

Fig. 5. Formulation diagram of the traditional mass flow model. 8

Energy Conversion and Management 206 (2020) 112491

D. Chen, et al.

for nodal temperature, the Jacobian matrix has a size of 33 × 33. A smaller Jacobian matrix requires fewer iterations. The size is significantly reduced from 33 × 33 in the traditional mass flow model to 26 × 26 in the proposed model. Therefore, a higher convergence speed can be achieved.

network optimization or system planning, in which the locations of the nodes remain constant at each iteration, whereas the connections between the nodes change. Therefore, choosing nodal parameters as the variables for network analysis is more suitable and can reduce the complexity of the algorithm. For example, when the topology changes from Fig. 3 to Fig. 6, the only change is the addition of the pipeline connecting H6 to H14. In the nodal pressure model, the number of unknown variables remains two at each node. The total number of unknown variables is 26. The exact changes in the model formulation are shown in Fig. 7. Therefore, in the analysis of the changed network, resetting the connection between the nodes will enable the model to be applicable again at each iteration. However, in the traditional mass flow model, the branch of the connection between the two nodes changes. Moreover, the number of variables increases with the branches, and a new loop pressure conservation equation is required. The exact changes in the model formulation are shown in Fig. 8. Therefore, in the analysis of the changed network, both resetting the connection between the nodes and reformulating the loop pressure drop conservation are necessary at each iteration.

(2) Convergence speed Numerical simulations are carried out to evaluate the convergence speed of the two model formulations by quantifying the number of iterations. The convergence process is shown in Fig. 11. The logarithmic iteration errors are applied due to the small values obtained after several iterations. From Fig. 11, it is clear that the proposed model converges faster than the mass flow model. If the stopping condition is that the relative errors between two adjacent iterations should be 10−12, pressure and temperature converge within nine iterations in the proposed model. In contrast, mass flow and temperature converge within 11 and 12 iterations, respectively, in the mass flow model. This is primarily due to the smaller size of the Jacobian matrix in the proposed model than in the mass flow model. The temperature convergence speed is quite lower than that of the mass flow. As the temperature drop is calculated by the mass flow in the traditional mass flow model and by nodal pressure in the proposed model, the temperature converges after the pressure and mass flow. The absolute errors of the heating flow calculation results between the two models are presented in Fig. 12. Both the errors of nodal temperature and pressure are minimal, which is suitable for most circumstances in practice. The minimal errors show the correctness of the proposed model.

5.1.3. Convergence rate comparison The most widely used numerical method for the energy flow calculation is the Newton algorithm [16,17]. This algorithm is used in this study to solve the heating flow equations. The convergence speed is measured by the number of iterations required to meet a stopping condition to avoid computing time errors caused by performance differences of different software platforms and computers. In the Newton algorithm, the number of iterations is highly dependent on the size of the Jacobian matrix, and the convergence rate comparison is carried out by comparing the Jacobian matrices of the two models.

5.1.4. Effectiveness of initial iteration value setting When the difference between the initial iteration values and the results increases, the convergence rate will decrease. Therefore, the effectiveness of the initial iteration value setting is evaluated by comparing the relative fluctuation between the initial iteration values and the results for the equations of the heating flow model. According to engineering experience, the relative variation of nodal temperature is small at each node in heating networks. Therefore, the initial iteration value for the temperature at each node is set to be the temperature of the slack node, that is, 170 °C. The pressure of the slack node (node 1) is set to 1.0 MPa. According to the engineering experience, the relative variation range of nodal pressure is small. When the normalized pressure is applied, as the denominator in Eq. (16) cannot be zero, the initial iteration value for pressure at node i is

(1) Jacobian matrix for two heating flow models For the heating network shown in Fig. 3 and the equation sets in Fig. 4, when the Newton algorithm is applied, the Jacobian matrix of the proposed model is shown in Fig. 9. Both the mass flow conservation equations F and the energy flow conservation equations G are functions of the pressure and temperature. Therefore, the partial derivatives ∂F/ ∂p, ∂F/∂T, ∂G/∂p, and ∂G/∂T are not equal to zero. As there are 13 nodal pressures and 13 nodal temperatures, the Jacobian matrix has a size of 26 × 26. The Jacobian matrix of the traditional mass flow model is shown in Fig. 10. The formulation consists of the mass flow conservation equations F, loop pressure drop equations Δp, and energy flow conservation equations G. As there are 20 variables for mass flow and 13 variables

pi = 1 − (i − 1) × 0.001

Fig. 6. Changed district heating network diagram. 9

(35)

Energy Conversion and Management 206 (2020) 112491

D. Chen, et al.

Fig. 7. Formulation variation for the proposed model.

Fig. 8. Formulation variation for traditional mass flow model.

initial iteration values and the results is not significant. Unlike nodal pressure and temperature, the maximum mass flow cannot be determined directly. After several attempts, the maximum mass flow is set to 10 kg/s. When the algorithm is run, the initial iteration data for the mass flow in each pipeline is set to 10 kg/s. Then, the normalized mass flow in each pipeline is shown in Fig. 13(b). From Fig. 13(b), it can be seen that the mass flow varies greatly from the minimum value of –0.35 to the maximum value of 0.98. The absolute maximum variation amplitude is 135%. By numerical simulation, the convergence regions for initial iteration values in [0, 20] are shown in Fig. 14. It can be found that the convergence regions are [2.4, 3.9], [5, 5.7], [9.7, 10.2], [12.5, 13.4], and [15, 17.5]. Thus, in the mass flow model, with a probability of 30.5%, the initial iteration value leads to convergence. By comparing the relative fluctuation between the two models, it is clear that the maximum fluctuation of the nodal pressure is quite smaller than the maximum fluctuation of the mass flow in the pipelines. Therefore, the proposed model is more effective than the traditional mass flow model in setting initial iteration data. In general, in heating networks, the fluctuation of state parameters (such as pressure) between different nodes is smaller than that of operation parameters (such as mass flow). Therefore, it is easier to set the proper initial iteration values for state parameters. The proposed model uses a state parameter as the iterative variable, and thus it is more

Fig. 9. Jacobian matrix for the nodal pressure model.

By applying the proposed model, the normalized nodal pressure and temperature are determined, as shown in Fig. 13(a), where it can be seen that the minimum normalized pressure is 0.799 at node 11, and the minimum normalized temperature is 0.839 at node 13. The absolute maximum variation amplitude for pressure is 20.1%, and for temperature, it is 16.1%. The maximum relative fluctuation between the 10

Energy Conversion and Management 206 (2020) 112491

D. Chen, et al.

Fig. 10. Jacobian matrix for the traditional mass flow model.

effective in setting initial iteration values than the traditional mass flow model, which uses an operation parameter. 5.2. Energy transmission loss in integrated electricity–gas–heat energy systems Energy transmission losses comprise active power loss in power networks, heat loss in heating networks, and pressure loss in both natural gas and heating networks. These are analyzed in this section. To calculate the load of energy conversion equipment, the user load is converted according to Eqs. (1)–(6) for calculating the load of energy conversion equipment. In this process, the energy conversion equations of heating, natural gas, and electrical devices are applied sequentially. 5.2.1. Active power loss analysis in power distribution systems In this part, the active power loss of adding the load of multi-energy conversion devices are analyzed. The load of CHPs, WPs, and GCs is shown in Table 3. Other parameters are shown in Appendixes D and E. The active power loss change by adding the loads is determined in standard IEEE-14; detailed data can be found in [37], and the results are shown in Fig. 15. It can be seen from Fig. 15 that the loss primarily concentrates on feeders 1–2, 1–5, and 2–3. When the loads vary at nodes 1, 2, and 3 owing to the use of multi-energy conversion devices, the total active power loss increases from 13.393 MW to 13.545 MW even if the CHP supplies power to the substation. However, the increase is quite small. Therefore, the active power loss of multi-energy devices is small and can be neglected.

Fig. 11. Logarithmic convergence speed comparison between two models.

5.2.2. Heat loss analysis in district heating systems Heat loss stems from heat exchange between hot water in pipelines and ambient soil. Based on the pressure and temperature shown in Fig. 13(a) and the branch mass flow shown in Fig. 13(b), the loss is determined by Eq. (34). The loss and loss rate (the ratio of heat loss and total transmission heat energy) are shown in Fig. 16. In Fig. 16, it is evident that the heat loss varies slightly except in pipeline 4–9 owing to its more significant length. However, the loss at each branch is not more than 100 kW even when the total transmission heat is as high as 2822.84 kW. Except for pipelines 12–13 and 13–14, the loss rate is low and not

Fig. 12. Absolute errors of the heating flow calculation results between two models.

11

Energy Conversion and Management 206 (2020) 112491

D. Chen, et al.

Fig. 13. Normalized nodal pressure, temperature, and branch mass flow.

Fig. 14. Convergence regions for initial iteration values in [0, 20].

Fig. 15. Active power loss.

Table 3 Electrical power load of CHP, WP, and GC. Device

Location (node)

Load (MW)

CHP WP GC EB EB

1 1 1 2 3

−4.34482 0.0224 0.1739 0.4348 0.9239

more than nearly 10%. From Eq. (34), it can be concluded that smaller mass flow and longer transmission distance result in more considerable heat loss. In practice, branches with such a high loss rate can be deleted. 5.2.3. Pressure loss analysis in district heating and natural gas systems (1) Pressure loss in district heating system Pressure loss is the other form of energy transmission losses in heating systems and is determined by Eq. (33). Based on the nodal pressure shown in Fig. 13(a) and the mass flow shown in Fig. 13(b), the pressure loss at each branch is shown in Fig. 17.

Fig. 16. Heating transmission loss.

12

Energy Conversion and Management 206 (2020) 112491

D. Chen, et al.

flow models, the proposed method has the following advantages: (1) The nodal pressure heating flow model can simultaneously handle mass and thermal energy transmission. (2) By integrating loop pressure drop conservation and nodal mass flow conservation, a more concise formulation and higher convergence speed can be obtained. The nodal pressure model has an efficient formulation and can be applied for system planning. Variables and equations in the model are the same types, and their number is constant, thus allowing the analysis of heating networks under topology changes. (3) By using nodal pressure to define independent variables, the initial iteration values can be more conveniently set because of the narrowed relative fluctuation between the initial iteration data and the result. (4) The nodal pressure in heating and gas networks can be analogized to nodal voltage in electricity networks. Therefore, a unified energy network flow model, based on nodal pressure (voltage and pressure) to describe the integrated energy network systems, such as power-gas-heating network systems, is obtained. Further studies may be concentrated on the following two categories: 1) More types of energy flow models should be developed to describe more and more complex integrated energy systems; 2) The application of the proposed energy flow models, such as dispatching, network optimization, and energy system planning.

Fig. 17. Pressure loss analysis in district heating network.

It can be seen from Fig. 17 that the maximum pressure loss is 0.733 kW in pipeline 1–2 and is quite smaller than the thermal loss of 66.164 kW. The maximum ratio of thermal loss is less than 1.4%. From the perspective of energy loss, pressure loss can be neglected compared with heat loss. However, the pressure drop is essential in energy system planning and dispatching, and it is necessary to check user pressure restrictions.

CRediT authorship contribution statement Dongwen Chen: Conceptualization, Methodology, Software, Writing - original draft. Xiao Hu: Conceptualization, Methodology. Yong Li: Conceptualization, Supervision, Writing - review & editing. Ruzhu Wang: Supervision. Zulkarnain Abbas: Writing - review & editing. Shunqi Zeng: Funding acquisition, Project administration. Li Wang: Funding acquisition, Project administration.

(2) Pressure drop energy loss in gas networks Pressure loss is a significant energy transmission loss in natural gas systems. Based on the nodal pressure shown in Fig. 18, the pressure loss at each branch is shown in Fig. 19. By comparing the loss in the heating systems, it is found that the losses in these two networks are at the same level even if the flow medium is gas rather than water. The pressure loss in branch 1–2 is quite more significant than in other pipelines. To reduce the loss, larger diameter for pipeline 1–2 is an appropriate choice.

Declaration of Competing Interest The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper.

6. Conclusions

Acknowledgment

This study proposed a nodal pressure-based heating flow model for analyzing integrated energy systems. Compared with previous heating

This work is supported by the National Key R&D Program of China (2016YFB0901300).

Appendixes A. Roughness of pipeline materials

Material

Surface roughness (mm)

Glass Welded steel New cast iron Old cast iron Corrosion Galvanized iron Steel for water Steel for steam

0–0.02 0.015–0.06 0.15–0.50 1–1.5 0.25 0.25 0.50 0.20

B. Pipeline parameters in heating network

Pipeline

Length (m)

Diameter (m)

1–2 1–5 2–3 2–4

2000 2000 2000 2000

0.150 0.100 0.125 0.125

13

Unit thermal conductivity (W/m·°C) 0.2077 0.1578 0.1830 0.1830 Fig. 18. Normalized nodal pressure in gas network.

Energy Conversion and Management 206 (2020) 112491

D. Chen, et al.

Fig. 19. Pressure drop energy loss analysis in gas network. 2–5 3–4 4–5 4–7 4–9 5–6 6–11 6–12 6–13 7–8 7–9 9–10 9–14 10–11 12–13 13–14

2100 2000 2000 2000 3000 2000 2000 2000 2200 2000 2200 2000 2100 2000 2000 2000

0.125 0.125 0.125 0.125 0.125 0.125 0.100 0.125 0.125 0.125 0.125 0.125 0.125 0.080 0.125 0.080

0.1830 0.1830 0.1830 0.1830 0.1830 0.1830 0.1578 0.1830 0.1830 0.1830 0.1830 0.1830 0.1830 0.1373 0.1830 0.1373

C. Natural gas and user heating loads

Node

Natural gas load (Nm3)

Heating load (kW)

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23

0.3550 0.1200 0.1380 0.1500 0.0432 0.2100 0.2160 0.2400 0.2820 0.2100 0.1500 0.2100 0.2100 0.1200 0.3120 0.1800 0.1500 0.1650 0.1500 0.1560 0.2400 0.1620 0.2250

\ −400 −850 900 −1350 600 600 870 760 930 1260 750 500 700 \ \ \ \ \ \ \ \ \

D. Pipeline parameters in gas network

Pipeline

Length (m)

Diameter (m)

1–2 1–5

1600 2100

0.50 0.40

14

Energy Conversion and Management 206 (2020) 112491

D. Chen, et al. 2–3 2–4 2–18 3–7 3–15 4–5 4–7 4–10 5–6 6–11 6–12 6–13 7–8 7–9 9–10 9–14 10–11 12–13 13–14 15–16 16–17 17–20 18–19 18–23 19–20 19–22 20–21

1000 2100 1600 1600 1600 1600 1000 1600 1600 1600 1600 1600 1600 1600 800 1600 1200 1600 1600 1600 1600 1000 1600 1000 1600 1000 1000

0.35 0.30 0.35 0.30 0.30 0.20 0.20 0.30 0.30 0.20 0.20 0.20 0.20 0.30 0.20 0.20 0.20 0.20 0.20 0.30 0.20 0.25 0.30 0.20 0.20 0.20 0.20

E. Fluid and device parameters

Items

Parameters

Value

Unit

Γ Cosφ Ρ Μ Cp Gas boiler CHP CHP EB WP GC

Calorific value Power coefficient Water density Kinematic viscosity Specific heat Efficiency Generating efficiency Heat recovery efficiency Efficiency Efficiency Adiabatic efficiency

34 0.95 853 1.7782 × 10−6 4.5476 0.92 0.36 0.85 0.92 0.65 0.85

MJ/Nm3 \ Kg/m3 m2/s kJ/(kg·°C) \ \ \ \ \ \

calculation by the Newton method. Appl Therm Eng 2013;61(2):163–70. [13] Ancona MA, Bianchi M, Branchini L, Melino F. District heating network design and analysis. Energy Procedia 2014;45:1225–34. [14] Coulson, Richardson JM, Backhurst JF, Harker JR, JH. Coulson & richardson's chemical engineering: fluid flow, heat transfer and mass transfer 1992;1, 47(2), 513–4. [15] Shabanpour-Haghighi A, Seifi AR. An integrated steady-state operation assessment of electrical, natural gas, and district heating networks. IEEE Trans Power Syst 2015;31(5):3636–47. [16] Liu X, Wu J, Jenkins N, Bagdanavicius A. Combined analysis of electricity and heat networks. Appl Energy 2016;162:1238–50. [17] Liu X, Mancarella P. Modelling, assessment and Sankey diagrams of integrated electricity-heat-gas networks in multi-vector district energy systems. Appl Energy 2016;167:336–52. [18] Chen Y, Xu Y, Li Z, Feng X. Optimally coordinated dispatch of combined-heat-andelectrical network with demand response. IET Gener Transm Distrib 2019;13(11):2216–25. [19] Juanwei C, Tao Y, Yue X, Xiaohua C, Bo Y, Baomin Z. Fast analytical method for reliability evaluation of electricity-gas integrated energy system considering dispatch strategies. Appl Energy 2019;242:260–72. [20] Zhang F, Cheng CS. A modified Newton method for radial distribution system power flow analysis. IEEE Trans Power Syst 1997;12(1):389–97. [21] Low SH. Convex relaxation of optimal power flow—part I: Formulations and equivalence. IEEE Trans Control Network Syst 2014;1(1):15–27. [22] Low SH. Convex relaxation of optimal power flow—part II: Exactness. IEEE Trans Control Network Syst 2014;1(2):177–89. [23] Tinney WF, Hart CE. Power flow solution by Newton's method. IEEE Trans Power Apparatus Syst 1967;11:1449–60. [24] Ma T, Wu J, Hao L. Energy flow modeling and optimal operation analysis of the micro energy grid based on energy hub. Energy Convers Manage 2017;133:292–306. [25] Chen S, Wei Z, Sun G, Sun Y, Zang H, Zhu Y. Optimal power and gas flow with a limited number of control actions. IEEE Trans Smart Grid 2017;9(5):5371–80. [26] Lu S, Gu W, Zhou J, Zhang X, Wu C. Coordinated dispatch of multi-energy system with district heating network: Modeling and solution strategy. Energy

References [1] Liu T, Zhang D, Wang S, Wu T. Standardized modelling and economic optimization of multi-carrier energy systems considering energy storage and demand response. Energy Convers Manage 2019;182:126–42. [2] Asl DK, Hamedi A, Seifi AR. Planning, operation and flexibility contribution of multi-carrier energy storage systems in integrated energy systems. IET Renew Power Gener 2019. [3] Narula K, Villasmil W, Patel MK. Simulation method for assessing hourly energy flows in district heating system with seasonal thermal energy storage. Renewable Energy 2019. [4] Hyman B, Ozalp N, Varbanov PS, Fan YV. Modelling energy flows in industry: general methodology to develop process step models. Energy Convers Manage 2019;181:528–43. [5] Qiao Z, Guo Q, Sun H, Pan Z, Liu Y, Xiong W. An interval gas flow analysis in natural gas and electricity coupled networks considering the uncertainty of wind power. Appl Energy 2017;201:343–53. [6] 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):2156–66. [7] Shao C, Wang X, Shahidehpour M, Wang X, Wang B. An MILP-based optimal power flow in multicarrier energy systems. IEEE Trans Sustainable Energy 2016;8(1):239–48. [8] Abeysekera M, Wu J, Jenkins N, Rees M. Steady state analysis of gas networks with distributed injection of alternative gas. Appl Energy 2016;164:991–1002. [9] 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 2016;32(3):1970–9. [10] Zlotnik A, Roald L, Backhaus S, Chertkov M, Andersson G. Coordinated scheduling for interdependent electric power and natural gas infrastructures. IEEE Trans Power Syst 2017;32(1):600–10. [11] Qiu J, Zhao J, Yang H, Dong ZY. Optimal scheduling for prosumers in coupled transactive power and gas systems. IEEE Trans Power Syst 2018;33(2):1970–80. [12] Wang W, Cheng X, Liang X. Optimization modeling of district heating networks and

15

Energy Conversion and Management 206 (2020) 112491

D. Chen, et al.

[32] Wang J, Hu Z, Xie S. Expansion planning model of multi-energy system with the integration of active distribution network. Appl Energy 2019;253:113517. [33] Housos EC, Irisarri G, Porter RM, Sasson AM. Steady state network equivalents for power system planning applications. IEEE Trans Power Apparatus Syst 1980;99(6):2113–20. [34] Munson BR, Okiishi TH, Huebsch WW, Rothmayer AP. Fluid mechanics. Singapore: Wiley; 2013. [35] Jiang XS, Jing ZX, Li YZ, Wu QH, Tang WH. Modelling and operation optimization of an integrated energy based direct district water-heating system. Energy 2014;64:375–88. [36] Bejan A. Advanced engineering thermodynamics. John Wiley & Sons; 2016. [37] Zia MF, Elbouchikhi E, Benbouzid M. Microgrids energy management systems: a critical review on methods, solutions, and prospects. Appl Energy 2018:222.

2018;152:358–70. [27] Liu X, Yan Z, Wu J. Optimal coordinated operation of a multi-energy community considering interactions between energy storage and conversion devices. Appl Energy 2019;248:256–73. [28] Chaudry M, Jenkins N, Qadrdan M, Wu J. Combined gas and electricity network expansion planning. Appl Energy 2014;113:1171–87. [29] Russell B, Seth B, Rene VHP, Conrado BS, Mehdi S. Joint electricity and natural gas transmission planning with endogenous market feedbacks. IEEE Trans Power Syst 2018;33(6):6397–409. [30] Ma T, Wu J, Hao L, Lee WJ, Yan H, Li D. The optimal structure planning and energy management strategies of smart multi energy systems. Energy 2018;160:122–41. [31] Li Y, Feng B, Li G, Qi J, Zhao D, Mu Y. Optimal distributed generation planning in active distribution networks considering integration of energy storage. Appl Energy 2018;210:1073–81.

16