Gas network simulation method and system based on electricity-gas dynamic analysis

Abstract

The invention relates to the technical field of gas network control, in particular to a gas network simulation method and system based on electricity-gas dynamic analysis, and the method comprises the steps: building a natural gas dynamic model based on a variable coefficient partial differential equation; calculating the gas pressure, the mass flow rate, the friction coefficient and the gas flow rate at the initial moment; determining boundary conditions of the natural gas dynamic model; and solving the natural gas dynamic model by adopting a three-stage leapfrog finite difference method to obtain the gas pressure, the mass flow rate, the friction coefficient and the gas flow rate of the natural gas at different moments in the pipeline. Accurate simulation of the natural gas network in the comprehensive electric power and natural gas system is achieved, so that natural gas resources can be distributed more reasonably, the energy utilization efficiency is improved, the operation cost of the gas network is reduced, and the reliability and economical efficiency of dispatching of the comprehensive electric power and natural gas system are improved.

Description

Technical Field

[0001] The present invention relates to the technical field of gas pipeline network control, and in particular to a gas network simulation method and system based on electro-gas dynamic analysis. Background Technique

[0002] The integrated electric and gas system (IEGS), as a representative of the integrated energy system, is a key strategy to promote the utilization of renewable energy and improve environmental and economic benefits. In recent years, many coal-fired power units have been replaced by cleaner and more flexible gas-fired power units, resulting in an increasingly strong coupling between the electric power network (EPN) and the natural gas network (NGN). However, the strong coupling relationship between the electric power network and the natural gas network significantly amplifies the system complexity, and inaccurate assessment of the operating state changes in the natural gas network may lead to serious consequences. Therefore, in order to ensure the reliability and economy of system scheduling, it is necessary to accurately simulate and model the natural gas flow rate.

[0003] Most studies on natural gas systems use the Weymouth equations (WEs) to approximate the gas flow model to reduce the computational complexity. This model relates gas pressure and mass flow rate (MFR) using constant coefficients, and this simplified model has been widely used in various fields of the energy system. However, the Weymouth equations ignore the dynamic diffusion process of the gas and assume that the mass flow rate in the gas network is equal at any position at any time, which does not conform to the actual situation. Therefore, it is necessary to introduce a natural gas dynamic model. The natural gas dynamic model is usually described by a set of partial differential equations (PDEs), which makes it complex to derive an analytical natural gas dynamic model.

[0004] To solve this problem, numerical methods have been used to approximately solve the natural gas dynamic model, and certain achievements have been made. For example, the finite difference method (FMD) based on the Euler format is used to discretize the simplified partial differential equations; the Wendroff finite difference method is used to convert the original partial differential equations into a set of equivalent algebraic equations; the finite difference method based on the central format is used to approximately model the gas flow or the Laplace transform is used to transform the gas dynamic model into an ordinary differential equation (ODE) with constant coefficients in the frequency domain to reduce the mathematical complexity, etc.

[0005] The above research has explored the natural gas dynamic model from various perspectives, but they are all based on the gas flow model and constant coefficient partial differential equations. Although it simplifies the computational complexity, it inevitably affects the accuracy of the model, reduces the simulation accuracy of the gas network, and ultimately affects the operation of the integrated power and natural gas system. Summary of the Invention

[0006] Therefore, the technical problem to be solved by the present invention is to overcome the problem that the simulation accuracy and efficiency of the natural gas dynamic model in the prior art are low, resulting in reduced reliability and economy of the integrated power and natural gas system scheduling.

[0007] To solve the above technical problems, the present invention provides a gas network simulation method based on electro-gas dynamic analysis, including:

[0008] Constructing a natural gas dynamic model based on variable coefficient partial differential equations, including:

[0009] Constructing an initial natural gas dynamic model based on the continuity equation, momentum equation and state equation of one-dimensional hydrodynamics; constructing a friction coefficient equation based on the mass flow rate, and constructing a gas velocity equation based on the mass flow rate and gas pressure; expressing the friction coefficient and gas velocity in the initial natural gas dynamic model as the friction coefficient equation and gas velocity equation to obtain the natural gas dynamic model;

[0010] Calculating the gas pressure, mass flow rate, friction coefficient and gas velocity at the initial moment;

[0011] Solving the natural gas dynamic model based on the boundary conditions of the natural gas dynamic model to obtain the gas pressure, mass flow rate, friction coefficient and gas velocity of natural gas at different moments in the pipeline.

[0012] Preferably, the friction coefficient equation is constructed based on the mass flow rate, and the formula is:

[0013]

[0014] Where represents the friction coefficient equation, D and A respectively represent the pipeline diameter and cross-sectional area, c1, c2 and c3 are empirical parameters, M represents the mass flow rate, λ and μ respectively represent the pipeline roughness and gas dynamic viscosity;

[0015] The gas velocity equation is constructed based on the mass flow rate and gas pressure, and the formula is:

[0016]

[0017] Where represents the gas velocity equation, P represents the gas pressure, and c represents the speed of sound;

[0018] The formula of the natural gas dynamic model is expressed as:

[0019]

[0020] P = c 2 ρ

[0021] where x and t represent the spatial coordinate and the time coordinate respectively, and ρ represents the gas density.

[0022] Preferably, calculating the gas pressure, mass flow rate, friction coefficient and gas flow velocity at the initial moment includes:

[0023] Initializing the values of the gas flow velocity and gas pressure at the initial moment; the mass flow rate at the initial moment is the function value of the mass flow rate change at the pipeline outlet at the initial moment; and calculating the friction coefficient at the initial moment according to the function value of the mass flow rate change at the pipeline outlet at the initial moment;

[0024] Iteratively updating the gas pressure at the initial moment according to the function value of the mass flow rate change at the pipeline outlet at the initial moment, the gas flow velocity at the initial moment and the friction coefficient, and then iteratively updating the gas flow velocity at the initial moment according to the updated gas pressure at the initial moment until the maximum number of iterations is reached or the change amount of the initial gas pressure and the initial gas flow velocity is less than the preset threshold.

[0025] Preferably, calculating the friction coefficient at the initial moment according to the function value of the mass flow rate change at the pipeline outlet at the initial moment, the formula is:

[0026]

[0027] where, represents the friction coefficient at the initial moment at the spatial coordinate i, D and A represent the diameter and cross-sectional area of the gas pipeline respectively, c1, c2 and c3 are empirical parameters, λ and μ represent the pipeline roughness and the gas dynamic viscosity respectively; M outlet (t0) represents the function value of the mass flow rate change at the pipeline outlet at the initial moment;

[0028] Iteratively updating the gas pressure at the initial moment according to the function value of the mass flow rate change at the pipeline outlet at the initial moment, the gas flow velocity at the initial moment and the friction coefficient, the formula is:

[0029] P j,0 = P j-1,0 + ε i,0 M Outlet (t0)

[0030]

[0031] where, P j,0 and P j-1,0respectively represent the gas pressures at the initial moment at spatial coordinate j and spatial coordinate j - 1, ε i,0 is a non - zero coefficient, represents the gas flow velocity at the initial moment at spatial coordinate i, and h represents the spatial difference step;

[0032] Iteratively update the gas flow velocity at the initial moment according to the updated gas pressure at the initial moment, and the formula is:

[0033]

