Electrical integrated energy system fault analysis method based on semi-implicit simulation

Abstract

The invention discloses an electrical integrated energy system fault analysis method based on semi-implicit simulation, and belongs to the field of integrated energy system simulation analysis. Specifically, in order to quantify the propagation velocity of the pipe explosion and leakage faults of the natural gas system and analyze the propagation mechanism of the faults in an electrical bidirectional coupling system, the method comprises the following steps: firstly, establishing a natural gas system differential algebraic equation model for describing the pipe explosion and leakage faults based on a characteristic line method; then, in consideration of the high-dimensional, nonlinear and rigid characteristics of the electrical bidirectional coupling system, a semi-implicit dynamic simulation method based on a rigid precise Rosenbrock format is provided; and finally, a key event positioning strategy based on a continuous Runge-Kutta algorithm is used to accurately position a key simulation event. According to the invention, high efficiency and robustness of simulation can be ensured. In combination with a key event positioning strategy, the fault analysis of the electrical integrated energy system can be efficiently and accurately realized.

Description

Technical Field

[0001] The present invention belongs to the field of integrated energy system simulation analysis, and specifically, relates to a fault analysis method for an electrical integrated energy system based on semi-implicit simulation. Background Art

[0002] As a clean alternative to fossil fuels represented by coal, natural gas has become the mainstream choice for power generation scenarios around the world. With gas turbines and power-to-gas (P2G) devices, the power system has obtained sufficient flexibility resources, so that electricity and natural gas are more closely coupled together.

[0003] While fully utilizing the advantages of electrical coupling operation, power equipment is also more susceptible to fluctuations in the state of natural gas. Faults such as pipe bursts and leaks will propagate at a faster speed in the natural gas network, causing equipment such as gas turbines to be cut off, affecting the normal operation of the power network. Therefore, it is necessary to quantify the propagation speed of faults such as pipe bursts and leaks through simulation methods, so that the control and protection devices can issue instructions in time to prevent the expansion of faults and cascading failures.

[0004] The simulation model of the electrical integrated energy system considering pipe burst and leakage faults has not yet been established. Considering that it is in the form of a set of high-dimensional, nonlinear, rigid partial differential-algebraic equations, the traditional iterative method based on finite difference / characteristic line method is inefficient, which seriously affects the efficiency of fault analysis. Summary of the invention

[0005] Technical problem: The technical problem to be solved by the present invention is to establish an accurate simulation model of pipe burst and leakage in natural gas network and provide an efficient simulation method suitable for fault analysis of electrical integrated energy system.

[0006] Technical solution: To solve the above technical problems, the present invention adopts the following technical solution:

[0007] A method for fault analysis of an electrical integrated energy system based on semi-implicit simulation comprises the following steps:

[0008] Step 1: Model natural gas pipe burst and leakage failures:

[0009] Step 101: Use a third-order WENO format to spatially discretize the isothermal natural gas dynamic equation to obtain an ordinary differential equation for pressure p and an ordinary differential equation group for flow rate q;

[0010] Step 102: Based on the ordinary differential equation about pressure p and the ordinary differential equation group about flow q obtained in step 101, the boundary conditions of pipeline burst fault and pipeline leakage fault are established respectively by using the characteristic line method, so as to obtain the upstream and downstream leakage flow-pressure algebraic equations of the natural gas pipeline burst fault position and the upstream and downstream leakage flow-pressure algebraic equations of the leakage position;

[0011] Step 2: Semi-implicit dynamic simulation solution based on rigid exact Rosenbrock format:

[0012] Using a semi-implicit dynamic simulation method based on a rigid and precise Rosenbrock format, jointly solve the ordinary differential equations about pressure p and the ordinary differential equations about flow q obtained in steps 101 and 102, the algebraic equations of leakage flow-pressure upstream and downstream of the natural gas pipeline burst fault location, and the algebraic equations of leakage flow-pressure upstream and downstream of the leakage location;

[0013] Step 3: Key event location strategy based on continuous Runge-Kutta:

