Zhang et al. (2021) proposed a two-layer quadratic curve optimal energy flow model. The optimal flow model of conical AC power was established in the upper level model. And the power generation cost of natural gas generator units was calculated based on the marginal price of natural gas. Jia et al. (2020) proposed a probabilistic energy flow convex optimization method for the integrated energy system. The piecewise linear approximation method was used to linearize the nonlinear objective function. Thereby the probabilistic energy flow calculation model was transformed into a linear programming model. The gas turbines were applied to improve power generation and consumption balance in short-term dispatching to achieve closer coordination between power system and natural gas system (Mirzaei et al., 2020). Zeng et al. (2019) proposed an integrated electricity and natural gas systems (IEGS) model that takes into account electricity to natural gas facilities and gas storage facilities. In addition, sequential Monte Carlo (SMC) method was used to evaluate the reliability of IEGS. Nasiri et al. (2023) proposed a two-layer dispatching model based on the multiple energy service providers (MESP). MESP minimizes the cost of purchasing electricity and natural gas by operating energy storage systems and demand response plans (DRPs). A Lagrangian relaxation multi-step strategy was proposed to realized collaborative optimization of IEGS (Faridpak et al., 2020).

Qi et al. (2017) proposed a flexible model to evaluate the performance of (IEGS) in extreme weather conditions, namely the energy assistance between the power system and the natural gas network. Shao et al. (2017) proposed a comprehensive planning algorithm for IEGS to enhance the resilience of the power system under extreme conditions. And a variable uncertainty set has been developed to describe the interaction between power system expansion status and extreme events. Correa-Posada and Sánchez-Martın (2014) proposed an optimal power flow calculation model, which considers the N-1 criteria for IEGS. And a fault analysis method for natural gas systems was established to realized rapid fault analysis based on the linear sensitivity factor method. Meng et al. (2019) proposed a sequential energy flow analysis method for IEGS based on multiple balance nodes. Ni et al. (2016) proposed a multi energy optimal power flow solution method for IEGS. Based on the concept of energy hubs, the coupling relationship between different energy systems was considered.

Qin et al. (2019) proposed a hybrid multi-objective optimization and game theory approach (HMOGTA) to realized optimized operation of an integrated energy system consisting of an electric and natural gas (E&G) utility network, multiple distributed energy stations (DESs), and multiple energy users (EU). Unified models and data-driven methods: However, due to the different characteristics of control mechanisms, network topology, and dynamic processes between power system and gas system, it is difficult to realized the collaborative optimization of IEGS. In order to solve this problem, a unified transient energy flow model was proposed based on the singular perturbation theory and dynamic model of IEGS (Huang et al., 2022). Wang et al. (2021) proposed a data-driven opportunity constrained method for calculating the optimal power flow of IEGS. Tabebordbar et al. (2023) proposed an optimization framework to control the optimal scale of cogeneration and P2G technology in IEGS. Gao et al. (2023) designed a rolling scheduling framework for IEGS. Based on the framework a partial differential equation mathematical model was proposed. The dynamics of natural gas pipelines and the operational characteristics of P2G facilities can be expressed accurately to maintain the authenticity of IEGS in the rolling scheduling process. Wang et al. (2020) proposed a two-stage low-carbon operational planning model based on bilateral trading mechanism and active demand side management (DSM), aiming to reduce carbon emissions. The carbon emission flow (CEF) model is used to track emissions and calculate carbon intensity, considering energy storage systems.

Wu et al. (2022) proposed a novel carbon-oriented expansion planning model of IEGS with FCSs to determine the optimal alternatives, locations, and sizes for ecofriendly candidate assets, including roof-top PV panels and fuel cell (FC) units in each FCS, as well as renewable energy units and carbon capture and storage (CCS) systems in IEGS. Salvatti et al. (2024) an optimized power management strategy for FCS with integrated battery energy storage systems (BESS) was proposed. The proposed strategy aimed to monitor the variation in AC voltage at the point of common coupling (PCC) and the state of charge (SOC) of the BESS, with the objective of establishing a deterministic formulation to find the optimal instantaneous BESS power level. Shi et al. (2022) proposed a novel sensitivity analysis-based FCS planning approach, which considers the voltage sensitivity of each sub-network in the distribution network and charging service availability for EV drivers in the transportation network.