[0034] where c represents the speed of sound.

[0035] Preferably, the boundary conditions of the natural gas dynamic model include: the gas pressure and mass flow rate at the inlet of the natural gas pipeline; the gas pressure and mass flow rate at the outlet of the natural gas pipeline; the gas pressure and mass flow rate at the connection of adjacent natural gas pipelines; the safety range of the gas pressure and mass flow rate inside the natural gas pipeline.

[0036] Preferably, the three - stage leap - frog finite - difference method is used to solve the natural gas dynamic model, including:

[0037] The friction coefficient and gas flow velocity are flow coefficients, and the gas pressure and mass flow rate are state variables; the solution region is discretized using a staggered grid to decouple the coupling relationship between the gas flow velocity, mass flow rate, and flow coefficient within the same time scale, including:

[0038] Obtain the spatial difference step h and the time difference step τ, and divide the solution region by two sets of parallel lines with a time interval of τ / 2 and a spatial interval of h / 2, where the spatial coordinates inside the natural gas pipeline are staggered - indexed by i and j = i + 1 / 2, and the time coordinates are staggered - indexed by n and m = n + 1 / 2;

[0039] Assume the state variables and flow coefficients at time level t n and time level t m are both known, then the solution process of the state variables and flow coefficients at time levels t n+1 and t m+1 includes:

[0040] Discretize the momentum equation using the central - difference derivative of the state variables, and solve for the mass flow rate at time level t m based on the gas pressure and flow coefficient at time level t n and the mass flow rate at time level t n+1 ;

[0041] Discretize the continuity equation using the central - difference derivative of the state variables, and solve for the mass flow rate at time level t n+1 based on the mass flow rate at time level t m and the gas pressure at time level tm+1 gas pressure;

[0042] Discretize the friction coefficient equation and the gas flow velocity equation at the node (x i+1 , t n+1 ) using the central difference derivatives of the friction coefficient, gas flow velocity, and gas pressure. Solve for the flow coefficient within the time layer t m+1 and time layer t m based on the gas pressure at time layer t n+1 , the mass flow rate at time layer t m , and the flow coefficient at time layer t m+1 .

[0043] Preferably, solve for the mass flow rate at time layer t m based on the gas pressure and flow coefficient at time layer t n and the mass flow rate at time layer t n+1 . The formula includes:

[0044] M i,n+1 = α i,m M i,n + β i,m (P j-1,m - P j,m ), 0 ≤ x i < L

[0045]

[0046] where M i,n+1 and M i,n represent the mass flow rates at time layer t n+1 and time layer t n at the spatial coordinate i, respectively. α i,m and β i,m are non-zero coefficients. P j-1,m and P j,m represent the gas pressures at time layer t m at the spatial coordinates j - 1 and j, respectively. x i represents the spatial coordinate at the i-th location, and L represents the total length of the solution domain; and represent the gas flow velocity and friction coefficient at time layer t m at the spatial coordinate i, respectively. D and A represent the pipe diameter and cross-sectional area; τ represents the time difference step size;

[0047] Solve for the gas pressure at time layer t n+1 based on the mass flow rate at time layer t m and the gas pressure at time layer t m+1 . The formula includes:

[0048] Pj,m+1 = P j,m + χ(M i,n+1 - M i+1,n+1 ), 0 ≤ x j < L

[0049]

[0050] where P j,m+1 and P j,m represent the gas pressures at time level t and time level t m+1 at spatial coordinate j respectively, χ represents a non - zero coefficient, M m and M i,n+1 represent the mass flows at time level t at spatial coordinate i and spatial coordinate i + 1 respectively i+1,n+1 , x n+1 represents the spatial coordinate at the j - th location, c represents the speed of sound; h represents the spatial difference step size; j Based on the gas pressures at time level t

[0051] and time level t m+1 and the mass flows at time level t m and the flow coefficient at time level t n+1 solve for the flow coefficient within time level t m , the formula includes: m+1 where

[0052]

[0053] where and represent the friction coefficients at time level t and time level t m+1 at spatial coordinate i + 1 respectively, M m represents the mass flow at time level t at spatial coordinate i + 1 i+1,n+1 , c1, c2 and c3 are empirical parameters, λ and μ represent the pipe roughness and the gas dynamic viscosity respectively; n+1 and and represent the gas flow velocities at time level t and time level t m+1 at spatial coordinate i + 1 respectively, P m and P j,m and P j+1,m represent the gas pressures at time level t at spatial coordinate j and spatial coordinate j + 1 respectively m , P j,m+1 and P j+1,m+1 represent the gas pressures at time level t at spatial coordinate j and spatial coordinate j + 1 respectively m+1 at time level t

[0054] Preferably, a multi-objective solution method is adopted to obtain the spatial difference step size h and the temporal difference step size τ. The formulas include:

[0055]

[0056] s.t. h min ≤ h ≤ h max

[0057] τ min ≤ τ ≤ τ max

[0058] ω1 + ω2 = 1

[0059]

[0060] where ω1 and ω2 represent weight coefficients, ζ(·) and respectively represent the complexity and error of the difference model, h min and h max respectively represent the lower and upper limits of the spatial difference step size, τ min and τ max respectively represent the lower and upper limits of the temporal difference step size; and represent the truncation error of the solution, and respectively represent the per-unit forms of the spatial step size and the temporal step size, and H and Γ respectively represent the reference spatial step size and the reference temporal step size.

[0061] Preferably, the global stability condition of the three-stage leapfrog finite difference method is:

[0062]

[0063] where h represents the spatial difference step size, τ represents the temporal difference step size; D represents the pipe diameter, c represents the speed of sound, and are respectively the gas flow velocity and the friction coefficient at the time layer t m at the spatial coordinate i.

[0064] The present invention also provides a gas network simulation system based on electro-gas dynamic analysis, including:

[0065] A model construction module for constructing a natural gas dynamic model based on variable coefficient partial differential equations, including: constructing a natural gas dynamic model based on the continuity equation, momentum equation, and state equation of one-dimensional hydrodynamics; constructing a friction coefficient equation based on the mass flow rate, and constructing a gas flow velocity equation based on the mass flow rate and gas pressure; expressing the friction coefficient and gas flow velocity in the natural gas dynamic model as the friction coefficient equation and the gas flow velocity equation to obtain the natural gas dynamic model;

[0066] An initialization module for calculating the gas pressure, mass flow rate, friction coefficient, and gas flow velocity at the initial moment;

[0067] A simulation solution module for solving the natural gas dynamic model based on the boundary conditions of the natural gas dynamic model to obtain the gas pressure, mass flow rate, friction coefficient, and gas flow velocity of natural gas in the pipeline at different moments.

[0068] The above technical solution of the present invention has the following beneficial effects compared with the prior art:

[0069] A gas network simulation method based on electro-gas dynamic analysis according to the present invention constructs a natural gas dynamic model based on a variable coefficient partial differential equation, considers the flow coefficient, that is, the dynamic changes of the friction coefficient and the gas flow velocity, avoids the loss of simulation accuracy caused by linearizing the fixed average gas velocity, realizes an accurate equivalent representation of the original partial differential equation, and improves the accuracy of the natural gas dynamic model. Furthermore, the present invention further proposes a three-stage leapfrog finite difference method to solve the natural gas dynamic model. By separating the difference units of the gas flow velocity, mass flow rate, and flow coefficient, direct solution of each discrete equation is realized, significantly improving the solution efficiency, and thus updating the flow coefficient and state variables in real time and efficiently during the difference process. By constructing and solving the natural gas dynamic model, the present invention can understand the gas pressure, mass flow rate, friction coefficient, and gas flow velocity of natural gas in the natural gas pipeline at different moments in real time, realizes accurate simulation of the gas network in the integrated power and natural gas system, can thus more reasonably allocate natural gas resources, improves energy utilization efficiency, reduces the operation cost of the gas network, and improves the reliability and economy of the dispatching of the integrated power and natural gas system.

[0070] Furthermore, the present invention analyzes the stability conditions of the three-stage leapfrog finite difference method and uses them as the constraint conditions for the spatial difference step size and the time difference step size to ensure the stability of the three-stage leapfrog finite difference method, thereby improving the accuracy and reliability of the gas network simulation results and ensuring the safe and stable operation of the gas network. Description of the Drawings

[0071] In order to make the content of the present invention easier to be clearly understood, the following further details the present invention according to the specific embodiments of the present invention in combination with the drawings, wherein:

[0072] Figure 1 is a flowchart of a gas network simulation method based on electro-gas dynamic analysis of the present invention;

[0073] Figure 2 is a comparison diagram of different energy flow models;

[0074] Figure 3 is a schematic diagram of the three-stage leapfrog finite difference method;

[0075] Figure 4 It is a schematic diagram for solving the three - stage leapfrog finite - difference method;

[0076] Figure 5 It is a diagram of the time scales and data exchange relationships between the power network and the natural gas network;

[0077] Figure 6 It is a flow chart for solving the optimal power flow model;

[0078] Figure 7 It is a schematic diagram of a natural gas pipeline;

[0079] Figure 8 It is a comparison diagram of the outlet P and inlet M curves of different models and the benchmark, where Figure 8 (a) in it represents the outlet P curve diagram, Figure 8 and (b) in it represents the inlet M curve diagram;

[0080] Figure 9 It is a diagram of the variation of the flow coefficient of each section of the pipeline with time, where Figure 9 (a) in it represents the fluctuation diagram of the gas flow velocity, Figure 9 and (b) in it represents the fluctuation diagram of the friction coefficient;

[0081] Figure 10 It is a diagram of the time - varying characteristics of the resistance term of each section of the pipeline, where Figure 10 and (a) in it represents the time - varying characteristic diagram of the resistance term, Figure 10 and (b) in it represents the time - varying characteristic diagram of the fluctuation of the resistance term. Detailed implementation mode

[0082] The present invention will be further described below in conjunction with the accompanying drawings and specific embodiments, so that those skilled in the art can better understand the present invention and be able to implement it, but the embodiments given do not limit the present invention.

[0083] Existing studies have explored natural gas dynamic models from various angles, but they are all based on gas flow models and constant - coefficient partial differential equations. Although they simplify the computational complexity, this inevitably affects the accuracy of the models, reduces the simulation accuracy of the gas network, and ultimately affects the operation of the integrated power and natural gas systems.

[0084] In addition, the current finite - difference method can only solve the constant - coefficient partial differential equations of existing natural gas dynamic models and is not applicable to the partial differential equations of gas flow with varying flow coefficients. Therefore, a large number of large - scale simultaneous equations need to be solved in the process of solving natural gas dynamic models, reducing the simulation efficiency of the gas network and making it difficult for natural gas dynamic models to be applied to system real - time scheduling.

[0085] The present invention provides a new method for constructing a natural gas dynamic model based on variable coefficient partial differential equations, and provides a three-stage leapfrog finite difference method for solving the natural gas dynamic model constructed based on variable coefficient partial differential equations.

[0086] Embodiment 1

[0087] Refer to Figure 1 As shown, the present invention proposes a gas network simulation method based on electro-gas dynamic analysis, including:

[0088] Step 1: Construct a natural gas dynamic model based on a variable coefficient partial differential equation (Variable coefficient Partial Differential Equation, VC-PDE), including:

[0089] Construct an initial natural gas dynamic model based on the continuity equation, momentum equation, and state equation of one-dimensional hydrodynamics; construct a friction coefficient equation based on the mass flow rate, and construct a gas velocity equation based on the mass flow rate and gas pressure; express the friction coefficient and gas velocity in the initial natural gas dynamic model as the friction coefficient equation and gas velocity equation to obtain the natural gas dynamic model;

[0090] Step 2: Calculate the gas pressure, mass flow rate, friction coefficient, and gas velocity at the initial moment;

[0091] Step 3: Determine the boundary conditions of the natural gas dynamic model;

[0092] Step 4: Solve the natural gas dynamic model using the three-stage leapfrog finite difference method (Three-stage leapfrog Finite Difference Method, TL-FDM) to obtain the gas pressure, mass flow rate, friction coefficient, and gas velocity of natural gas in the pipeline at different times.

[0093] In the said Step 1, constructing the natural gas dynamic model specifically includes the following steps:

[0094] Refer to the comparison chart of different energy flow models Figure 2 As shown. Usually, the state variables of gas flowing along the pipeline can be described by one-dimensional hydrodynamics, including three basic equations: the continuity equation, the momentum equation, and the state equation. However, even under one-dimensional conditions, the dynamic equations controlling gas flow are very complex, and past research has shown that ignoring the convective term in the momentum equation has the least impact on the model accuracy. Therefore, the initial natural gas dynamic model can be represented by a set of partial differential equations of mass flow rate, namely Formulas (1)-(3):

[0095]

[0096] P = c 2 ρ (3)

[0097] M = ρvA (4)

[0098] Wherein, P represents the gas pressure, M represents the mass flow rate, f represents the friction coefficient, v represents the gas flow velocity, x and t represent the spatial coordinate and the time coordinate respectively; D and A represent the diameter and the cross-sectional area of the gas pipeline respectively, c represents the speed of sound, and ρ represents the gas density. Among them, the friction coefficient and the gas flow velocity are flow coefficients, and the gas pressure and the mass flow rate are state variables.

[0099] Since the resistance term in formula (2) contains the square term of the flow velocity, it is non-linear. Due to the frequent changes of the gas flow velocity and the friction coefficient with the state variables, simply using the average gas velocity in the pipeline as a constant coefficient for substitution may cause the simulation results to be distorted. Therefore, the dynamic gas flow model described by (1)-(3) basically constitutes a set of VC-PDEs.

[0100] To more accurately simulate the dynamic gas flow, the real-time update of the flow coefficient value during the gas transmission process is crucial. The flow coefficient consists of the friction coefficient and the gas flow velocity , and their values depend on the fixed pipeline parameters and the time-varying state variables. The friction coefficient is a dimensionless quantity that is non-linearly negatively correlated with the Reynolds number Re, and it is usually described by the Hofer equation (5) applicable to the pipeline gas flow:

[0101]

[0102] Wherein, c1, c2 and c3 are empirical parameters, and λ and μ represent the pipeline roughness and the gas dynamic viscosity respectively.

[0103] Combining formulas (4), (5) and (6), a friction coefficient equation (7) is constructed based on the mass flow rate; combining formulas (3) and (4), a gas flow velocity equation (8) is constructed based on the mass flow rate and the gas pressure; then substituting formulas (4), (7) and (8) into formula (2), the equivalent VC-PDE form (9) of the momentum equation can be obtained as follows:

[0104]

[0105] Wherein, represents the friction coefficient equation, represents the gas flow velocity equation.

[0106] Therefore, the natural gas dynamic model can be represented by formulas (1), (9) and (3).

[0107] However, the variable flow coefficient in Equation (9) and are determined by the unknowns M and P, making it difficult to solve the VC-PDE using existing finite difference methods. To address this issue, the present invention proposes a new TL-FDM algorithm that can update the state variables and flow coefficients simultaneously.

[0108] In step 4, solving the natural gas dynamic model using the three-stage leapfrog finite difference method specifically includes the following steps:

[0109] Referring to Figure 3 as shown, the solution domain is discretized using a staggered grid to decouple the coupling relationship between gas velocity, mass flow rate, and flow coefficient within the same time scale, including:

[0110] Obtain the spatial difference step size h and the temporal difference step size τ, and divide the solution domain by two sets of parallel lines with a time interval of τ / 2 and a spatial interval of h / 2, where the spatial coordinates inside the natural gas pipeline are staggered-indexed by i and j = i + 1 / 2, and the time coordinates are staggered-indexed by n and m = n + 1 / 2.

[0111] Construct three sets of node clusters at the intersection points of the grid lines. Referring to Equation (10), and represent the mass flow rate, gas pressure, and flow coefficient respectively:

[0112]

[0113] At this time, the mass flow rate, gas pressure, and flow coefficient are discretized and stored in independent difference cells, so it is allowed to solve the unknown variables using only known quantities.

[0114] Referring to Figure 4 as shown, the proposed TL-FDM includes three stages. Assume that the state variables and flow coefficients at time level t n and time level t m are both known. Then, the solution process for the state variables and flow coefficients at time levels t n+1 and t m+1 includes the following steps:

[0115] Stage 1: Discretize the momentum equation using the central difference derivatives (11)-(12) of the state variables, and solve for all the mass flow rates at time level t m based on the gas pressure and flow coefficient at time level t n and the mass flow rate at time level t n+1 :

[0116]

[0117] Mi,n+1 = α i,m M i,n + β i,m (P j-1,m - P j,m ), 0 ≤ x i < L (14)

[0118]

[0119] where M i,n+1 and M i,n represent the mass flow rates at time level t at spatial coordinate i n+1 and time level t n respectively, α i,m and β i,m are non - zero coefficients, P j-1,m and P j,m represent the gas pressures at time level t at spatial coordinates j - 1 and j respectively m x i represents the spatial coordinate at the i - th location, and L represents the total length of the solution domain, i.e., the total length of the pipeline; and represent the gas flow velocity and the friction coefficient at time level t at spatial coordinate i m respectively, and τ represents the time difference step size.

[0120] Stage 2: Discretize the continuity equation using the central - difference derivatives (17)-(18) of the state variables, and solve for all the gas pressures within time level t n+1 based on the mass flow rate at time level t m and the gas pressure at time level t m+1 :

[0121]

[0122] P j,m+1 = P j,m + X(M i,n+1 - M i+1,n+1 ), 0 ≤ x j < L (19)

[0123]

[0124] where P j,m+1 and P j,m represent the gas pressures at time level t at spatial coordinate j m+1 and time level t m respectively, X represents a non - zero coefficient, M i,n+1 and M i+1,n+1 represent the mass flow rates at time level t at spatial coordinates i and i + 1 respectively n+1The mass flow rate, obtained from Stage 1; x j represents the spatial coordinate at the j-th location; h represents the spatial difference step size.

[0125] Stage 3: Discretize the friction coefficient equation and the gas velocity equation at the node (x i+1 , t n+1 ) using the central difference derivatives (21)-(23) of the friction coefficient, gas velocity, and gas pressure. Solve for the flow coefficient within the time layer t m+1 and the time layer t m based on the gas pressure at the time layer t n+1 , the mass flow rate at the time layer t m , and the flow coefficient at the time layer t m+1 :

[0126]

[0127]

[0128] where and represent the friction coefficients at the spatial coordinate i + 1 for the time layers t m+1 and t m respectively, and M i+1,n+1 represents the mass flow rate at the spatial coordinate i + 1 for the time layer t n+1 ; and represent the gas velocities at the spatial coordinate i + 1 for the time layers t m+1 and t m respectively, and P j,m and P j+1,m represent the gas pressures at the spatial coordinates j and j + 1 for the time layer t m respectively; P j,m+1 and P j+1,m+1 represent the gas pressures at the spatial coordinates j and j + 1 for the time layer t m+1 , obtained from Stage 2.

[0129] Therefore, the TL-FDM of the natural gas dynamic model proposed by the present invention can be expressed by formulas (14), (19), (24), and (25), where each discrete equation contains only one unknown variable, allowing for direct solution and having significant advantages in terms of computational efficiency. In addition, the flow coefficient can be calculated and stored based on the values of the state variables changing within each differential cell, providing a more realistic flow coefficient value for the solution of the state variables, thereby greatly improving the accuracy of the simulation results.

[0130] To obtain an accurate solution for the dynamic energy flow of natural gas in a pipeline, boundary conditions must be specified first. In step 3, determining the boundary conditions of the natural gas dynamic model specifically includes:

[0131] Assume that the pipeline inlet P 0,t (x = 0) and the outlet M L,t (x = L) are known in formulas (26) and (27) respectively. Substituting them into formula (19) can obtain the boundary conditions (28)-(29) at the pipeline inlet M 0,t and the outlet P L,t . The formulas are as follows:

[0132] P 0,t = P Inlet (t), 0 ≤ t ≤ T (26)

[0133] M L,t = M Outlet (t), 0 ≤ t ≤ T (27)

[0134]

[0135] P L,m+1 = P L,m + X(M L-1,n+1 - M Outlet (t n+1 )) (29)

[0136] Among them, P 0,t represents the gas pressure at the pipeline inlet, and P Inlet (t) represents the function of the gas pressure change at the pipeline inlet; M L,t represents the mass flow rate at the pipeline outlet, and M Outlet (t) represents the function of the mass flow rate change of the gas at the pipeline outlet; M 0,n+1 and M 1,n+1 respectively represent the mass flow rates at the pipeline inlet and at the spatial coordinate 1 at time layer t n+1 ; P L,m+1 and P L,m respectively represent the gas pressures at the pipeline outlet at time layer t m ; M L-1,n+1 represents the mass flow rate at the spatial coordinate L - 1 at time layer t n+1 ; L and T respectively represent the length and time of the solution region.

[0137] Therefore, the boundary conditions for the gas pressure and mass flow rate at the natural gas pipeline inlet are formulas (26) and (28); the boundary conditions for the gas pressure and mass flow rate at the natural gas pipeline outlet are formulas (27) and (29).

[0138] The boundary pressures at the connections between adjacent pipes should be equal, and according to the law of conservation of mass, the mass fluxes flowing in and out at the junction should be consistent, as expressed by the following formulas:

[0139] P Pipe-1,m = P Pipe-2,m = P Pipe-3,m … (30)

[0140] M Pipe-1,n + M Pipe-2,n + M Pipe-3,n … = 0 (31)