[0014] In the process of solving the problem using the semi-implicit dynamic simulation method based on the rigid exact Rosenbrock format as described in step 2, the continuous Runge-Kutta method is used to locate the propagation of the limit pressure caused by the pipeline burst or leakage failure to obtain the time and location of the pipeline burst failure and / or pipeline leakage failure propagation.

[0015] Preferably, in step 102, the upstream and downstream leakage flow-pressure algebraic equations of the natural gas pipeline burst fault location are as follows:

[0016]

[0017] Where: j b represents the discrete location corresponding to the natural gas pipeline burst fault location, denoted as the burst location, j b The value of is [0, M], where M represents the total number of discrete locations of the natural gas pipeline; 爆管,上游 represents the burst pipe location j b The upstream leakage flow, q 爆管，下游 represents the burst pipe location j b The downstream leakage flow rate, S represents the pipe cross-sectional area, c represents the speed of sound, represents the discrete site j b -1 pressure; represents the discrete site j b -1 flow rate; represents the discrete site j b -2 pressure; represents the discrete site j b -2 flow rate; represents the discrete site j b The pressure; represents the discrete site j b +2 pressure; represents the discrete site j b +2 traffic; represents the discrete site j b +1 pressure; represents the discrete site j b +1 traffic;j b +1 indicates the burst pipe location j b The next discrete site, j b -1 represents the burst pipe location j b The previous discrete position of b +2 indicates the burst pipe location j b The last two discrete positions of j b -2 indicates the burst pipe location j b The first two discrete positions of .

[0018] Preferably, the upstream and downstream leakage flow-pressure algebraic equations of the leakage position of the natural gas pipeline are:

[0019]

[0020] q 泄露,上游 =q 泄露 +q 泄露，下游

[0021] Where: j x The discrete location corresponding to the leakage position of the natural gas pipeline is denoted as the leakage location, j x The value of is [0, M], where M represents the total number of discrete locations of the natural gas pipeline; 泄露,上游 Indicates the leakage site j x The upstream leakage flow, q 泄露,下游 Indicates the leakage site j x The downstream leakage flow rate; S represents the pipe cross-sectional area, c represents the speed of sound, represents the discrete site j x +2 pressure; represents the discrete site j x +2 traffic; represents the discrete site j x +1 pressure; represents the discrete site j x +1 traffic; represents the discrete site j x The pressure; represents the discrete site j x -2 pressure; represents the discrete site jx -2 flow rate; represents the discrete site j x -1 pressure; represents the discrete site j x -1 flow; j x +1 indicates leakage site j x The next point, j x -1 represents leakage site j x The previous point, j x +2 indicates leakage site j x The last two discrete positions of j x -2 indicates leakage site j x The distance between two adjacent discrete locations on the natural gas pipeline is Δx;

[0022] If the leaked site j x The pressure is less than or equal to the critical pressure p cr , leakage flow q at the leakage location 泄露 for:

[0023]

[0024] If the leaked site j x The pressure is greater than the critical pressure p cr , leakage flow q at the leakage location 泄露 for:

[0025]

[0026] Where k is the adiabatic coefficient, M is the molar mass of natural gas, T = 273.15 Kelvin, Z is the natural gas compressibility coefficient, R is the ideal gas constant, S h is the leakage area;

[0027] Critical pressure P cr for:

[0028]

[0029] Where, Hcr = 1.3; P a Preferably, in step three, the continuous Runge-Kutta method is used to locate the ultimate pressure propagation caused by pipeline burst or leakage failure, specifically, the following equation is constructed:

[0030]

[0031] Among them, p i (t) represents the pressure p of the gas grid node i in the electrical integrated energy system i As a function of time t, p criticalIt indicates the ultimate pressure after a natural gas pipeline bursts or leaks. P2G (t) represents the output flow at the node where the power-to-gas equipment P2G is located in the electrical integrated energy system, q max Indicates the maximum flow rate that can be output at the node where the power-to-gas device P2G is located.

[0032] The semi-implicit dynamic simulation method based on the rigid exact Rosenbrock format described in step 2 specifically constructs a solution based on the embedded Runge-Kutta method in the time interval [t0, t1]:

[0033]

