Comprehensive energy system load recovery method based on risk assessment and scene matching

Abstract

The invention discloses an integrated energy system load recovery method based on risk assessment and scene matching. The method comprises the following steps: firstly, generating a typical fault scene set under a typhoon disaster based on a Monte Carlo mode method; then constructing a comprehensive energy system load recovery model; after disasters occur and the energy system breaks down, a recovery strategy is generated regularly by using the comprehensive energy system load recovery model, and the recovery strategy is input to carry out load recovery; when the recovery strategy is generated, matching is carried out on the typical fault scene set, then the weight in the comprehensive energy system load recovery model is obtained according to the matching result and historical data, and then the recovery strategy is obtained based on the comprehensive energy system load recovery model. According to the method, the weight is optimized through dynamic scene matching, a multi-dimensional risk evaluation system is introduced, thermal and natural gas dynamic constraints are considered, and the load recovery effect and the system stability are effectively improved.

Description

Technical Field

[0001] The present invention relates to the technical field of power system dispatching, and particularly relates to a method for load restoration of an integrated energy system based on risk assessment and scenario matching. Background Art

[0002] With the wide application of coupling devices such as combined heat and power units, gas turbines, and electric boilers, traditional distribution networks have gradually evolved into integrated energy systems with multi-energy complementarity of electricity, heat, and gas. Energy conversion and coordinated dispatching are achieved through coupling devices, significantly improving energy utilization efficiency. However, the deep coupling of multi-energy systems also brings new challenges: failures in a single subsystem may spread to other systems through coupling devices, triggering cascading failures. Especially under extreme natural disasters such as typhoons, large-scale damage to key components of the system (such as power lines and gas pipelines) will lead to low load restoration efficiency and even cause long-term and large-scale energy supply interruptions.

[0003] Although the existing technologies for risk assessment and restoration strategies of power systems have been relatively mature, there are still the following limitations:

[0004] 1. Most existing load restoration models generate strategies based on a fixed set of fault scenarios, while the fault scenarios after typhoon disasters are dynamically changing. Moreover, there are dynamic constraints and multi-energy coupling characteristics in the thermal and natural gas parts of the integrated energy system. The existing methods cannot adapt to the changing and complex system states after disasters, resulting in insufficient matching between the restoration strategy and the real-time disaster situation and poor load restoration effect.

[0005] 2. Most existing index systems start from the global risk of the system layer and lack refined evaluation of local weak links such as key nodes and coupling devices, resulting in unsatisfactory operational stability of the system after restoration. Summary of the Invention

[0006] The present invention proposes a method for load restoration of an integrated energy system based on risk assessment and scenario matching, and its purpose is to solve the problems of poor load restoration effect and unsatisfactory stability.

[0007] The technical solution of the present invention is as follows:

[0008] A method for load restoration of an integrated energy system based on risk assessment and scenario matching, the steps include:

[0009] Step 1: Generate a set of typical fault scenarios under typhoon disasters based on the Monte Carlo method;

[0010] Step 2: Construct a load restoration model for the integrated energy system. The objective function of the load restoration model for the integrated energy system includes weights corresponding to different energy types respectively, and the constraint conditions of the load restoration model for the integrated energy system include multi-dimensional risk constraints;

[0011] Step 3: After a disaster occurs and the energy system fails, regularly use the integrated energy system load restoration model to generate a restoration strategy and implement the restoration strategy for load restoration;

[0012] When generating the restoration strategy, first match the typical fault scenario closest to the current situation from the typical fault scenario set, then optimize the weights in the integrated energy system load restoration model based on the historical data corresponding to the typical fault scenario, and finally obtain the restoration strategy based on the integrated energy system load restoration model.

[0013] As a further improvement of the integrated energy system load restoration method based on risk assessment and scenario matching, Step 1 specifically includes:

[0014] Step 1-1: Build a typhoon disaster simulation model for calculating the real-time wind speed at different locations in the affected area according to typhoon information;

[0015] Step 1-2: Calculate the probability distribution of the starting time of power line faults under the influence of typhoon disasters;

[0016] Step 1-3: Obtain a fault scenario set based on the Monte Carlo method and the probability distribution of the starting time of each power branch fault, and then cluster the fault scenario set to obtain a typical fault scenario set.

[0017] As a further improvement of the integrated energy system load restoration method based on risk assessment and scenario matching:

[0018] In Step 1-1, the calculation process of the typhoon disaster simulation model is as follows:

[0019] Step 1-1-1: Calculate the pressure difference between the typhoon center and the outside at time t:

[0020] ΔH(t) = 0.75ΔH 0 - 0.508[1 + sinθ]t;

[0021] In the formula: ΔH 0 is the initial pressure difference between the typhoon center and the outside; θ is the angle between the moving direction of the typhoon center at time t and the horizontal axis of the coordinate system;

[0022] Step 1-1-2: Calculate the gradient wind speed and maximum wind speed radius of the typhoon at each moment according to the pressure difference between the typhoon center and the outside:

[0023]

[0024] R max (t) = exp(5.104 - 1.124ΔH(t) 0.6 );

[0025] v gw(t) is the gradient wind speed at time t, R max (t) is the maximum wind speed radius at time t;

[0026] Step 1-1-3: Calculate the maximum wind speed v within the typhoon influence range at time t according to the gradient wind speed and the typhoon movement speed: Rmax (t):

[0027] v Rmax (t) = 0.865v gw (t) + 0.5v T ;

[0028] Step 1-1-4: Calculate the wind speed at any position within the integrated energy system:

[0029]

[0030] In the formula: v ra (t) is the real-time wind speed at a certain position at time t, and d(t) is the real-time distance between the typhoon center position and this position at time t.

[0031] As a further improvement of the integrated energy system load restoration method based on risk assessment and scenario matching:

[0032] In Step 1-2, for the i-th power branch, according to the wind speed v i (t) borne at its location at time t, calculate its real-time probability h i (t) of failure at time t in the following way:

[0033]

[0034] In the formula, v des is the designed wind speed of this power branch; is the maximum wind speed that this power branch can withstand; is a preset model parameter;

[0035] The probability p i,t that the failure start time of this power branch is at time t is:

[0036]