[0141] where P Pipe-1,m , P Pipe-2,m and P Pipe-3,m represent the gas pressures of the first pipe, the second pipe, and the third pipe at time level t m respectively; M Pipe-1,n , M Pipe-2,n and M Pipe-3,n represent the mass fluxes of the first pipe, the second pipe, and the third pipe at time level t n respectively.

[0142] Therefore, the boundary conditions for the gas pressure and mass flux at the connections between adjacent natural gas pipes are given by formulas (30) and (31).

[0143] Since the boundary calculations for the mass flux and gas pressure are independent of the flow coefficient, the boundary conditions for the flow coefficient are not derived.

[0144] In practical applications, the safety ranges and constraint conditions for the gas pressure and mass flux inside the natural gas pipe can be expressed by formula (32):

[0145]

[0146] where P j,m represents the gas pressure at time level t m at spatial coordinate j, and M i,n represents the mass flux at time level t n at spatial coordinate i; P j and represent the minimum and maximum gas pressure limits respectively, M i and represent the minimum and maximum mass flux limits respectively.

[0147] In step 2, the process of calculating the initial moment conditions includes:

[0148] Assume that the system is stable in the initial period. Then, the partial derivatives of the mass flow rate and gas pressure with respect to time satisfy equations (33)-(34). Combining equations (1) and (11), the distribution of the mass flow rate at the initial moment can be obtained, namely equations (35) and (36). The latter indicates that the mass flow rate in the pipeline is equal everywhere during the first two time periods.

[0149]

[0150] M i+1,0 =M i,0 (35)

[0151] M i,n =M Outlet (t0), n ≤ 2 (36)

[0152] where M i+1,0 and M i,0 represent the mass flow rates at the initial moment at spatial coordinates i + 1 and i, respectively; M outlet (t0) represents the value of the mass flow rate change function at the pipeline outlet at the initial moment.

[0153] Therefore, the mass flow rate at the initial moment is the value of the mass flow rate change function at the pipeline outlet at the initial moment.

[0154] Combining equations (9), (12), (13), (34), and (36), the gas pressure distribution at the initial moment can be obtained:

[0155] P j,0 =P j-1,0 + ε i,0 M Outlet (t0) (37)

[0156]

[0157] where P j,0 and P j-1,0 represent the gas pressures at the initial moment at spatial coordinates j and j - 1, respectively, ε i,0 is a non-zero coefficient, represents the gas flow velocity at the initial moment at spatial coordinate i, represents the friction coefficient at the initial moment at spatial coordinate i. The value of P j,0 is defined as the gas well pressure in the natural gas system.

[0158] Furthermore, by combining equations (7) and (8) with equation (36) respectively, the distributions of the friction coefficient and gas flow velocity in the initial period can be calculated:

[0159]

[0160] Since the state variables and the flow coefficient are coupled with each other in the initial period, for the problem of initializing the gas flow rate, the process of calculating the gas pressure, mass flow rate, friction coefficient and gas flow velocity at the initial moment in the present invention specifically includes:

[0161] Initializing the values of the gas flow velocity and gas pressure at the initial moment; the mass flow rate at the initial moment is the function value of the mass flow rate change at the pipeline outlet at the initial moment; and calculating the friction coefficient at the initial moment according to the function value of the mass flow rate change at the pipeline outlet at the initial moment;

[0162] Iteratively updating the gas pressure at the initial moment according to the function value of the mass flow rate change at the pipeline outlet at the initial moment, the gas flow velocity at the initial moment and the friction coefficient, and then iteratively updating the gas flow velocity at the initial moment according to the updated gas pressure at the initial moment until the maximum number of iterations is reached or the change amounts of the initial gas pressure and the initial gas flow velocity are less than the preset threshold.

[0163] The pseudo-code of the gas flow rate initialization algorithm is as follows:

[0164] 1. Input: Convergence tolerance δ and maximum number of iterations K;

[0165] 2. Initialization: Number of iterations k ← 1. Initialize the values of and; set the values of P(t) and M(t), t = 0; and ; set P Inlet (t) and M Outlet (t), t = 0;

[0166] 3. Obtain by solving formula (39)

[0167] 4. For the fixed coefficients and Obtain by solving the coefficient equation (38) Obtain by solving the initial condition (37)

[0168] 5. Obtain by solving the coefficient equation (40)

[0169] 6. Let the iteration deviation

[0170] 7. Make

[0171] 8. k ← k + 1

[0172] 9. Return to step 4 for iteration until Δσ ≤ δ or k > K + 1.

[0173] The algorithm is executed at most K times, or until Δσ is less than δ. Usually, the number of iterations for the algorithm to converge does not exceed four, and the calculation time is less than one second.

[0174] Furthermore, the spatial difference step size h and the temporal difference step size τ must satisfy certain conditions to ensure the stability of the TL-FDM proposed in the present invention. Therefore, the stability of the TL-FDM is analyzed in this embodiment.

[0175] Using the frozen coefficient method, and are fixed as constants for the stability analysis. By using the discrete Fourier transform, equations (14) and (19) are re-expressed as equation (41):

[0176]

[0177] where, and represent the frequency functions of the mass flow rate and the gas pressure respectively, represents the amplification factor matrix, which can be expressed as equation (42):

[0178]

[0179] where, i ′ is the imaginary unit, α and β represent the coefficients in the algebraic equations obtained from the finite differences of the natural gas mass flow rate, and ξ represents the unit in the frequency domain.

[0180] According to the von Neumann stability condition, the spectral radius of must be less than or equal to 1 to ensure the stability of the difference scheme. According to Euler's formula, its eigenvalue γ is:

[0181]

[0182] Substituting equations (16) and (20) into equation (43), the eigenvalue can be re-expressed as equation (44):

[0183]

[0184] If then γ is a real number. In this case, when ξ = ±π, the right side of the above inequality takes the maximum value. Therefore, it is required to always maintain The stability condition max|γ| ≤ 1 is equal to:

[0185]

[0186] Since and All are greater than 0. Therefore, formula (45) always holds for any value of ξ. Thus, the TL-FDM proposed by the present invention is unconditionally stable under the following circumstances:

[0187]

[0188] If then γ is a complex number, and the stability condition max|γ|≤1 is equal to:

[0189]

[0190] The root terms in the above formula are imaginary numbers. Therefore, formula (47) can be reformulated as:

[0191]

[0192] All terms in the above formula are positive. Therefore, when ξ = ±π, the right-hand side term reaches its maximum value. Thus, the stability condition under this condition can be reformulated as:

[0193]

[0194] Because the right-hand side term of formula (46) is greater than the right-hand side term of formula (49). The stability condition of the TL-FDM proposed by the present invention is

[0195] This condition only ensures the local stability at the coefficient freezing point. To ensure global stability, the set of all local stability conditions should be considered. Therefore, the global stability condition of the TL-FDM proposed by the present invention is:

[0196]

[0197] Where and are the gas flow velocity and friction coefficient at time layer t m at spatial coordinate i, respectively.

[0198] Example Two

[0199] Based on Example One, an application example of a gas network simulation method based on electro-gas dynamic analysis in the dynamic optimal power flow of an integrated power and natural gas system is provided.

[0200] I. Limitations of the Power Network