The above literature mainly studies the issue of FCS integration into IEGS from a planning perspective. But the operation guidance for IEGS is not provided. Li et al. (2021) proposed a bi-level optimal scheduling model for community integrated energy system (CIES) with an charging station in multi-stakeholder scenarios, where an integrated demand response program comprising a dynamic pricing mechanism was designed. Weng et al. (2020) proposed a joint planning model of active distribution network and transportation network including electricity, gas, heat, and traffic loads. Zhang et al. (2024) proposed a two-stage robust operation method of electricity-gas-heat integrated MEMGs considering heterogeneous uncertainties was proposed to coordinate multiple energy carriers. Jiang et al. (2023) has conducted research on the coordinated scheduling of IEGS from the perspective of carbon emissions. And the committed carbon emission operation region (CCEOR) was proposed. Zheng et al. (2021) peoposed a data-driven integrated electricity-gas system stochastic co-optimisation model to co-optimise these two energy systems for day-ahead market clearance. For the public transport in IEGS, Zhu et al. (2023) proposed coordination planning of IEGS and charging station considering carbon emission reduction in this article that responds effectively to emerging concerns with energy and transportation integration inefficiencies and is less environmentally friendly.

Although the above research has conducted a systematic study on IEGS, it has ignored the load characteristics of DN. More and more electric vehicles (EVs) will be connected to distribution network (DN), especially the connection of fast charging stations (FCSs), which will have a negative impact on the operation of DN. In response to the shortcomings of the above research, the main contributions of this paper are as follows:

(1) Based on DTALite simulation, the travel status of vehicles is simulated. The process of EV travel and charging is simulated through hourly and timely rolling to achieve multi time scale charging load modeling. The charging price follows the principle of hourly rolling based on TOU. The changes of travel path will be simulated due to charging. And the charging load model is established.

In this paper, the model of EVs can be divided into three parts for modeling, namely the initial power model, the travel power model, and the charging power model.

(1) The initial power model

The initial power model represents the remaining battery capacity of EVs when the EVs departs. The initial power model is shown in Equation 1. S EV = S EV ini (1) where S EV ini is the initial battery remaining power of EV; S EV is the battery remaining power of EV.

(2) The travel power model

In this paper, EVs will choose nearby fast charging stations for charging. Therefore, the charging load is related to the number of vehicles going to FCS for charging. The charging load model is shown in Equations 4–5. P FCS = ∑ i FCS = 1 n FCS μ i FCS p pile (4) μ i FCS = 1 Charging pile used 0 Charging station used (5) where P FCS is the charging load of FCS; μ i FCS is a binary variable representing the occupancy of charging pile; i FCS is the ID of the fast charging pile; n FCS is the number of charging piles within FCS.

The Weymouth equation is an equation used to calculate the flow characteristics of fluids in pipelines under flow pressure and temperature conditions. This equation considers the influence of fluid compressibility on flow measurement, especially under high pressure and high flow conditions. The Weymouth equation can help engineers predict and calculate fluid flow in pipelines more accurately. In natural gas supply systems, the Weymouth equation is usually used to describe the pipeline retention state of the natural gas network which is described in Equation 6. G m n , ave , t G m n , ave , t = C m n 2 π m , t 2 − π n , t 2 (6) where G m n , ave , t is the average flow of natural gas in pipeline m n at time t ; π m , t and π n , t are the gas pressures at node m and node n of the natural gas pipeline respectively; C m n is the Weymouth constant, which is related to physical properties such as pipeline temperature, pipeline length, pipeline diameter, and pipeline friction.

The flow balance constraints of the natural gas network are as follows: ∑ n ∈ R P n G m n , t e − G m n , t h + ∑ s ∈ R s n G s , t o u t − G s , t i n + ∑ c ∈ R C n G c , t C = ∑ l ∈ R G L n G l , t L + ∑ g ∈ R g n G g , t G T (7) In Equation 7, R P n is the set of natural gas pipelines connected to node n ; R s n is the set of natural gas storge connected to node n ; R C n is the set of natural valve stations connected to node n ; R G L n is the set of natural gas load connected to node n ; R g n is the set of gas turbine connected to node n ; G m n , t e and G m n , t h are the gas flow at the end and beginning of the natural gas pipeline m n at time t respectively; G s , t o u t and G s , t i n are the release and intake volume of the natural gas storage at time t ; G c , t C is the natural gas flow of the natural valve station at time t ; G l , t L is the gas flow of the natural gas load connected to node n at time t; G g , t G T is the gas flow of the gas turbine connected to the node n at time t .

Meanwhile, a linear storage model is applied to describe the dynamic characteristics of natural gas pipelines. R m n , t = L m n π m n , ave , t (8) R m n , t = R m n , t − 1 + G m n , t h − G m n , t e (9) π m n , ave , t = π m , t + π n , t 2 (10) G m n , ave , t = G m n , t e + G m n , t h 2 (11) where R m n , t is the natural gas reserve of the pipeline m n at time t ; π m n , ave , t is the average pressure of natural gas pipeline m n ; L m n is the length of the pipeline.