[0034] Among them, b m (t) represents the polynomial with given coefficients determined by the embedded Runge–Kutta method; p node (t0) represents the value of the node pressure at time t0; q P2G (t0) represents the value of the flow rate output at the node where the power-to-gas device P2G is located at time t0; represents the intermediate coefficient of the nodal pressure calculated by the rigid exact Rosenbrock scheme described in step 2; represents the intermediate level coefficient of the flow rate output at the node where the power-to-gas device P2G is located calculated by the rigid and accurate Rosenbrock format described in step 2;

[0035] The node pressure t node (t) and the output flow q at the node where the power-to-gas device P2G is located P2G Substituting the expression of (t) into formula (2) and solving a univariate polynomial equation about t, we can obtain the time t required for the ultimate pressure to propagate to each node and for the P2G device to produce natural gas flow saturation cr .

[0036] Preferably, in step 101, the isothermal natural gas dynamic equation is as follows:

[0037]

[0038] Where p represents the pressure at any position on the natural gas pipeline, q represents the flow rate at any position on the natural gas pipeline; c represents the speed of sound; S represents the cross-sectional area of ​​the natural gas pipeline; λ represents the pipeline friction coefficient of the natural gas pipeline; D represents the pipeline diameter of the natural gas pipeline; t represents time, and x represents the length of the natural gas pipeline

[0039] Beneficial effects: Compared with the prior art, the present invention has the following beneficial effects: an accurate dynamic model for pipe burst and leakage simulation of the natural gas network is established; the semi-implicit dynamic simulation method based on the rigid and precise Rosenbrock format has the characteristics of non-iterative, efficient and accurate, and can realize robust simulation of fault dynamics; the simulation event location technology based on the continuous Runge-Kutta method can accurately locate simulation events, and can provide guidance for the formulation of protection and control strategies. BRIEF DESCRIPTION OF THE DRAWINGS

[0040] Figure 1 It is a flow chart of the electrical integrated energy system fault analysis method based on semi-implicit simulation described in the present invention.

[0041] Figure 2 It is the partial differential equation difference template used in the present invention.

[0042] Figure 3 It is a schematic diagram of pipe burst fault modeling in the present invention.

[0043] Figure 4 It is a schematic diagram of leakage fault modeling in the present invention.

[0044] Figure 5 This is the simulation verification result of Example 1 of the present invention.

[0045] Figure 6 It is a curve diagram showing the change of the limit removal time with the leakage position under different leakage holes in Example 2.

[0046] Figure 7 It is a schematic diagram of the electrical integrated energy system described in Example 3.

[0047] Figure 8 This is the mass flow diagram injected into node 1 in the electrical integrated energy system described in Example 3.

[0048] Fig. 9 It is the simulation result of each node in the electrical integrated energy system described in Example 3.

[0049] Fig.10 It is a schematic diagram of the cascading fault in the electrical integrated energy system described in Example 3 captured using the fault analysis method described in the present invention. DETAILED DESCRIPTION

[0050] In order to make the purpose, technical solution and advantages of the present invention more clearly understood, the present invention is described in detail below in conjunction with the accompanying drawings and implementation examples. It should be understood that the specific implementation examples described herein are only used to explain the present invention and are not used to limit the invention.

[0051] like Figure 1As shown, the present invention proposes a method for fault analysis of an electric integrated energy system based on semi-implicit simulation, which specifically includes the following steps:

[0052] Step 1: Model natural gas pipe burst and leakage failures:

[0053] Step 101: Use a third-order WENO format to spatially discretize the isothermal natural gas dynamic equation to obtain an ordinary differential equation for pressure p and an ordinary differential equation group for flow rate q;

[0054] Among them: The isothermal natural gas dynamic equation is a partial differential equation that takes into account the spatiotemporal dynamics of natural gas network transportation, as follows:

[0055]

[0056] Where p and q are the distribution of pressure (unit: Pa) and flow (unit: kg / s) in the natural gas pipeline in terms of time and space; c represents the speed of sound (m / s); S is the cross-sectional area of ​​the natural gas pipeline (m 2 ), λ is the pipeline friction coefficient of the natural gas pipeline, D is the pipeline diameter of the natural gas pipeline (m), t represents time, and x represents the length of the natural gas pipeline.