[0037] As a further improvement of the integrated energy system load restoration method based on risk assessment and scenario matching, Step 1-3 includes:

[0038] Step 1-3-1: Randomly construct multiple different typhoon information, where the typhoon information includes the typhoon center trajectory, the typhoon center moving speed, and the initial pressure difference between the typhoon center and the outside; at the same time, randomly construct multiple different system states for the integrated energy system, where the system state includes the load demands of each node in the system at different times, the output of each coupling device, and the operating state of each coupling device;

[0039] Randomly combine the constructed typhoon information and system states to obtain multiple operation scenarios;

[0040] Step 1-3-2: For each operation scenario, calculate the probability distribution of the start time of all power line faults respectively, and then use the Monte Carlo simulation method to randomly generate fault samples based on this probability distribution. The fault sample refers to the fault scenario where a certain power branch starts to fail at a certain time under the corresponding operation scenario, so as to obtain a fault scenario set composed of multiple fault scenarios; each fault scenario corresponds to a feature vector, and the feature vector includes the corresponding typhoon information, system state, as well as the fault start time and the power branch number where the fault occurs;

[0041] Step 1-3-3: Cluster the fault scenario set according to the feature vectors of each fault scenario, and the cluster centers in each obtained cluster are the typical fault scenarios, and the obtained typical fault scenarios form a typical fault scenario set.

[0042] As a further improvement of the integrated energy system load restoration method based on risk assessment and scenario matching, the objective function of the integrated energy system load restoration model is:

[0043]

[0044] where, Ω L is the multi-energy load set of the regional integrated energy system, including electric load, thermal load and natural gas load; and are respectively the reduction amounts of electric power, natural gas, and thermal energy supply to the load node i at the moment t in the formula when the integrated energy system implements the restoration strategy next time, which are decision variables; α e 、α g and α h are respectively the weights of the three energy types of electric power, natural gas, and thermal energy.

[0045] As a further improvement of the integrated energy system load restoration method based on risk assessment and scenario matching, the constraint conditions of the integrated energy system load restoration model include:

[0046] (2.1) Power system constraints:

[0047] For any load node j:

[0048]

[0049] Wherein, i, j, and k are the numbers of load nodes, ij represents the power branch with the starting point at load node i and the ending point at load node j, and jk represents the power branch with the starting point at load node j and the ending point at load node k; Ω e is the set of load nodes in the topology after the distribution network fails; represents Ω e all the power branches in which the starting point is load node j, represents Ω e all the power branches in which the ending point is load node j; P jk,t , Q jk,t are respectively the active power transmission power and the reactive power transmission power of branch jk at time t, P ij,t , Q tj,t are respectively the active power transmission power and the reactive power transmission power of branch ij at time t, which are decision variables; I ij,t is the current of branch ij at time t; r ij , x ij are respectively the resistance and reactance of branch ij; are respectively the active and reactive power outputs of load node j as an energy supply device at time t, which are decision variables; are respectively the active and reactive power loads of load node j at time t, which are decision variables; are respectively the active and reactive power losses of load node j at time t;

[0050] For any load node j and all power branches ij with the ending point at load node j:

[0051]

[0052] For any load node i:

[0053]

[0054] Wherein, V i,t , V j,t are respectively the voltage values of load node i and load node j at time t; z ij,t is the on-off state of branch ij at time t. If the branch ij is disconnected due to a typhoon, it is represented as 0, otherwise it is 1; M is a preset constant; V 0 is the reference voltage value; are respectively the upper and lower limits of the active power transmission power of branch ij; are respectively the upper and lower limits of the reactive power transmission power of branch ij; are respectively the upper and lower limits of the voltage of load node i;

[0055] (2.2) Natural gas system dynamic constraints:

[0056] For any load node i:

[0057]

[0058] In the formula, are the natural gas flow rates flowing out and into pipeline ij connecting load node i and load node j at time t, respectively; is the set of all natural gas nodes connected to load node i through pipelines and with load node i as the gas outlet end; is the gas supply volume from the gas source node to load node i at time t; is the gas load of load node i at time t, which is a decision variable; is the set of energy conversion devices using the natural gas of load node i as the energy source; is the gas load of energy conversion device k at time t, which is a decision variable; is the reduction amount of natural gas energy supply of load node i, which is a decision variable;

[0059]

[0060] In the formula, G ij,t is the average pipeline flow rate of pipeline ij at time t, K ij is the preset pipeline parameter of pipeline ij, D ij is the natural gas flow direction of pipeline ij, μ i,t and μ j,t are the gas pressures of load node i and load node j at time t, respectively, which are known quantities;

[0061] For any pipeline ij:

[0062]

[0063] For any load node i:

[0064]

[0065] μ min ≤ μ i,t ≤ μ max ;

[0066] In the formula, are the upper and lower limits of the natural gas pipeline flow rate, G s,max and G s,min are the upper and lower limits of the gas source output, μ max and μ min are the upper and lower limits of the node gas pressure, respectively;

[0067] (2.3) Thermal system dynamic constraint:

[0068]

[0069] Wherein, are the heating powers of the electric boiler and the CHP unit in the integrated energy system at time t, respectively, which are decision variables; c p is the specific heat capacity of water; m ij,t is the sectional flow rate of the heating pipeline ij at time t; T sj,t and T ri,t are the temperatures at the end of the heating pipeline ij and the beginning of the corresponding return water pipeline at time t, respectively, which are known quantities; is the heat load loss at the node at time t, which is a decision variable;

[0070]

[0071] Wherein, is the heat load at the node at time t;

[0072] For any load node h:

[0073]

[0074] Wherein: are the sets of heating pipelines with load node h as the first head and the sets of heating pipelines with load node h as the first end, respectively; respectively represent the outlet temperature of the heating part and the outlet temperature of the heat regeneration part of the heating pipeline k; are the heating temperature and the heat regeneration temperature of load node g, respectively; are the flow rates of the heating part and the heat regeneration part of pipeline k, respectively;

[0075] For the heating temperature T of any load node at time t s,t and the heat regeneration temperature T of any load node r,t and the sectional flow rate m of any heating pipeline ij ij,t :

[0076]