[0201] Since the power flow in the power grid transmits at the speed of light and the dynamic time scale is in the millisecond level, the dynamic process of the EPN can be ignored. Its operating constraints include: node active power balance constraint, generation capacity constraint, generation ramp constraint, and transmission line active power flow limit, expressed as formulas (51)-(54):

[0202]

[0203] Among them, and represent the minimum and maximum power generations respectively, and S l represent the ramping capacity and the thermal limit respectively, and represent the active powers of the generator set and the load respectively; and represent the active powers of the kth coal-fired unit, gas-fired unit, and renewable energy unit in the time period t e respectively; and represent the sets of coal-fired units, gas-fired units, and renewable energy units connected to the bus i respectively; represents the set of power system buses; represents the active power through the power line, and represent the phase angles at the start and end points of the power line respectively, b k represents the admittance of the power line k, t e represents the time step of the power network.

[0204] Since the natural gas network uses a slower time step to describe its dynamic process, the time steps between the EPN and the NGN are mismatched. Therefore, this embodiment adopts a multi-rate method. As Figure 5 shown, whenever the EPN simulates one step, the NGN completes N = Γ / τ steps in parallel, and exchanges data between the two networks, including the operation information and safety status of the gas-fired units. Therefore, the gas consumption of the GFU in the natural gas network is a constant value within each time step of the power grid, which can be expressed as:

[0205]

[0206] Among them, η represents the energy conversion efficiency, represents the time index t of the EPN e corresponding to the gas consumption of the kth GFU at the nth time index of the NGN network, N represents the multiple of the length of a power system time period to the length of a natural gas system time period, and the time index t e of the power system corresponds to the time index of the natural gas system as t g = Nt e .

[0207] II. Selection of the differential step size

[0208] Increasing the density of the difference step size can improve the simulation accuracy, but at the same time it increases the dimension of the constraint set, resulting in an increase in the computational amount. Therefore, the selection of the difference step size needs to consider the balance relationship between the model accuracy and the computational efficiency. In this embodiment, a multi-objective solution method, that is, formula (56), is used to select the optimal combination of the difference step size, and the size of the constraint set is used to quantify the computational complexity, while the error of the solution is represented by the truncation error, that is, formula (57):

[0209]

[0210] Among them, ω1 and ω2 represent the weight coefficients, ζ(·) and respectively represent the complexity and error of the difference model, h min and h max respectively represent the lower and upper limits of the spatial difference step size, τ min and τ max respectively represent the lower and upper limits of the time difference step size; and represent the truncation error of the solution, and respectively represent the per-unit forms of the spatial step size and the time step size, H and Γ respectively represent the reference spatial step size and the reference time step size.

[0211] Among them, the selection of the weight coefficient depends on the scale of the system. Large IEGSs usually give priority to high computational efficiency, while smaller IEGSs tend to focus on the accuracy of the results.

[0212] III. Gas flow initialization

[0213] This step is initialized using the gas flow initialization algorithm in Embodiment 1.

[0214] IV. Optimal power flow model

[0215] The objective of the optimal power flow model is to minimize the total fuel cost of the IEGS. This cost includes the power generation costs of coal-fired units and gas-fired units. The cost coefficient of the renewable energy unit is set to 0 to maximize its output, and the cost coefficient of the coal-fired unit is set higher than that of the gas turbine to reduce fossil fuel consumption. Generally, it can be represented by formula (59):

[0216]

[0217] Among them, c c and c g respectively represent the cost coefficients of the coal-fired unit and the gas well, t g represents the time index (in seconds) of the natural gas system, t e represents the time index (in hours) of the power system, Denote the mass flow rate of natural gas supplied by the k-th natural gas unit in period t (kg / s). Multiply it by the time step τ of this period, which is the total amount of natural gas supplied by the k-th natural gas unit in period t. g As shown in g it is the total amount of natural gas supplied by the k-th natural gas unit in period t.

[0218] Refer to Figure 6 As shown, the proposed optimal power flow model includes three stages:

[0219] In the first stage, based on the requirements of computational efficiency and simulation accuracy, obtain the optimal combination of the spatial difference step h and the time difference step τ under the global stability condition of TL-FDM, and incorporate the difference stability constraints of h and τ into the constraint set. Since the value of the flow coefficient is currently unknown, first use the reference value of historical data to replace in the global stability condition and verify it in the third stage; use the gas flow initialization algorithm to solve the initial values of gas pressure and gas velocity in the natural gas network;

[0220] In the second stage, solve the optimal power flow model of IEGS according to the results of the first stage;

[0221] In the third stage, verify the feasibility of the reference value of the flow coefficient. If then it meets the stability condition of the model, and output the current system operation results; otherwise, the model may not converge, resulting in unreliable results. Then use of the current system operation results to update and return to the first stage.

[0222] Embodiment 3

[0223] This embodiment will take a single natural gas pipeline as an example to verify the effectiveness of the above high-precision gas network model. The pipeline schematic diagram is as shown in Figure 7 The pipeline inlet is connected to a 55 bar constant pressure gas well, and the outlet is connected to a gas load that changes with time. This part simulates three models to analyze the influence of different flow coefficients on the gas flow distribution, including the WEs static model (S1), the traditional CC-PDE dynamic model (S2), and the natural gas dynamic model (S3) based on variable coefficient partial differential equations proposed by the present invention.

[0224] Analyze the dynamic gas flow in the pipeline. The jump at the outlet M causes changes in the gas flow in the pipeline. The curves of the outlet P and the inlet M of each model are as shown in Figure 8 where Figure 8 (a) represents the curve of the outlet P, Figure 8In (b), it represents the inlet M curve graph. Due to the slow dynamic characteristics of the natural gas in the pipeline, the diffusion process needs to be extended to reach a new steady state. However, S1 ignores the time delay of gas transmission and assumes that both the outlet P and the inlet M can reach the steady state instantaneously, which contradicts the physical laws. Therefore, S1 can only calculate the state variables of gas flow at the steady state and cannot simulate the actual gas flow process, while the dynamic models S2 and S3 can describe the physical process of gas diffusion. The natural gas produced by the Wells gas well cannot be immediately transported to the load. In fact, the natural gas demand is supported by the gas stored in the pipeline package. Therefore, the pipeline outlet P gradually decreases as the gas storage in the pipeline package is consumed and reaches a new steady state after the gas supplemented by the Wells meets the gas load. However, frequent load fluctuations may affect the gas network to reach the steady state. Therefore, the accuracy of the dynamic gas flow model is crucial for ensuring the reliability of gas network operation.

[0225] Flow coefficient during gas diffusion and fluctuations are as Figure 9 shown, where Figure 9 in (a) represents the fluctuation graph of gas flow velocity, Figure 9 in (b) represents the fluctuation graph of friction coefficient. Relative to its initial value and the maximum fluctuations are 70.48% and 5.59% respectively. The gray curves represent the and values of different sections of the pipeline, and the length of each section is the value of the spatial step. Therefore, and fluctuations lead to a significant change in the resistance term in the momentum equation, which in turn has a considerable impact on the state variables in the pipeline.