[0057] A third-order WENO format is used to spatially discretize the isothermal natural gas dynamic equations. Specifically, Figure 2 As shown, firstly, points 0, 1, 2 to M are selected at intervals of Δx on a natural gas pipeline of length L. Then, a third-order WENO format is used to discretize the spatial derivatives of the above isothermal natural gas dynamic equations, and the ordinary differential equations about the pressure p and the ordinary differential equations about the flow q can be obtained respectively.

[0058] The ordinary differential equation for pressure p is as follows:

[0059]

[0060] The ordinary differential equation for flow rate q is as follows:

[0061]

[0062] Where j represents a discrete point, and the numerical flux at points j+1 / 2 and j-1 / 2 is and Respectively meet:

[0063]

[0064] Below is the negative flux at point j+1 / 2 and the positive flux at point j-1 / 2 For example, give its expression: You can do this by Subtract 1 from the discrete sites in the expression to obtain, You can do this by The discrete sites in the expression of are obtained by adding 1.

[0065]

[0066] The differential template 0 and differential template 1 at point j+1 / 2 respectively satisfy:

[0067]

[0068] Where ω0 and ω1 are The weight factor of the differential template satisfies:

[0069]

[0070] Among them, the template selection coefficients α0 and α1 satisfy

[0071]

[0072] Among them, the smoothness characterization coefficients β0 and β1 satisfy

[0073]

[0074] in and for The weight factor of the differential template satisfies

[0075]

[0076] Among them, the template selection coefficient satisfy

[0077]

[0078] The differential template 0 and differential template 1 at point j-1 / 2 respectively satisfy

[0079]

[0080] In the above formula, for the ordinary differential equation about pressure p, the value of and Satisfy

[0081]

[0082] For the ordinary differential equation for pressure q, the pressure at point j is and Satisfy

[0083]

[0084] Step 102: Based on the ordinary differential equation about pressure p and the ordinary differential equation group about flow q obtained in step 101, the boundary conditions of pipeline burst failure and pipeline leakage failure are established respectively using the characteristic line method, so as to obtain the upstream and downstream leakage flow-pressure algebraic equations of the natural gas pipeline burst failure location and the upstream and downstream leakage flow-pressure algebraic equations of the leakage location.

[0085] Specifically, Figure 3 As shown in the figure, the characteristic line method is used to establish the boundary conditions of pipeline burst failure. 爆管,上游 With q 爆管，下游 They represent the upstream and downstream leakage flow rates after the pipe burst, Π air Represents the ambient pressure. The upstream and downstream leakage flow-pressure algebraic equations of the natural gas pipeline burst fault location are as follows:

[0086]

[0087] In the formula: In the formula: j b represents the discrete location corresponding to the natural gas pipeline burst fault location, denoted as the burst location, j b The value of is [0, M], where M represents the total number of discrete locations of the natural gas pipeline; 爆管,上游 represents the burst pipe location j b The upstream leakage flow, q 爆管，下游 represents the burst pipe location j b The downstream leakage flow rate, S represents the pipe cross-sectional area, c represents the speed of sound, represents the discrete site j b -1 pressure; represents the discrete site j b -1 flow rate; represents the discrete site j b -2 pressure; represents the discrete site j b -2 flow rate; represents the discrete site j b The pressure; represents the discrete site j b +2 pressure; represents the discrete site j b +2 traffic; represents the discrete site j b +1 pressure; represents the discrete site j b +1 traffic;j b +1 indicates the burst pipe location j b The next discrete site, j b-1 represents the burst pipe location j b The previous discrete position of b +2 indicates the burst pipe location j b The last two discrete positions of j b -2 indicates the burst pipe location j b The first two discrete positions of .

[0088] like Figure 4 As shown in the figure, the characteristic line method is used to establish the boundary conditions of pipeline leakage failure. j1 With q j2 They represent the flow rates upstream and downstream of the leakage location, q leak Represents the leakage flow. The upstream and downstream leakage flow-pressure algebraic equations of the leakage position on the natural gas pipeline are:

[0089]

[0090] q 泄露,上游 =q 泄露 +q 泄露，下游

[0091] Where: j x The discrete location corresponding to the leakage position of the natural gas pipeline is denoted as the leakage location, j x The value of is [0, M], where M represents the total number of discrete locations of the natural gas pipeline; 泄露,上游 Indicates the leakage site j x The upstream leakage flow, q 泄露,下游 Indicates the leakage site j x The downstream leakage flow rate; S represents the pipe cross-sectional area, c represents the speed of sound, represents the discrete site j x +2 pressure; represents the discrete site j x +2 traffic; represents the discrete site j x +1 pressure; represents the discrete site j x +1 traffic; represents the discrete site j x The pressure; represents the discrete site j x -2 pressure; represents the discrete site j x -2 flow rate; represents the discrete site j x -1 pressure; represents the discrete site j x -1 flow; j x +1 indicates leakage site j x The next point, jx -1 represents leakage site j x The previous point, j x +2 indicates leakage site j x The last two discrete positions of j x -2 indicates leakage site j x The first two discrete locations of the natural gas pipeline are Δx.

[0092] If the leaked site j x The pressure is less than or equal to the critical pressure p cr , leakage flow q at the leakage location 泄露 for:

[0093]

[0094] If the leaked site j x The pressure p jx Greater than the critical pressure p cr , leakage flow q at the leakage location 泄露 for:

[0095]

[0096] Where k is the adiabatic coefficient, M is the molar mass of natural gas, T = 273.15 Kelvin, Z is the natural gas compressibility coefficient, R is the ideal gas constant, S h is the leakage area;

[0097] Critical pressure P cr for:

[0098]

[0099] Where, Hcr = 1.3; P a Indicates standard atmospheric pressure.

[0100] Step 2: Semi-implicit dynamic simulation solution based on rigid exact Rosenbrock format:

[0101] Using a semi-implicit dynamic simulation method based on a rigid and precise Rosenbrock format, jointly solve the ordinary differential equations about pressure p and the ordinary differential equations about flow q obtained in steps 101 and 102, the algebraic equations of leakage flow-pressure upstream and downstream of the natural gas pipeline burst fault location, and the algebraic equations of leakage flow-pressure upstream and downstream of the leakage location;

[0102] In this step, a semi-implicit dynamic simulation method based on a fourth-order rigid exact Rosenbrock format is used to solve the above-mentioned natural gas network dynamic model and the nonlinear time-domain simulation model of the power system.

[0103] Specifically, for a differential algebraic equation of the form

[0104]

[0105] Where M is the singular mass matrix of the differential algebraic equation system, y(t) is the vector of variables to be solved in the differential algebraic equation system, and F represents the differential algebraic equation expression. In any time interval [t0, t1], use the following format to solve it:

[0106]

[0107] Where h represents the simulation step size,

[0108]

[0109] t0 is the initial time, y0 represents the initial value of the variable, u m With k m represents the coefficients of each level obtained by recursive solution, The other parameters are: γ = 0.25, α 21 =0.386, α 31 =0.146,α 32 =0.0639, α 41 =-0.331, α 42 =0.711,α 43 =0.250,α 51 =-4.55, α 52 =1.71,α 53 =4.01,α 54 =-1.72, α 61 =2.43,α 62 =-3.83, α 63 =-1.86, α 64 =0.560, α 65 =0.250,γ 21 =-0.354,γ 31 =-0.134,γ 32 =-0.0129,γ 41 =1.527,γ 42 =-0.534,γ 43 =-1.279,γ 51 =6.981,γ 52 =-2.093,γ 53 =-5.87,γ 54 =0.732,γ 61 =-2.08,γ 62 =0.596,γ63 =1.702,γ 64 =-0.089,γ 65 =-0.379, b1=0.348, b2=0.213, b3=-0.154, b4=0.471, b5=-0.129, b6=0.25.

[0110] The embedded Runge-Kutta method is used to estimate the truncation error and adaptively adjust the simulation step size. Specifically, let

[0111]