[0077] m min ≤ m ij,t ≤ m max ;

[0078] Wherein, are the upper and lower limits of the temperatures of each node in the water supply network, respectively, are the upper and lower limits of the temperatures of each node in the return water network, m max and m minis the upper and lower limits of the flow rate of the heating pipeline cross-section ij;

[0079] (2.4) Output constraint of the coupling device:

[0080]

[0081] In the formula, is the output electric power of the gas turbine at time t, which is a decision variable; is the natural gas consumed by the gas turbine at time t, which is a decision variable; L NG is the calorific value of natural gas, which is 9.78 J / kg; η GT is the gas-to-electricity conversion efficiency of the gas turbine; θ GT is the start-stop state of the gas turbine, which is a 0-1 variable, 0 for the off state and 1 for the on state, and is a decision variable;

[0082]

[0083] In the formula, are the electric and heat powers output by the CHP unit at time t, respectively, which are decision variables; θ CHP is the start-stop state of the CHP unit, which is a 0-1 variable, 0 for the off state and 1 for the on state, and is a decision variable; is the natural gas consumed by the CHP unit at time t, which is a decision variable; is the heat loss rate of the CHP unit; η GTE 、η GTH are the gas-electricity conversion efficiency and gas-heat conversion efficiency of the CHP unit, respectively;

[0084]

[0085] In the formula, is the heat power output by the electric boiler at time t, which is a decision variable; is the electric power consumed by the electric boiler at time t, which is a decision variable; η EB is the electricity-to-heat conversion efficiency of the electric boiler; is the heat power loss rate of the electric boiler; θ EB is the start-stop state of the electric boiler, which is a 0-1 variable, 0 for the off state and 1 for the on state, and is a decision variable;

[0086] (2.5) Output constraint of the energy supply equipment, i.e., the coupling device:

[0087]

[0088] In the formula, P eq,t is the power output of the coupling device at time t, which is a decision variable; are the upper and lower limits of the power output of the coupling device, respectively.

[0089] As a further improvement of the above-mentioned integrated energy system load restoration method based on risk assessment and scenario matching: the multi-dimensional risk constraints specifically include system-level risk indicators, node-level risk indicators, and coupled equipment-level risk indicators; when solving the integrated energy system load restoration model, each indicator in the multi-dimensional risk constraints must be within a preset range.

[0090] As a further improvement of the above-mentioned integrated energy system load restoration method based on risk assessment and scenario matching:

[0091] The calculation method of the system-level risk indicator is:

[0092]

[0093] In the formula, R sys,t is the average load shedding rate of the system at time t; α e , α g , α h are the weights of electricity, gas, and heat energy respectively; and are the shedding amounts of electricity, gas, and heat loads of the system at time t respectively, which are decision variables; and are the total demands of electricity, gas, and heat loads at time t respectively, which are known quantities;

[0094] The calculation method of the node-level risk indicator is:

[0095] For any load node node:

[0096]

[0097] In the formula, ΔP node,t is the shedding amount of the system on the load node node at time t, which is a decision variable; L node,t is the load demand of the load node node at time t, which is a known quantity;

[0098] The calculation method of the coupled equipment-level risk indicator is:

[0099] For any coupled equipment j:

[0100] R eq,j,t =R sys,t -X eq,j,t ;

[0101]

[0102] In the formula: X eq,j,t is the average load shedding rate of the system at time t when the coupled equipment j fails; and The power, gas, and heat load reduction amounts of the system at time t when the coupling device j fails.

[0103] As a further improvement of the comprehensive energy system load restoration method based on risk assessment and scenario matching, in step 3, the specific steps for generating and implementing the restoration strategy are as follows:

[0104] Step 3-1: Collect typhoon information, system status, fault start time, and power branch numbers where faults occur in the current actual fault scenario, construct the current feature vector, then calculate the distances between the current feature vector and the feature vectors of each typical fault scenario, find the typical fault scenario corresponding to the minimum distance, and then obtain M sets of historical data of the comprehensive energy system during the post-disaster restoration process. The historical data includes the reduction amounts of power, natural gas, and heat energy supply for each load node.

[0105] Step 3-2: Construct the weight solution objective function:

[0106]

[0107] In the formula, α e , α g , and α h are the weights of the three energy types of power, natural gas, and heat respectively, and are the current values to be determined; and are the reduction amounts of power, natural gas, and heat energy supply for load node i by the comprehensive energy system at time t respectively;

[0108] Step 3-3: Set initial values for the weights α e , α g , and α h , and initialize the iteration count to 0;

[0109] Step 3-4: Substitute the M sets of historical data into the weight solution objective function, denote the objective function value obtained from the m-th set of historical data as obj m , calculate the average value and increment the iteration count by 1;

[0110] Step 3-5: Determine whether the current iteration count has reached the preset value. If it has, take the current weights α e , α g , α h as the optimal weight values for this type of fault scenario set, and jump to step 3-6. Otherwise, with the minimum of the average value as the goal, update the weight values α e , α g , α h using the gradient descent method, and then jump to step 3-4;

[0111] Step 3-6: Substitute the obtained weights α e , α g , α h into the integrated energy system load restoration model for solution to obtain a restoration strategy composed of decision variables; then implement the restoration strategy, and the control system performs load restoration until the next moment for calculating the restoration strategy arrives.

[0112] Compared with the prior art, the present invention has the following beneficial effects:

[0113] 1. The present invention establishes a dynamic scenario matching mechanism, adaptively optimizes the weights in the integrated energy system load restoration model according to the historical data of the successfully matched scenarios, making the restoration strategy output by the model closer to the actual scenario and improving the restoration effect.

[0114] 2. A multi-dimensional risk assessment index system is introduced into the integrated energy system load restoration model. The system-level index restricts the global energy supply loss, the node-level index diversifies the allocation of load shedding priorities, and the coupled equipment-level index inhibits the cross-system fault propagation, solving the problem that the traditional global assessment ignores local weak links, avoiding introducing new risks during the restoration process, and improving the stability of system operation.

[0115] 3. The dynamic constraints between the thermal and natural gas parts and the multi-energy coupling characteristics of the system are fully considered in the integrated energy system load restoration model, further improving the load restoration effect and the stability of the system. Specific Embodiments