[0226] Table 1 gives the simulation results of the flow coefficient, resistance term and state variables at the end of the pipeline. The resistance term in the pipeline increases with the increase of M, which will lead to an increase in the pressure drop in the pipeline, and vice versa. In addition, the difference between the resistance terms at the inlet and outlet of the pipeline also increases with M, as Figure 10 shown, where Figure 10 in (a) represents the time-varying characteristic graph of the resistance term, Figure 10 in (b) represents the time-varying characteristic graph of the fluctuation of the resistance term. This difference will further affect the gas flow distribution in the pipeline. The S2 model ignores the increase of the resistance term, and the outlet P is 4.63 bar higher than the reference, with an error rate of 13.53%; the S3 model considers the real-time change of the flow coefficient, and the simulation results are in good agreement with the reference values. The relative error between each model and the reference is shown in Table 2. The solution accuracy of S3 is more than 33 times that of S2, the calculation accuracy of S3 is significantly higher than that of S2, and due to avoiding solving complex simultaneous equations, the calculation speed of S3 is more than 9 times faster than that of S2.

[0227] Table 1. Outlet values of flow coefficient, resistance term and state variable at different times

[0228]

[0229] Table 2. Comparison of simulation accuracies MRE and LRE of S1, S2 and S3

[0230]

[0231] Therefore, the simulation results of this embodiment show that the high-precision gas network model proposed by the present invention has good performance in terms of simulation accuracy and solution efficiency. At the same time, the model can be further applied to the solution of the integrated electrical energy system to achieve accurate perception of the natural gas system and improve the overall simulation accuracy and solution efficiency.

[0232] Embodiment 4

[0233] Based on the gas network simulation method based on electro-gas dynamic analysis described in Embodiment 1, this embodiment further provides a gas network simulation system based on electro-gas dynamic analysis, including:

[0234] A model construction module for constructing a natural gas dynamic model based on variable coefficient partial differential equations, including: constructing a natural gas dynamic model based on the continuity equation, momentum equation and state equation of one-dimensional fluid dynamics; constructing a friction coefficient equation based on mass flow, and constructing a gas velocity equation based on mass flow and gas pressure; expressing the friction coefficient and gas velocity in the natural gas dynamic model as the friction coefficient equation and the gas velocity equation to obtain the natural gas dynamic model;

[0235] An initialization module for calculating the gas pressure, mass flow, friction coefficient and gas velocity at the initial moment;

[0236] A simulation solution module for solving the natural gas dynamic model based on the boundary conditions of the natural gas dynamic model to obtain the gas pressure, mass flow, friction coefficient and gas velocity of natural gas in the pipeline at different times.

[0237] In summary, for the gas network simulation method based on electro-gas dynamic analysis described in the present invention, a natural gas dynamic model is constructed based on variable coefficient partial differential equations, taking into account the dynamic changes of the flow coefficient, i.e., the friction coefficient and the gas flow velocity, avoiding the loss of simulation accuracy caused by linearizing the fixed average gas velocity, achieving an accurate equivalent representation of the original partial differential equation, and improving the accuracy of the natural gas dynamic model. Moreover, the present invention further proposes a three-stage leapfrog finite difference method to solve the natural gas dynamic model. By separating the difference units of the gas flow velocity, mass flow rate, and flow coefficient, the direct solution of each discrete equation is realized, significantly improving the solution efficiency, and thus updating the flow coefficient and state variables in real time and efficiently during the difference process. By constructing and solving the natural gas dynamic model, the present invention can understand the gas pressure, mass flow rate, friction coefficient, and gas flow velocity of natural gas in the natural gas pipeline at different times in real time, achieving an accurate simulation of the natural gas network in the integrated power and natural gas system, thereby enabling a more reasonable allocation of natural gas resources, improving the energy utilization efficiency, reducing the operating cost of the gas network, and improving the reliability and economy of the integrated power and natural gas system scheduling.

[0238] Furthermore, the present invention analyzes the stability conditions of the three-stage leapfrog finite difference method and uses them as the constraint conditions for the spatial difference step size and the time difference step size to ensure the stability of the three-stage leapfrog finite difference method, thereby improving the accuracy and reliability of the natural gas network simulation results and ensuring the safe and stable operation of the gas network.

[0239] Obviously, the above embodiments are merely examples given for clear illustration and are not limitations on the implementation manners. For those of ordinary skill in the art, other different forms of changes or modifications can be made based on the above description. It is not necessary and impossible to list all the implementation manners here. And the obvious changes or modifications derived therefrom are still within the protection scope of the present invention.

Claims

1. A gas network simulation method based on electric-gas dynamic analysis, characterized in that: include: The natural gas dynamic model is constructed based on variable coefficient partial differential equations, including: An initial natural gas dynamic model is constructed based on the continuity equation, momentum equation and state equation of one-dimensional fluid dynamics; A friction coefficient equation is constructed based on mass flow, and a gas flow rate equation is constructed based on mass flow and gas pressure; the friction coefficient and gas flow rate in the initial natural gas dynamic model are expressed as a friction coefficient equation and a gas flow rate equation to obtain a natural gas dynamic model; Calculate the gas pressure, mass flow rate, friction coefficient and gas flow rate at the initial moment; The natural gas dynamic model is solved based on its boundary conditions to obtain the gas pressure, mass flow rate, friction coefficient and gas flow rate of the natural gas in the pipeline at different times.

2. A gas network simulation method based on electric-gas dynamic analysis according to claim 1, characterized in that: The friction coefficient equation is constructed based on the mass flow rate, and the formula is: in, represents the friction coefficient equation, D and A represent the pipe diameter and cross-sectional area respectively, c1, c2 and c3 are empirical parameters, M represents the mass flow rate, λ and μ represent the pipe roughness and gas dynamic viscosity respectively; The gas flow rate equation is constructed based on the mass flow rate and gas pressure. The formula is: in, represents the gas flow rate equation, P represents the gas pressure, and c represents the speed of sound; The formula of natural gas dynamic model is expressed as: P6c 2 ρ Among them, x and t represent the spatial coordinate and time coordinate respectively, and ρ represents the gas density.

3. A gas network simulation method based on electric-gas dynamic analysis according to claim 1, characterized in that: Calculate the gas pressure, mass flow rate, friction coefficient and gas flow rate at the initial moment, including: Initialize the values ​​of the gas flow rate and gas pressure at the initial moment; the mass flow at the initial moment is the function value of the mass flow change at the pipeline outlet at the initial moment; and calculate the friction coefficient at the initial moment according to the function value of the mass flow change at the pipeline outlet at the initial moment; The gas pressure at the initial moment is iteratively updated according to the function value of the pipeline outlet mass flow rate change at the initial moment, the gas flow rate at the initial moment and the friction coefficient, and then the gas flow rate at the initial moment is iteratively updated according to the updated gas pressure at the initial moment, until the maximum number of iterations is reached or the changes in the initial gas pressure and the initial gas flow rate are less than a preset threshold.