[0112] in, The error constructed satisfies Embedded solution.

[0113] Use the expression:

[0114]

[0115] Make an estimate of the local error. And use:

[0116]

[0117] to adjust the step size, where Atol and Rtol are the preset absolute and relative error limits respectively, h represents the simulation step size, and q represents the order of the format.

[0118] Step 3: Key event location strategy based on continuous Runge-Kutta:

[0119] In the process of solving the problem using the semi-implicit dynamic simulation method based on the rigid exact Rosenbrock format as described in step 2, the continuous Runge-Kutta method is used to locate the propagation of the limit pressure caused by the pipeline burst or leakage failure to obtain the time and location of the pipeline burst failure and / or pipeline leakage failure propagation.

[0120] Specifically, in order to obtain the time required for the ultimate pressure to propagate to each gas grid node in the electrical integrated energy system and the node where the power-to-gas equipment P2G is located to produce natural gas flow saturation after a pipe burst or leakage failure occurs, the following equation is constructed:

[0121]

[0122] Among them, p i (t) represents the pressure p of the gas grid node i in the electrical integrated energy system i As a function of time t, p critical It indicates the ultimate pressure after a natural gas pipeline bursts or leaks. P2G(t) represents the output flow at the node where the power-to-gas equipment P2G is located in the electrical integrated energy system, q max Indicates the maximum flow rate that can be output at the node where the power-to-gas device P2G is located.

[0123] In step 2, a semi-implicit dynamic simulation method based on the rigid exact Rosenbrock format is used to construct a solution based on the embedded Runge-Kutta method in the time interval [t0, t1]:

[0124]

[0125] Among them, b m (t) represents a polynomial with given coefficients determined by the embedded Runge-Kutta method; p node (t0) represents the value of the node pressure at time t0, that is, y0 in the format corresponds to p node part;q P2G (t0) represents the value of the output flow at the node where the power-to-gas device P2G is located at time t0, that is, y0 in the format corresponds to p P2G part of; represents k calculated by the rigid exact Rosenbrock scheme described in step 2 m The corresponding p node part; represents k calculated by the rigid exact Rosenbrock scheme described in step 2 m The corresponding p P2G part;

[0126] The p node (t) and q P2G Substituting the expression of (t) into (2) and solving a univariate polynomial equation about t, we can obtain the time t required for the ultimate pressure to propagate to each node and for the P2G device to produce natural gas flow saturation cr .

[0127] Three embodiments are given below.

[0128] The present invention is compared with finite difference and characteristic line methods.

[0129] The format for finite differences is:

[0130]

[0131] The format of the characteristic line method is:

[0132]

[0133] Among them, Δt and Δx are the time and space difference steps respectively.

[0134] Example 1

[0135] Take a natural gas pipeline with a length of 51 kilometers, a diameter of 0.5901 meters, and a friction coefficient of 0.5901 as an example. The spatial difference step size is 100 meters. A pipe burst occurred at the center of the pipeline. The result of the characteristic line method is used as the benchmark value. Figure 5 The simulation verification results are shown. The finite difference method shows strong numerical oscillation, while the present invention does not produce numerical oscillation and can accurately reflect the dynamic changes of leakage flow.

[0136] Example 2

[0137] Add leakage faults at different locations of the pipeline and calculate the time required for the pipeline outlet pressure to drop to 2.8MPa. When the pipeline outlet pressure drops to 2.8MPa, the outlet equipment needs to be cut off, so this time is called the limit cut-off time. Figure 6 The figure shows the change of the limit removal time with the change of the leakage position under different leakage hole sizes. It can be seen from the figure that the limit removal time gradually increases with the leakage hole principle pipeline outlet. And the larger the leakage hole, the smaller the limit removal time. When d / D=1, the limit removal time of the leakage model is similar to the limit removal time of the burst pipe model, which further verifies the accuracy of the fault model established by the present invention.

[0138] Example 3