[0116] The technical solution of the present invention will be described in detail below. Obviously, the described embodiments are only a part of the embodiments of the present invention, rather than all of the embodiments.

[0117] An integrated energy system load restoration method based on risk assessment and scenario matching, the steps include:

[0118] Step 1: Generate a set of typical fault scenarios under typhoon disasters based on the Monte Carlo method.

[0119] Due to the different characteristics of various natural disasters, there are significant differences in the failure probabilities of key components such as distribution network lines, heating pipelines, and gas transmission pipelines. Among them, typhoon is a typical extreme natural disaster with strong destructive power and wide influence range.

[0120] The specific steps include:

[0121] Step 1-1: Build a typhoon disaster simulation model for calculating the real-time wind speed at different positions in the affected area according to typhoon information.

[0122] The calculation process of the typhoon disaster simulation model is as follows:

[0123] Step 1-1-1: During the movement of the typhoon, the pressure difference between the typhoon center and the outside gradually decays with time t. Calculate the pressure difference between the typhoon center and the outside at time t:

[0124] ΔH(t) = 0.75ΔH 0 - 0.508[1 + sinθ]t;

[0125] Where: ΔH 0 is the initial pressure difference between the typhoon center and the outside; θ is the angle between the moving direction of the typhoon center at time t and the horizontal axis of the coordinate system.

[0126] Step 1-1-2: Calculate the gradient wind speed and the maximum wind speed radius of the typhoon at each moment according to the pressure difference between the typhoon center and the outside:

[0127]

[0128] R max (t) = exp(5.104 - 1.124ΔH(t) 0.6 ).

[0129] v gw (t) is the gradient wind speed at time t, and R max (t) is the maximum wind speed radius at time t.

[0130] Step 1-1-3: Calculate the maximum wind speed v Rmax (t) within the typhoon influence range at time t according to the gradient wind speed and the typhoon moving speed:

[0131] v Rmax (t) = 0.865v gw (t) + 0.5v T ;

[0132] Step 1-1-4: Calculate the wind speed at any position within the integrated energy system:

[0133]

[0134] Where: v ra (t) is the real-time wind speed at a certain position at time t, and d(t) is the real-time distance between the typhoon center position and this position at time t.

[0135] Step 1-2: Calculate the probability distribution of the starting time of power line faults under the influence of typhoon disasters.

[0136] Under extreme typhoon weather, power branches are extremely vulnerable to damage due to excessive wind speed. Therefore, it is necessary to determine the failure probability according to the wind speed borne by the power branches.

[0137] For the i-th power branch, according to the wind speed v i (t) borne at its location at time t, calculate the real-time probability h i (t) of it having a fault at time t in the following way:

[0138]

[0139] where v des is the designed wind speed of this power branch; is the maximum wind speed that this power branch can withstand; is a model parameter, which is set to 0.4 in this embodiment.

[0140] The probability p i,t that the start time of the fault of this power branch is time t is:

[0141]

[0142] Thus, the calculation relationship between the probability distribution of all power branches in the integrated energy system with different times as the start time of the fault and the typhoon information is obtained.

[0143] Step 1-3: Based on the Monte Carlo method and the probability distribution of the start time of the fault of each power branch, obtain a fault scenario set, and then cluster the fault scenario set to obtain a typical fault scenario set.

[0144] Step 1-3-1: Randomly construct multiple different typhoon information, where the typhoon information includes the typhoon center trajectory, the moving speed of the typhoon center, and the initial pressure difference between the typhoon center and the outside. At the same time, randomly construct multiple different system states for the integrated energy system, where the system state includes the load demands of each node at different times in the system, the output of each coupling device, and the operating state of each coupling device.

[0145] Randomly combine the constructed typhoon information and system state to obtain multiple operation scenarios.

[0146] Step 1-3-2: For each operation scenario, calculate the probability distribution of the start time of the fault of all power lines respectively, and then use the Monte Carlo simulation method to randomly generate fault samples based on this probability distribution. The fault sample refers to the fault scenario where a certain power branch starts to have a fault at a certain time under the corresponding operation scenario, so as to obtain a fault scenario set composed of multiple fault scenarios. Each fault scenario corresponds to a feature vector, and the feature vector includes the corresponding typhoon information, system state, as well as the start time of the fault and the number of the power branch where the fault occurs.

[0147] Step 1-3-3: Cluster the fault scenario set according to the feature vectors of each fault scenario. The cluster centers in each obtained cluster are the typical fault scenarios, and the obtained typical fault scenarios form the typical fault scenario set.

[0148] In this embodiment, the K-means clustering algorithm is used for data clustering. Its core idea is to iteratively optimize and assign sample points to the nearest cluster center. Specifically, first randomly initialize K cluster centers, then calculate the Euclidean distance of the feature vectors between each sample point and these centers, and assign it to the nearest cluster center. Then, update the position of the cluster center to the mean value of all sample points in the current cluster. Through multiple iterations, until the cluster center no longer changes significantly or reaches the preset number of iterations, the clustering division is finally realized. The purpose of this method is to make the sample points within the same cluster have high similarity, while the sample points between different clusters have large differences.

[0149] Step 2: Construct a comprehensive energy system load restoration model. The objective function of the comprehensive energy system load restoration model contains weights corresponding to different energy types respectively, and the constraint conditions of the comprehensive energy system load restoration model contain multi-dimensional risk constraints.

[0150] (1) Objective function.

[0151] In this model, when the failure causes insufficient energy supply, the regional comprehensive energy system will optimize the system energy supply with the goal of minimizing the total reduction of multi-energy loads. Therefore, the objective function of the comprehensive energy system load restoration model is:

[0152]

[0153] Among them, Ω L is the multi-energy load set of the regional comprehensive energy system, including electric power load, heat load and natural gas load; and are the reduction amounts of electric power, natural gas, and heat energy supply to load node i at the moment t in the next implementation of the restoration strategy in the comprehensive energy system, that is, in the formula, which are decision variables. If a load node only has the energy supply of some energy types, the reduction amounts of other types are defaulted to 0. α e 、α g and α h are the weights of the three energy types of electric power, natural gas, and heat respectively.