The pressure constraints at natural gas pipeline nodes, pipeline operation constraints, and gas supply constraints at natural gas valve stations are defined in Equations 12–14: π t , ⁡ min ≤ π t ≤ π t , ⁡ max (12) R m n , t , ⁡ min ≤ R m n , t ≤ R m n , t , ⁡ max (13) G c , t , ⁡ min C ≤ G c , t C ≤ G c , t , ⁡ max C (14) where π t , ⁡ max and π t , ⁡ min are the upper and lower limits of the node gas pressure respectively; R m n , t , ⁡ max and R m n , t , ⁡ min are the upper and lower limits of the pipeline retention respectively; G c , t , ⁡ max C and G c , t , ⁡ min C are the upper and lower limits of natural valve station gas supply volume.

The constraints of natural gas storage devices are defined in Equations 15–20: R s , t = R s , t − 1 − w s , t o u t G s , t o u t η s , o u t + w s , t i n G s , t i n η s , i n (15) G s , t . ⁡ min o u t ≤ G s , t o u t ≤ G s , t , ⁡ max o u t (16) G s , t . ⁡ min i n ≤ G s , t i n ≤ G s , t , ⁡ max i n (17) w s , t o u t = 0 Not releasing gas 1 Releasing gas (18) w s , t i n = 0 Not inflating gas 1  Inflating gas (19) 0 ≤ w s , t i n + w s , t o u t ≤ 1 (20) where R s , t is the storage capacity of the natural gas storage device at time t ; η s , i n and η s , o u t are the inflation efficiency and release efficiency of the gas storage device respectively; w s , t o u t and w s , t i n are binary variables that represent the inflation and deflation status of the gas storage device; G s , t , ⁡ max o u t and G s , t . ⁡ min o u t are the maximum and minimum release capacities of the gas storage device respectively; G s , t , ⁡ max i n and G s , t . ⁡ min i n are the maximum and minimum inflation capacities of the gas storage device respectively.

In this paper, a collaborative optimization dispatching strategy is proposed for IEGS with FCSs to realize safe and economical operation of the system.

Objection function:

The objection function is divided into two parts: the operating costs of the power system and the natural gas network. And the specific expression is as follows: min ∑ t = 1 T ∑ c = 1 N c λ c , t G c , t C + ∑ t = 1 T λ P , t P t G r i d (21) In Equation 21, λ c , t is the gas supply price of the natural valve station at time t; λ P , t is the TOU of DN; P t G r i d is the purchase electricity of DN; T is the number of time periods.

Distribution network constraints are defined in Equations 22–27. U min 2 ≤ V i , t ≤ U max 2 (22) 0 ≤ l i j , t ≤ l i j , ⁡ max (23) P i j , t + ∑ i G = 1 N G κ j P t , i G = r i j l i j , t + p j , t + ∑ j : j → k P j k , t + ∑ i FCS = 1 N FCS υ j P t , i FCS (24) Q i j , t = x i j l i j , t + q j , t + ∑ j : j → k Q j k , t (25) V i , t + l i j , t r i j 2 + x i j 2 = 2 r i j P i j , t + x i j Q i j , t + V j , t (26) 2 P i j , t 2 Q i j , t V i , t − l i j , t 2 ≤ V i , t + l i j , t (27) where R is the set of DN nodes; P i j , t is the active power of branch i j ; Q i j , t is the reactive power of branch i j ; V i , t is the square of the voltage amplitude; l i j , t is the square of the current amplitude; U min and U max are the minimum and maximum values of node voltage amplitude respectively; r i j is the resistance value of branch i j ; x i j is the reactance value of branch i j ; i G is the ID of gas turbine; N G is the number of gas turbine; P t , i G is the active power output of gas turbine; κ j is a binary variable that indicates whether the gas turbine is installed at node j ; p j , t is the active load of node j ; q j , t is the reactive load of node j ; υ j is a binary variable that indicates whether FCS is installed at node j .

Natural gas network constraints:

The gas supply system model is defined in this paper. The constraint conditions for the natural gas network are (3)–(8).

At the same time, the adjustment constraints for natural gas valve stations include natural gas supply flow constraints, supply volume adjustment constraints, supply volume adjustment frequency constraints, and supply flow limit constraints. The adjustment constraints for natural gas valve stations are defined in Equations 28–34. w t + + w t − Δ G c , t C ≤ G c , t C − G c , t − 1 C ≤ w t + + w t − Δ G c , t C (28) 0 ≤ w t + + w t + 1 − ≤ 1 (29) 0 ≤ w t − + w t + 1 + ≤ 1 (30) 0 ≤ w t + + w t + 1 + ≤ 1 (31) 0 ≤ w t − + w t + 1 − ≤ 1 (32) 0 ≤ w t + + w t − ≤ 1 (33) 0 ≤ ∑ t T w t + + w t − ≤ N max C (34) where Δ G c , t C is the amount of natural gas valve station climbing and landslide; w t + is a binary variable, when w t + = 1 , the natural gas valve station will increase the gas supply; when w t + = 0 , the natural gas valve station will not increase the gas supply; w t − is a binary variable, when w t − = 1 , the natural gas valve station will reduce the gas supply; when w t − = 0 , the natural gas valve station will not reduce the gas supply; N max C is the maximum number of times a natural gas valve station can adjust its gas supply.