[0139] exist Figure 7 The simulation is carried out in the electrical integrated energy system shown in FIG. The power grid and the natural gas grid are bidirectionally coupled to the gas turbine through the power-to-gas device. At t = 300 seconds, a leakage fault occurred in the pipeline between gas grid node 2 and gas grid node 3. The diameter of the leakage hole was 60% of the pipeline diameter. After the output of the power-to-gas device was saturated, it was converted from a constant pressure node to a constant flow node. The spatial difference step size Δx = 100 meters, and the total simulation time was 16 hours.

[0140] Figure 8-10 The simulation results of Example 3 are shown. Using the key event location strategy, the present invention accurately locates the output saturation moment of the power-to-gas equipment, i.e., t=32423.00 seconds, and successfully calculates the limit removal moment of each natural gas node and the gas consumption curve of each gas turbine. In addition, at t=52561.62, the No. 2 gas turbine GT2 was removed from the system; due to the transfer of the power system load under the action of the governor, at t=52566.73, the No. 1 gas turbine GT1 was removed, and the method successfully captured the cascading failure of the system.

[0141] Table 1 compares the time consumed for simulation of the present invention and the characteristic line method, and it can be seen that the present invention has a significant efficiency advantage.

[0142] Table 1 Comparison of calculation time consumption in Example 3

[0143]

[0144] It can be seen from the above three embodiments that the present invention can realize efficient and accurate fault analysis of the electrical integrated energy system and provide guidance for the formulation of protection and control strategies.

Claims

1. A method for fault analysis of an electrical integrated energy system based on semi-implicit simulation, characterized in that: The steps include: Step 1: Model natural gas pipe burst and leakage failures: Step 101: Use a third-order WENO format to spatially discretize the isothermal natural gas dynamic equation to obtain an ordinary differential equation for pressure p and an ordinary differential equation group for flow rate q; Step 102: Based on the ordinary differential equation about pressure p and the ordinary differential equation group about flow q obtained in step 101, the boundary conditions of pipeline burst fault and pipeline leakage fault are established respectively by using the characteristic line method, so as to obtain the upstream and downstream leakage flow-pressure algebraic equations of the natural gas pipeline burst fault position and the upstream and downstream leakage flow-pressure algebraic equations of the leakage position; Step 2: Semi-implicit dynamic simulation solution based on rigid exact Rosenbrock format: Using a semi-implicit dynamic simulation method based on a rigid and precise Rosenbrock format, jointly solve the ordinary differential equations about pressure p and the ordinary differential equations about flow q obtained in steps 101 and 102, the algebraic equations of leakage flow-pressure upstream and downstream of the natural gas pipeline burst fault location, and the algebraic equations of leakage flow-pressure upstream and downstream of the leakage location; Step 3: Key event location strategy based on continuous Runge-Kutta: In the process of solving the problem using the semi-implicit dynamic simulation method based on the rigid exact Rosenbrock format as described in step 2, the continuous Runge-Kutta method is used to locate the propagation of the limit pressure caused by the pipeline burst or leakage failure to obtain the time and location of the pipeline burst failure and / or pipeline leakage failure propagation.

2. According to claim 1, a method for fault analysis of an electrical integrated energy system based on semi-implicit simulation is characterized in that: In step 102, the upstream and downstream leakage flow-pressure algebraic equations of the natural gas pipeline burst fault location are as follows: Where: j b represents the discrete location corresponding to the natural gas pipeline burst fault location, denoted as the burst location, j b The value of is [0, M], where M represents the total number of discrete locations of the natural gas pipeline; 爆管,上游 represents the burst pipe location j b The upstream leakage flow, q 爆管，下游 represents the burst pipe location j b The downstream leakage flow rate, S represents the pipe cross-sectional area, c represents the speed of sound, represents the discrete site j b -1 pressure; represents the discrete site j b -1 flow rate; represents the discrete site j b -2 pressure; represents the discrete site j b -2 flow rate; represents the discrete site j b The pressure; represents the discrete site j b +2 pressure; represents the discrete site j b +2 traffic; represents the discrete site j b +1 pressure; represents the discrete site j b +1 traffic;j b +1 indicates the burst pipe location j b The next discrete site, j b -1 represents the burst pipe location j b The previous discrete position of b +2 indicates the burst pipe location j b The last two discrete positions of j b -2 indicates the burst pipe location j b The first two discrete positions of .