[0154] (2) Constraint conditions.

[0155] The constraint conditions of the comprehensive energy system load restoration model include:

[0156] (2.1) Power system constraints:

[0157] For any load node j:

[0158]

[0159] Wherein, i, j, and k are the numbers of load nodes, ij represents the power branch with the starting point at load node i and the ending point at load node j, and jk represents the power branch with the starting point at load node j and the ending point at load node k; Ω e is the set of load nodes in the topology after the distribution network fails; denotes Ω e all the power branches in which the starting point is load node j, denotes Ω e all the power branches in which the ending point is load node j; P jk,t , Q jk,t are respectively the active power transmission power and the reactive power transmission power of branch jk at time t, P ij,t , Q ij,t are respectively the active power transmission power and the reactive power transmission power of branch ij at time t, which are decision variables; I ij,t is the current of branch ij at time t; r ij , x ij are respectively the resistance and reactance of branch ij; are respectively the active and reactive power outputs of load node j as an energy supply device at time t, which are decision variables; are respectively the active and reactive power loads of load node j at time t, which are decision variables; are respectively the active and reactive power losses of load node j at time t.

[0160] For any load node j and all power branches ij with the ending point at load node j:

[0161]

[0162] For any load node i:

[0163]

[0164] Wherein, V i,t , V j,t are respectively the voltage values of load node i and load node j at time t; z ij,t is the on-off state of branch ij at time t. If the branch ij is disconnected due to a typhoon, it is represented as 0, otherwise it is 1; M is a preset constant; V 0 is the reference voltage value; are respectively the upper and lower limits of the active power transmission power of branch ij; are respectively the upper and lower limits of the reactive power transmission power of branch ij. are the upper and lower limits of the voltage of load node i, respectively.

[0165] (2.2) Dynamic constraints of the natural gas system:

[0166] For any load node i:

[0167]

[0168] In the formula, are the natural gas flow rates flowing out and flowing into pipeline ij connecting load node i and load node j at time t, respectively; is the set of all natural gas nodes connected to load node i through pipelines and with load node i as the gas outlet end; is the gas supply volume from the gas source node to load node i at time t; is the gas load of load node i at time t and is a decision variable; is the set of energy conversion devices using the natural gas of load node i as the energy source; is the gas load of energy conversion device k at time t and is a decision variable; is the reduction amount of the natural gas energy supply of load node i and is a decision variable.

[0169]

[0170] In the formula, G ij,t is the average pipeline flow rate of pipeline ij at time t, K ij is the preset pipeline parameter of pipeline ij, D ij is the natural gas flow direction of pipeline ij, μ i,t and μ j,t are the gas pressures of load node i and load node j at time t, respectively, and are known quantities.

[0171] For any pipeline ij:

[0172]

[0173] For any load node i:

[0174]

[0175] μ min ≤ μ i,t ≤ μ max ;

[0176] In the formula, are the upper and lower limits of the natural gas pipeline flow rate, G s,max and G s,min are the upper and lower limits of the gas source output, μ max and μ min are the upper and lower limits of the node gas pressure, respectively.

[0177] (2.3) Thermal system dynamic constraint:

[0178]

[0179] In the formula, are the heating powers of the electric boiler and the CHP unit in the integrated energy system at time t, respectively, which are decision variables; c p is the specific heat capacity of water; m ij,t is the sectional flow rate of the heating pipeline ij at time t; T sj,t and T ri,t are the temperatures at the end of the heating pipeline ij and the beginning of the corresponding return water pipeline at time t, respectively, which are known quantities; is the heat load loss at the node at time t, which is a decision variable.

[0180]

[0181] In the formula, is the heat load at the node at time t.

[0182] For any load node h:

[0183]

[0184] In the formula: are the heating pipeline sets with load node h as the first head and the heating pipeline sets with load node h as the first end, respectively; represent the outlet temperature of the heating part and the outlet temperature of the heat regeneration part of the heating pipeline k, respectively; are the heating temperature and the heat regeneration temperature of load node h, respectively; are the flow rate of the heating part and the flow rate of the heat regeneration part of pipeline k, respectively.

[0185] For the heating temperature T of any load node at time t s,t , the heat regeneration temperature T of any load node r,t and the sectional flow rate m of any heating pipeline ij ij,t :

[0186]

[0187] m min ≤m ij,t ≤m max ;

[0188] In the formula, are the upper and lower limits of the temperatures of each node in the water supply network, are the upper and lower limits of the temperatures of each node in the return water network, m max, m min is the upper and lower limits of the flow rate of the heating pipeline ij section.

[0189] (2.4) Output constraint of the coupling equipment:

[0190]

[0191] In the formula, is the output electric power of the gas turbine at time t, which is a decision variable; is the natural gas consumed by the gas turbine at time t, which is a decision variable; L NG is the calorific value of natural gas, which is 9.78 J / kg; η GT is the gas-to-electricity conversion efficiency of the gas turbine; θ GT is the start-stop state of the gas turbine, which is a 0-1 variable, 0 for the off state and 1 for the on state, and is a decision variable.

[0192]

[0193] In the formula, are the electric and heat powers output by the CHP unit at time t respectively, which are decision variables; θ CHP is the start-stop state of the CHP unit, which is a 0-1 variable, 0 for the off state and 1 for the on state, and is a decision variable; is the natural gas consumed by the CHP unit at time t, which is a decision variable; is the heat loss rate of the CHP unit; η GTE , η GTH are the gas-to-electricity conversion efficiency and gas-to-heat conversion efficiency of the CHP unit respectively.

[0194]

[0195] In the formula, is the heat power output by the electric boiler at time t, which is a decision variable; is the electric power consumed by the electric boiler at time t, which is a decision variable; η eB is the electric-to-heat conversion efficiency of the electric boiler; is the heat power loss rate of the electric boiler; θ EB is the start-stop state of the electric boiler, which is a 0-1 variable, 0 for the off state and 1 for the on state, and is a decision variable.

[0196] (2.5) Output constraint of the energy supply equipment, i.e., the coupling equipment:

[0197]

[0198] In the formula, P eq,t is the power output of the coupling equipment at time t, which is a decision variable; They are the upper and lower limits of the power output of the coupling device respectively.