In terms of charging load simulation calculation, based on simulating vehicle travel, it is assumed that EVs will go to the nearest charging station to complete charging. Therefore, the simulation process of vehicle travel and vehicle charging is shown in the Figure 1. And DTALite is applied for simulation of vehicle travel.

In terms of optimizing the calculation method of dispatching strategy, the calculation process is shown in Figure 2. In Figure 2, the calculation process is divided into hour time part and real time part. In hour time part, the vehicles travel plan for each hour within a day will be counted. In real time part, the travel of vehicles will be simulated. And the FCSs selection of EVs will be simulated every 15 min. And the process of generating travel data is shown in Figure 3.

The steps of the algorithm are as follows:

Step 1: Calculate the charging load based on the proposed charging load calculation method. The charging plan for EVs will not be changed. And EVs will go to the fixed FCSs for charging.

Step 2: (1) Input distribution network related data: load data, distribution network topology data.

(2) Input gas network related data: gas load data, gas network topology data.

Step 3: Optimize calculations: The optimization objective is to minimize the operation cost of IEGS. The constraints include DN operation constraints, GN operation constraints, FCSs constraints, and gas turbine constraints.

Objection function: (21)

Step 4: Output the dispatching results of IEGS. The dispatching results of IEGS include the operating status of DN and the operation costs of IEGS.

The voltage distribution of FCSs nodes is shown in Figures 9, 10.

For FCS 1 node, when DN is coupled with the gas network and equipped with energy and gas storage devices, the voltage of FCS 1 node does not exceed the limit. At the same time, when DN is coupled with the gas grid without energy storage and gas storage devices, the voltage of FCS 1 node also does not exceed the limit. When the distribution network is not coupled with the gas network, according to scenario 2 and scenario 3, it can be seen that even with energy storage installed, it is difficult to meet the safe and stable operation requirements of DN. And the voltage exceeds the lower limit. From scenarios 3 and 4, it can be seen that equipping energy storage devices in DN can to some extent alleviate the problem of voltage exceeding the lower limit caused by charging loads, but this problem has not been solved.

For FCS 2 node, the voltage does not exceed the lower limit. The voltage distribution in scenario 1 and scenario 4 is better than that in scenario 2 and scenario 3. When DN is coupled with the gas network, it can improve the operation of DN. With the support of gas turbines, the operation pressure of DN is alleviated.

From the voltage distribution, it can be seen that equipping a certain capacity of ES within DN cannot solve the voltage limit problem caused by charging loads. When the gas turbines are applied to realize system coupling between DN and GN, the voltage will be maintained within a good range. Therefore, the coupled operation of DN and the GN can better support the stable operation of DN.

The comparison of electricity purchase costs for DN is shown in Figure 11. The comparison of gas purchase costs for gas networks is shown in Figure 12. The comparison of energy purchase costs for IEGS is shown in Figure 13.

In Figure 11, when DN and gas network are coupled through gas turbines, the gas network can provide energy to DN to reducing the amount of electricity purchased of DN. The purchase of electricity has been reduced. When DN and gas network operate independently, DN lacks energy support from gas network, so more electricity needs to be purchased to meet the load demand. The purchase of electricity has been increased.

For Figure 12, the gas network does not need to provide energy support to DN in scenario 2 and scenario 5. It only needs to meet the gas load demand of the gas network, so the cost of purchasing gas is relatively low. When the energy coupling between the gas network and the distribution network is realized through gas turbines, the gas network can support the energy to DN, so the gas purchase cost of the gas network will be increased.

For Figure 13, due to the independent operation of the distribution network and gas network in scenario 2 and scenario 5, there is a lack of energy support between the two systems. The gas network cannot provide energy support to DN based on TOU, which increases the overall operating cost of IEGS.

When the gas turbines are applied to realize system coupling between DN and GN, the electricity purchase cost of DN will be reduced significantly. But the natural gas purchase cost of DN will be increased significantly. For the energy purchase cost of IEGS, when the stable state of the systems is ignored, the energy purchase cost of IEGS will be decrease significantly owing to the support to DN from GN.

The comparative distribution of time-sharing active power loss and total active power loss for DN is shown in Figures 14, 15.

In Figure 14, the active power loss of DN coupling with the gas network is much smaller than the active power loss of the independent operation of DN. In Figure 15, the total active power loss of DN coupling with the gas network is much smaller than the active power loss of the independent operation of DN. At the same time, when energy storage and gas storage equipment are equipped in the integrated energy system of electricity and gas coupling, it reduces the active power loss of the distribution network to a certain extent.