3. The method for fault analysis of an electrical integrated energy system based on semi-implicit simulation according to claim 1 is characterized in that: The upstream and downstream leakage flow-pressure algebraic equations of the leakage position of the natural gas pipeline are: q 泄露,上游 =q 泄露 +q 泄露，下游 Where: j x The discrete location corresponding to the leakage position of the natural gas pipeline is denoted as the leakage location, j x The value of is [0, M], where M represents the total number of discrete locations of the natural gas pipeline; 泄露,上游 Indicates the leakage site j x The upstream leakage flow, q 泄露,下游 Indicates the leakage site j x Downstream leakage flow; S represents the cross-sectional area of ​​the pipe, c represents the speed of sound, represents the discrete site j x +2 pressure; represents the discrete site j x +2 traffic; represents the discrete site j x +1 pressure; represents the discrete site j x +1 traffic; represents the discrete site j x The pressure; represents the discrete site j x -2 pressure; represents the discrete site j x -2 flow rate; represents the discrete site j x -1 pressure; represents the discrete site j x -1 flow; j x +1 indicates leakage site j x The next point, j x -1 represents leakage site j x The previous point, j x +2 indicates leakage site j x The last two discrete positions of j x -2 indicates leakage site j x The distance between two adjacent discrete locations on the natural gas pipeline is Δx; If the leaked site j x The pressure is less than or equal to the critical pressure p cr , leakage flow q at the leakage location 泄露 for: If the leaked site j x The pressure is greater than the critical pressure p cr , leakage flow q at the leakage location 泄露 for: Where k is the adiabatic coefficient, M is the molar mass of natural gas, T = 273.15 Kelvin, Z is the natural gas compressibility coefficient, R is the ideal gas constant, S h is the leakage area; Critical pressure P cr for: Where, Hcr = 1.3; P a Indicates standard atmospheric pressure.

4. The method for fault analysis of an electrical integrated energy system based on semi-implicit simulation according to claim 1 is characterized in that: In step three, the continuous Runge-Kutta method is used to locate the ultimate pressure propagation caused by pipeline burst or leakage failure. Specifically, the following equation is constructed: Among them, p node (t) represents the function of the gas grid node pressure in the electrical integrated energy system with respect to time t, p critical It indicates the ultimate pressure after a natural gas pipeline bursts or leaks. P2G (t) represents the function of the output flow rate at the node where the power-to-gas device P2G in the electrical integrated energy system is located with respect to time, q max Indicates the maximum flow rate that can be output at the node where the power-to-gas device P2G is located. The semi-implicit dynamic simulation method based on the rigid exact Rosenbrock format described in step 2 specifically constructs a solution based on the embedded Runge-Kutta method in the time interval [t0, t1]: Among them, b m (t) represents a polynomial with given coefficients determined by the embedded Runge-Kutta method; p node (t0) represents the value of the node pressure at time t0; q P2G (t0) represents the value of the flow rate output at the node where the power-to-gas device P2G is located at time t0; represents the intermediate coefficient of the nodal pressure calculated by the rigid exact Rosenbrock scheme described in step 2; represents the intermediate level coefficient of the flow rate output at the node where the power-to-gas device P2G is located calculated by the rigid and accurate Rosenbrock format described in step 2; The p node (t) and q P2G Substituting the expression of (t) into (2) and solving a univariate polynomial equation about t, we can obtain the time t required for the ultimate pressure to propagate to each node and for the P2G device to produce natural gas flow saturation cr .

5. The method for fault analysis of an electrical integrated energy system based on semi-implicit simulation according to claim 1 is characterized in that: In step 101, the isothermal natural gas dynamic equation is as follows: Among them, p represents the pressure at any position on the natural gas pipeline, q represents the flow rate at any position on the natural gas pipeline; c represents the speed of sound; S is the cross-sectional area of ​​the natural gas pipeline; λ is the pipeline friction coefficient of the natural gas pipeline; D is the pipeline diameter of the natural gas pipeline; t represents time, and x represents the length of the natural gas pipeline.

Images