[0199] (2.6) Multi-dimensional risk constraints:

[0200] According to the simulated real-time disaster information, accurately evaluate the overall risk status of the system. The present invention constructs a risk assessment index system for the integrated energy system from three dimensions: the system, nodes, and coupling devices.

[0201] When solving the load restoration model of the integrated energy system, each index in the multi-dimensional risk constraints must be within the preset range.

[0202] (2.6.1) System layer risk index:

[0203] Adopt the average load shedding rate of the system as the system layer restoration ability index. This index quantifies the real-time risk assessment of the system during the whole process of typhoon impact from a global perspective. By actively adjusting and reducing the load, formulate the load restoration strategy at multiple time scales. The specific calculation method is as follows:

[0204]

[0205] In the formula: R sys,t is the average load shedding rate of the system at time t; α e , α g , α h are the weights of electricity, gas, and heat energy respectively; and are the reduction amounts of electricity, gas, and heat loads of the system at time t respectively, which are decision variables; and are the total demands of electricity, gas, and heat loads at time t respectively, which are known quantities.

[0206] (2.6.2) Node layer risk index:

[0207] The movement process of the typhoon has spatio-temporal characteristics. Accurately identifying the important nodes of the system during each period of typhoon impact plays an important role in reducing the fault loss. Therefore, the present invention proposes a node importance index to identify the weak links of the system, which is quantified by the average load loss of the node in each period. For any load node node, the calculation method is as follows:

[0208]

[0209] In the formula: ΔP node,t is the reduction amount of the system to the load node node at time t, which is a decision variable; L node,t is the load demand of the load node node at time t, which is a known quantity.

[0210] (2.6.3) Coupling device layer risk index:

[0211] Multi - energy complementarity is one of the important characteristics of the integrated energy system, and the coupling device provides a conversion path for energy complementary substitution. The present invention proposes an importance index for the coupling device to quantitatively analyze the influence degree of the coupling device on the real - time load reduction rate, which is expressed as when the coupling device fails, it passively causes load reduction, resulting in the change amount of the system average load reduction rate. For any coupling device j, the calculation method is as follows:

[0212] R eq,j,t =R sys,t -X eq,j,t ;

[0213]

[0214] In the formula: X eq,j,t is the average load reduction rate of the system at time t when the coupling device j fails; and are the electricity, gas, and heat load reduction amounts of the system at time t when the coupling device j fails.

[0215] Step 3: After a disaster occurs and the energy system fails, regularly use the integrated energy system load recovery model to generate a recovery strategy and put the recovery strategy into operation for load recovery.

[0216] When generating the recovery strategy, first match the typical fault scenario closest to the current scenario from the typical fault scenario set, then optimize the weights in the integrated energy system load recovery model according to the historical data corresponding to the typical fault scenario, and then obtain the recovery strategy based on the integrated energy system load recovery model.

[0217] The specific steps for generating and putting the recovery strategy into operation are as follows:

[0218] Step 3 - 1: Collect the typhoon information, system status, fault start time, and the power branch numbers where the fault occurs in the current actual fault scenario, construct the current feature vector, then calculate the distances between the current feature vector and the feature vectors of each typical fault scenario, find the typical fault scenario corresponding to the minimum distance, and then obtain M groups of historical data of the integrated energy system during the post - disaster recovery process corresponding to the typical fault scenario. The historical data includes the reduction amounts of electricity, natural gas, and heat supply of each load node.

[0219] Step 3 - 2: Construct the weight - solving objective function:

[0220]

[0221] In the formula, α e 、α g and α hThey are the weights of the three energy types of electricity, natural gas, and heat respectively, and are the values to be decided currently; and They are the reduction amounts of electricity, natural gas, and heat energy supply of the integrated energy system to load node i at time t respectively.

[0222] Step 3-3: Set initial values for the weights α e , α g and α h , and initialize the iteration count to 0.

[0223] Step 3-4: Substitute the M groups of historical data into the objective function for solving the weights. Denote the objective function value obtained from the m-th group of historical data as obj m , calculate the average value , and increment the iteration count by 1.

[0224] Step 3-5: Determine whether the current iteration count has reached the preset value. If it has, take the current weights α e , α g , α h as the optimal weight values of this type of fault scenario set, and jump to Step 3-6. Otherwise, with the goal of minimizing the average value , update the weight values α e , α g , α h using the gradient descent method, and then jump to Step 3-4.

[0225] Step 3-6: Substitute the obtained weights α e , α g , α h into the integrated energy system load restoration model for solution to obtain the restoration strategy composed of decision variables. Then put into the restoration strategy and control the system to perform load restoration until the next moment for calculating the restoration strategy arrives.

[0226] It should be noted that for those skilled in the art, it is obvious that the present invention is not limited to the details of the above exemplary embodiments, and can be implemented in other specific forms without departing from the spirit or basic characteristics of the present invention. The scope of the present invention is defined by the claims rather than the above description.

Claims

1. A method for load restoration of an integrated energy system based on risk assessment and scenario matching, characterized in that the steps include: Step 1: Generate a typical fault scenario set under typhoon disaster based on the Monte Carlo model method; Step 2: Construct a comprehensive energy system load recovery model. The objective function of the comprehensive energy system load recovery model includes weights corresponding to different energy types, and the constraints of the comprehensive energy system load recovery model include multi-dimensional risk constraints. Step 3: After a disaster occurs and the energy system fails, the integrated energy system load recovery model is used regularly to generate a recovery strategy, and the recovery strategy is put into use to restore the load; When generating a recovery strategy, we first match the typical fault scenario closest to the current scenario from the typical fault scenario set, then optimize the weights in the integrated energy system load recovery model based on the historical data corresponding to the typical fault scenario, and then obtain the recovery strategy based on the integrated energy system load recovery model.

2. The integrated energy system load restoration method based on risk assessment and scenario matching according to claim 1, characterized in that: Step 1 specifically includes: Step 1-1, build a typhoon disaster simulation model to calculate the real-time wind speed at different locations in the disaster-stricken area based on typhoon information; Step 1-2, calculate the probability distribution of the start time of the power line fault under the influence of the typhoon disaster; Step 1-3, obtain the fault scenario set based on the Monte Carlo model method and the probability distribution of the fault start time of each power branch, and then cluster the fault scenario set to obtain the typical fault scenario set.