4. A gas network simulation method based on electric-gas dynamic analysis according to claim 3, characterized in that: The friction coefficient at the initial moment is calculated according to the function value of the mass flow rate change at the pipeline outlet at the initial moment. The formula is: in, represents the friction coefficient at the initial moment at the spatial coordinate i, D and A represent the diameter and cross-sectional area of ​​the gas pipeline respectively, c1, c2 and c3 are empirical parameters, λ and μ represent the pipeline roughness and gas dynamic viscosity respectively; M outlet (t0) represents the value of the pipeline outlet mass flow rate change function at the initial moment; The gas pressure at the initial moment is iteratively updated according to the value of the mass flow rate change function at the pipeline outlet at the initial moment, the gas flow rate at the initial moment, and the friction coefficient. The formula is: P j,0 =P j-1,0 +ε i,0 M Outlet (t0) Among them, P j,0 and P j-1,0 Respectively represent the gas pressure at the initial moment at the space coordinate j and the space coordinate j-1, ε i,0 is a non-zero coefficient, represents the gas flow rate at the initial moment at the spatial coordinate i, and h represents the spatial difference step size; The gas flow rate at the initial moment is iteratively updated according to the updated gas pressure at the initial moment. The formula is: Here, c represents the speed of sound.

5. A gas network simulation method based on electric-gas dynamic analysis according to claim 1, characterized in that: The boundary conditions of the natural gas dynamic model include: gas pressure and mass flow at the inlet of the natural gas pipeline; gas pressure and mass flow at the outlet of the natural gas pipeline; gas pressure and mass flow at the connection of adjacent natural gas pipelines; and the safe range of gas pressure and mass flow inside the natural gas pipeline.

6. A gas network simulation method based on electric-gas dynamic analysis according to claim 1, characterized in that: The three-stage frog leaping finite difference method is used to solve the natural gas dynamic model, including: The friction coefficient and gas flow rate are flow coefficients, and the gas pressure and mass flow rate are state variables. The solution area is discretized using staggered grids to decouple the coupling relationship between gas flow rate, mass flow rate and flow coefficient within the same time scale, including: Obtain the spatial difference step size h and the time difference step size τ, and divide the solution area by two sets of parallel straight lines with a time interval of τ / 2 and a spatial interval of h / 2, where the spatial coordinates inside the natural gas pipeline are indexed by i and j=i+1 / 2, and the time coordinates are indexed by n and m=n+1 / 2; Assume that the time layer t n and time layer t m The state variables and flow coefficient of are known, then the time layer t n+1 and t m+1 The process of solving the state variables and flow coefficient in includes: The momentum equation is discretized using the central difference derivative of the state variable, based on the time layer t m Gas pressure and flow coefficient and time layer t n The mass flow rate is solved for the time layer t n+1 Mass flow rate; The continuity equation is discretized using the central difference derivative of the state variable, based on the time layer t n+1 The mass flow rate and time layer t m The gas pressure is solved for the time layer t m+1 Gas pressure; The friction coefficient, gas velocity and gas pressure are used to calculate the central differential derivative of the node (x i+1 , t n+1 ) is discretized into the friction coefficient equation and the gas velocity equation, and based on the time layer t m+1 and time layer t m Gas pressure, time layer t n+1 The mass flow rate and time layer t m The flow coefficient of the time layer t is solved m+1 The flow coefficient inside.

7. A gas network simulation method based on electric-gas dynamic analysis according to claim 6, characterized in that: Based on the time layer t m Gas pressure and flow coefficient and time layer t n The mass flow rate is solved for the time layer t n+1 The mass flow rate is given by include: M i,n+1 =a i,m M i,n +b i,m (P j-1,m -P j,m ),0≤x i <L Among them, M i,n+1 and M i,n Respectively represent the time layer t at the spatial coordinate i n+1 and time layer t n The mass flow rate, α i,m and β i,m is a non-zero coefficient, P j-1,m and P j,m Respectively represent the time layer t at the spatial coordinate j-1 and the spatial coordinate j m Gas pressure, x i represents the spatial coordinate of the i-th location, and L represents the total length of the solution area; and Respectively represent the time layer t at the spatial coordinate i m The gas flow rate and friction coefficient, D and A represent the pipe diameter and cross-sectional area respectively; τ represents the time difference step; Based on the time layer t n+1 The mass flow rate and time layer t m The gas pressure is solved for the time layer t m+1 The gas pressure is given by: P j,m+1 =P j,m +χ(M i,n+1 -M i+1,n+1 ),0≤x j <L Among them, P j,m+1 and P j,m Respectively represent the time layer t at the spatial coordinate j m+1 and time layer t m The gas pressure, χ represents a non-zero coefficient, M i,n+1 and M i+1,n+1 Respectively represent the time layer t at spatial coordinate i and spatial coordinate i+1 n+1 The mass flow rate, x j represents the spatial coordinate at the jth location, c represents the speed of sound; h represents the spatial difference step length; Based on the time layer t m+1 and time layer t m Gas pressure, time layer t n+1 The mass flow rate and time layer t m The flow coefficient of the time layer t is solved m+1 The flow coefficient in the formula includes: in, and Respectively represent the time layer t at the spatial coordinate i+1 m+1 and time layer t m The friction coefficient, M i+1,n+1 Represents the time layer t at the spatial coordinate i+1 n+1 The mass flow rate, c1, c2 and c3 are empirical parameters, λ and μ represent the pipeline roughness and gas dynamic viscosity respectively; and Respectively represent the time layer t at the spatial coordinate i+1 m+1 and time layer t m The gas flow rate, P j,m and P j+1,m Respectively represent the time layer t at the spatial coordinate j and the spatial coordinate j+1 m The gas pressure, P j,m+1 and P j+1,m+1 Respectively represent the time layer t at the spatial coordinate j and the spatial coordinate j+1 m+1 of gas pressure.

8. A gas network simulation method based on electric-gas dynamic analysis according to claim 6, characterized in that: In the integrated power and natural gas system, a multi-objective solution method is used to obtain the spatial difference step h and the time difference step τ. The formulas include: s.t.h min ≤h≤h max t min ≤τ≤τ max ω1+ω2=1 Among them, ω1 and ω2 represent weight coefficients, ζ(·) and Respectively represent the complexity and error of the differential model, h min and h max Respectively represent the lower and upper limits of the spatial difference step size, τ min and τ max Respectively represent the lower and upper limits of the time difference step; and represents the truncation error of the solution, and denote the per-unit form of the spatial step and time step, respectively, H and Γ denote the reference spatial step and reference time step, respectively.

9. A gas network simulation method based on electric-gas dynamic analysis according to claim 6, characterized in that: The global stability condition of the three-stage frog leaping finite difference method is: Where h represents the spatial difference step, τ represents the time difference step; D represents the pipe diameter, c represents the sound speed, and are respectively the time layer t at spatial coordinate i m Gas flow rate and friction coefficient.

10. A gas grid simulation system based on electric-gas dynamic analysis, characterized in that: include: Model building module, used to build natural gas dynamic model based on variable coefficient partial differential equations, including: Construct a natural gas dynamic model based on the continuity equation, momentum equation and state equation of one-dimensional fluid dynamics; A friction coefficient equation is constructed based on mass flow, and a gas flow rate equation is constructed based on mass flow and gas pressure; the friction coefficient and gas flow rate in the natural gas dynamic model are expressed as a friction coefficient equation and a gas flow rate equation to obtain a natural gas dynamic model; Initialization module, used to calculate the gas pressure, mass flow, friction coefficient and gas at the initial moment Flow rate; The simulation solution module is used to solve the natural gas dynamic model based on the boundary conditions of the natural gas dynamic model. The gas pressure, mass flow rate, friction coefficient and gas flow rate of natural gas in the pipeline at different times are obtained.

Images