3. The integrated energy system load restoration method based on risk assessment and scenario matching according to claim 2 is characterized in that: In step 1-1, the calculation process of the typhoon disaster simulation model is: Step 1-1-1, calculate the pressure difference between the typhoon center and the outside at time t: ΔH(t)=0.75ΔH0-0.508[1+sinθ]t; Where: ΔH0 is the initial pressure difference between the typhoon center and the outside; θ is the angle between the moving direction of the typhoon center and the horizontal axis of the coordinate system at time t; Step 1-1-2: Calculate the typhoon's gradient wind speed and maximum wind speed radius at each moment based on the pressure difference between the typhoon center and the outside: R max (t)=exp(5.104-1.124ΔH(t) 0.6 ); v gw (t) is the gradient wind speed at time t, R max (t) is the maximum wind speed radius at time t; Step 1-1-3: Calculate the maximum wind speed v within the typhoon's influence range at time t based on the gradient wind speed and typhoon moving speed. Rmax (t): v Rmax (t)=0.865v gw (t)+0.5v T ; Step 1-1-4: Calculate the wind speed at any location in the integrated energy system: Where: v ra (t) is the real-time wind speed at a certain location at time t, and d(t) is the real-time distance between the typhoon center and the location at time t.

4. The integrated energy system load restoration method based on risk assessment and scenario matching according to claim 3 is characterized in that: In step 1-2, for the i-th power branch, according to the wind speed v at its location at time t i (t) Calculate the real-time probability h of failure at time t i (t) is as follows: In the formula, v des is the design wind speed of the power branch circuit; is the maximum wind speed that the power branch circuit can withstand; are the preset model parameters; The probability p that the power branch fault starts at time t i,t for:

5. The integrated energy system load restoration method based on risk assessment and scenario matching according to claim 2, characterized in that: Steps 1-3 include: Step 1-3-1, randomly construct multiple different typhoon information, the typhoon information includes the typhoon center trajectory, the typhoon center moving speed, and the initial pressure difference between the typhoon center and the outside; at the same time, randomly construct multiple different system states for the integrated energy system, the system state includes the load demand of each node in the system at different times, the output of each coupling device and the operating status of each coupling device; The constructed typhoon information and system status are randomly combined to obtain multiple operation scenarios; Step 1-3-2: for each operation scenario, respectively calculate the probability distribution of the start time of all power line faults, and then randomly generate fault samples based on the probability distribution using the Monte Carlo simulation method. The fault sample refers to a fault scenario in which a certain power branch starts to fail at a certain time under the corresponding operation scenario, thereby obtaining a fault scenario set consisting of multiple fault scenarios; each fault scenario corresponds to a feature vector, and the feature vector includes the corresponding typhoon information, system status, fault start time, and the number of the power branch where the fault occurs; Step 1-3-3: cluster the fault scenario set according to the feature vector of each fault scenario, the cluster center in each cluster obtained is the typical fault scenario, and the obtained typical fault scenarios constitute the typical fault scenario set.

6. The integrated energy system load restoration method based on risk assessment and scenario matching according to claim 1, characterized in that: The objective function of the integrated energy system load recovery model is: Among them, Ω L It is a collection of multi-energy loads of the regional integrated energy system, including electricity load, thermal load and natural gas load; and are the reductions in electricity, natural gas and thermal energy supply to load node i by the integrated energy system at the next time the recovery strategy is implemented, i.e., time t in the formula, which is the decision-making amount; α r , α g and α h They are the weights of the three energy types: electricity, natural gas and heat.

7. The integrated energy system load restoration method based on risk assessment and scenario matching according to claim 6, characterized in that: The constraints of the integrated energy system load recovery model include: (2.1) Power system constraints: For any load node j: Where i, j, and k are the numbers of the load nodes, ij represents the power branch with the starting point at load node i and the end point at load node j, and jk represents the power branch with the starting point at load node j and the end point at load node k; Ω e is the set of load nodes in the topology after the distribution network fails; Represents Ω e All power branches starting from load node j, Represents Ω e All power branches whose terminal is load node j; P jk,t , Q jk,t are the active transmission power and reactive transmission power of branch jk at time t, P ij,t , Q ij,t are respectively the active transmission power and reactive transmission power of branch ij at time t, which are the decision quantities; I ij,t is the current of branch ij at time t; r ij 、x ij are the resistance and reactance of branch ij respectively; are respectively the active and reactive outputs at time t when load node j is the energy supply equipment, which are the decision quantities; are respectively the active and reactive loads of load node j at time t, which are the decision quantities; are respectively the active and reactive losses of load node j at time t; For any load node j and all power branches ij whose end point is load node j: For any load node i: Where V i,t 、V j,t are the voltage values ​​of load nodes i and j at time t respectively; z ij,t is the on / off state of branch ij at time t. If the typhoon causes branch ij to be disconnected, it is expressed as 0, otherwise it is 1; M is a preset constant; V0 is the reference voltage value; are the upper and lower limits of active transmission power of branch ij respectively; They are the upper and lower limits of reactive transmission power of branch ij respectively; are the upper and lower limits of the voltage at load node i respectively; (2.2) Dynamic constraints of natural gas system: For any load node i: In the formula, are the natural gas flows outflowing and inflowing from pipeline ij connecting load node i and load node j at time t respectively; is the set of all natural gas nodes connected to load node i through pipelines and where load node i is the gas outlet; is the gas supply from the gas source node to the load node i at time t; is the gas load of load node i at time t, and is the decision quantity; is a collection of energy conversion devices that use natural gas at load node i as energy source; is the gas load of energy conversion equipment k at time t, which is the decision quantity; is the reduction amount of natural gas energy supply at load node i, which is the decision amount; In the formula, G ij,t is the average pipeline flow rate of pipeline ij at time t, K ij is the preset pipeline parameter of pipeline ij, D ij is the natural gas flow direction of pipeline ij, μ i,t , μ j,t are the air pressures at load nodes i and j at time t, respectively, which are known quantities; For any pipeline ij: For any load node i: m min ≤μ i,t ≤μ max ; In the formula, are the upper and lower limits of the natural gas pipeline flow, G s,max , G s,min are the upper and lower limits of the gas source output, μ max , μ min are the upper and lower limits of node pressure respectively; (2.3) Dynamic constraints of thermal system: In the formula, are the heating power of the electric boiler and CHP unit in the comprehensive energy system at time t, respectively, and are the decision quantities; c p is the specific heat capacity of water; m ij,t is the cross-sectional flow rate of the heating pipe ij at time t; T sj,t , T ri,t are the temperatures at the end of the heating pipe ij and the beginning of the corresponding return pipe at time t, which are known quantities; is the heat load loss of the node at time t, and is the decision quantity; In the formula, is the heat load of the node at time t; For any load node h: Where: They are respectively a set of heating pipes with load node h as the first head end and a set of heating pipes with load node h as the first end end; They are respectively represented as the outlet temperature of the heat supply part and the outlet temperature of the heat recovery part of the heat supply pipe k; are the heating temperature and reheating temperature of load node h respectively; are the flow rates of the heating part and the heat recovery part of pipeline k respectively; For any load node heating temperature T at time t s,t , the reheating temperature T of any load node r,t and the cross-sectional flow rate m of any heating pipe ij ij,t : m min ≤m ij,t ≤m max ; In the formula, are the upper and lower limits of the temperature of each node in the water supply network, is the upper and lower limits of the temperature of each node in the return water network, m max 、m min The upper and lower limits of the flow rate of the heating pipeline ij section; (2.4) Coupling equipment output constraints: In the formula, is the output power of the gas turbine at time t, which is the decision quantity; is the natural gas consumed by the gas turbine at time t, is the decision quantity; L NG is the calorific value of natural gas, which is 9.78 J / kg; η GT is the gas-to-electricity efficiency of the gas turbine; θ GT is the start / stop state of the gas turbine, which is a 0-1 variable, 0 is the shutdown state, 1 is the start state, and is the decision quantity; In the formula, are the electrical and thermal power output of the CHP unit at time t, respectively, which are the decision quantities; θ CHP is the start / stop state of the CHP unit, which is a 0-1 variable, 0 is the off state, 1 is the on state, and is the decision quantity; is the natural gas consumed by the CHP unit at time t, which is the decision quantity; is the heat loss rate of the CHP unit; η GTE , η GTH They are the gas-to-electricity conversion efficiency and gas-to-heat conversion efficiency of the CHP unit respectively; In the formula, is the thermal power output of the electric boiler at time t, which is the decision quantity; is the electric power consumed by the electric boiler at time t, which is the decision-making quantity; η EB is the electricity-to-heat efficiency of the electric boiler; is the thermal power loss rate of the electric boiler; θ EB is the start / stop state of the electric boiler, which is a 0-1 variable, 0 is the off state, 1 is the on state, and is the decision quantity; (2.5) Output constraints of energy supply equipment, i.e. coupling equipment: Where P eq,t is the power output of the coupling device at time t, and is the decision quantity; They are the upper and lower limits of the power output of the coupling device respectively.

8. The integrated energy system load restoration method based on risk assessment and scenario matching according to claim 6, characterized in that: The multidimensional risk constraints specifically include system-level risk indicators, node-level risk indicators, and coupled equipment-level risk indicators. When solving the integrated energy system load recovery model, each indicator in the multidimensional risk constraints must be within a preset range.

9. The integrated energy system load restoration method based on risk assessment and scenario matching according to claim 8, characterized in that: The calculation method of system-level risk indicators is: In the formula, R sys,t is the average load reduction rate of the system at time t; α e , α g , α h These are the weights of electricity, gas and heat respectively; and are the reductions in the electricity, gas and heat loads of the system at time t, respectively, which are the decision quantities; and are the total demands for electricity, gas and heat load at time t, which are known quantities; The calculation method of node layer risk index is: For any load node: In the formula, ΔP node,t is the reduction amount of the system on the load node node at time t, which is the decision amount; L node,t is the load demand of the load node node at time t, which is a known quantity; The calculation method of the coupling device layer risk index is: For any coupled device j: R eq,j,t =R sys,t -X eq,j,t ; Where: X eq,j,t is the average load reduction rate of the system at time t when a failure occurs in coupling device j; and is the reduction in the system's electricity, gas and heat loads at time t when coupling device j fails.

10. The integrated energy system load restoration method based on risk assessment and scenario matching according to any one of claims 1 to 9, characterized in that: In step 3, the specific steps to generate and implement the recovery strategy are: Step 3-1, collect typhoon information, system status, fault start time and fault branch number in the current actual fault scenario, construct a current feature vector, then calculate the distance between the current feature vector and the feature vector of each typical fault scenario, find the typical fault scenario corresponding to the minimum distance, and then obtain M groups of historical data of the comprehensive energy system in the post-disaster recovery process of the typical fault scenario, wherein the historical data includes the reduction amount of electricity, natural gas and thermal energy supply of each load node; Step 3-2: Construct weighted solution objective function: In the formula, α e , α g and α h are the weights of the three energy types of electricity, natural gas, and heat, which are the current values ​​to be decided; and are the reductions in electricity, natural gas and thermal energy supply to load node i by the integrated energy system at time t, respectively; Step 3-3, weight α e , α g and α h Set the initial value and initialize the number of iterations to 0; Step 3-4, substitute the M groups of historical data into the weighted objective function, and record the objective function value obtained from the mth group of historical data as obj m , find the average value And increase the number of iterations by 1; Step 3-5: Determine whether the current number of iterations has reached the preset value. If it has reached the preset value, the current weight α e , α g , α h As the optimal weight value of this type of fault scenario set, jump to steps 3-6, otherwise use the average value Minimum is the goal, and the weight value α is updated through the gradient descent method. e , α g , α h , then jump to step 3-4; Step 3-6: The obtained weight α e , α g , α h Substitute it into the integrated energy system load recovery model for solution, and obtain the recovery strategy composed of decision quantities; Then the recovery strategy is put into use and the control system performs load recovery until the next time the recovery strategy is calculated.

Images