CN111898816A - Dynamic state estimation method for comprehensive energy system - Google Patents

Dynamic state estimation method for comprehensive energy system Download PDF

Info

Publication number
CN111898816A
CN111898816A CN202010704522.5A CN202010704522A CN111898816A CN 111898816 A CN111898816 A CN 111898816A CN 202010704522 A CN202010704522 A CN 202010704522A CN 111898816 A CN111898816 A CN 111898816A
Authority
CN
China
Prior art keywords
node
pipeline
network
matrix
gas
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Pending
Application number
CN202010704522.5A
Other languages
Chinese (zh)
Inventor
陈亮
贺伟
丁宇
李扬
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Nanjing University of Information Science and Technology
Original Assignee
Nanjing University of Information Science and Technology
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by Nanjing University of Information Science and Technology filed Critical Nanjing University of Information Science and Technology
Priority to CN202010704522.5A priority Critical patent/CN111898816A/en
Publication of CN111898816A publication Critical patent/CN111898816A/en
Pending legal-status Critical Current

Links

Images

Classifications

    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06Q—INFORMATION AND COMMUNICATION TECHNOLOGY [ICT] SPECIALLY ADAPTED FOR ADMINISTRATIVE, COMMERCIAL, FINANCIAL, MANAGERIAL OR SUPERVISORY PURPOSES; SYSTEMS OR METHODS SPECIALLY ADAPTED FOR ADMINISTRATIVE, COMMERCIAL, FINANCIAL, MANAGERIAL OR SUPERVISORY PURPOSES, NOT OTHERWISE PROVIDED FOR
    • G06Q10/00—Administration; Management
    • G06Q10/04—Forecasting or optimisation specially adapted for administrative or management purposes, e.g. linear programming or "cutting stock problem"
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06Q—INFORMATION AND COMMUNICATION TECHNOLOGY [ICT] SPECIALLY ADAPTED FOR ADMINISTRATIVE, COMMERCIAL, FINANCIAL, MANAGERIAL OR SUPERVISORY PURPOSES; SYSTEMS OR METHODS SPECIALLY ADAPTED FOR ADMINISTRATIVE, COMMERCIAL, FINANCIAL, MANAGERIAL OR SUPERVISORY PURPOSES, NOT OTHERWISE PROVIDED FOR
    • G06Q10/00—Administration; Management
    • G06Q10/06—Resources, workflows, human or project management; Enterprise or organisation planning; Enterprise or organisation modelling
    • G06Q10/063—Operations research, analysis or management
    • G06Q10/0639—Performance analysis of employees; Performance analysis of enterprise or organisation operations
    • G06Q10/06393—Score-carding, benchmarking or key performance indicator [KPI] analysis
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06Q—INFORMATION AND COMMUNICATION TECHNOLOGY [ICT] SPECIALLY ADAPTED FOR ADMINISTRATIVE, COMMERCIAL, FINANCIAL, MANAGERIAL OR SUPERVISORY PURPOSES; SYSTEMS OR METHODS SPECIALLY ADAPTED FOR ADMINISTRATIVE, COMMERCIAL, FINANCIAL, MANAGERIAL OR SUPERVISORY PURPOSES, NOT OTHERWISE PROVIDED FOR
    • G06Q50/00—Systems or methods specially adapted for specific business sectors, e.g. utilities or tourism
    • G06Q50/06—Electricity, gas or water supply

Abstract

The invention discloses a dynamic state estimation method of a comprehensive energy system, which comprises the following steps: defining a specific state quantity form, determining the measurement quantity of the comprehensive energy system, calculating a step length, and forming a forecast error variance matrix and a measurement error variance matrix; forming a system state transition matrix and a measurement matrix according to the topological structure and the element parameters of the comprehensive energy system and the power balance relationship between the gas turbine unit and the cogeneration unit; forecasting the state quantity according to the state transition matrix to obtaintThe forecast value of the state quantity at +1 moment and a forecast error variance matrix; form atMeasuring phasor at +1 moment, calculating filter gain, and calculatingtState quantity estimated value and estimation error variance matrix at +1 moment; and judging whether the calculation is finished or not, outputting an estimation result, and otherwise, estimating at the next moment. The invention can not only be applied to a comprehensive energy systemThe dynamic process is forecasted, and meanwhile, the state quantity can be estimated, so that accurate dynamic information is obtained.

Description

Dynamic state estimation method for comprehensive energy system
Technical Field
The invention relates to a dynamic state estimation method for an integrated energy system, and belongs to the technical field of integrated energy system monitoring.
Background
With the increasing environmental problem, the scale of new energy utilization such as wind power and photovoltaic is continuously enlarged. Because the new energy has stronger random fluctuation, the safe and stable operation of the system is influenced to a certain extent by accessing the large-scale new energy into the power grid. The flexible turndown capability of the gas turbine and the complementary nature of the district heating grid and the electric power system through the cogeneration unit provide a new means of control over the operation of the electric power system. Therefore, the problem of monitoring the operation state of the comprehensive energy system consisting of the power system, the natural gas network and the centralized heating network is particularly important. Can monitor comprehensive energy system running state through data acquisition device and communication system, however, random error and bad data appear in the data are inevitable in collection and transmission process, consequently, need carry out state estimation to the data of gathering and just can use.
The existing state estimation method of the comprehensive energy system is mainly static state estimation, namely state estimation aiming at the condition that each state quantity does not change when the system is in stable operation. The methods cannot acquire the dynamic process information of the comprehensive energy system and do not have the state forecasting function. In actual operation, due to the constant changes of power supply and load, the compression of gas in the natural gas pipeline and the heat transfer characteristic of the centralized heating network, the integrated energy system is usually in a dynamic process which changes constantly. The real-time accurate acquisition of the dynamic process information and the prediction are more important for the safe and stable operation and the optimal control of the comprehensive energy system.
The comprehensive energy system dynamics composed of an electric power system, a natural gas network and a centralized heating network is a multi-time scale dynamic process formed by the interaction of various physical characteristic dynamic behaviors. The electromechanical transient process of the power system is a dynamic process of the change of the rotating speed of a generator rotor caused by the imbalance of the mechanical power and the electromagnetic power of the generator, and the duration time is from several seconds to tens of seconds; the dynamic change process of the power flow of the power system caused by the load change is usually from several minutes to several hours; the variations in gas pressure and flow caused by the compression of natural gas in pipelines are typically of the order of hours; the dynamic process of the central heating network caused by the change of the heating load usually lasts for hours, and the flow of water in the pipeline is usually constant due to the regulating action of the controller.
Disclosure of Invention
The invention aims to overcome the defects of the prior art and provides a method for estimating the dynamic state of an integrated energy system. The invention adopts the Kalman filtering algorithm to carry out dynamic state estimation on the comprehensive energy system. For the electric power system, the node voltage phasor is used as a state quantity, and an exponential smoothing model is adopted for forecasting; for a natural gas network, gas pressure and flow are used as state quantities, and algebraic equations and partial differential equations meeting mass balance and momentum conservation constraints are adopted to forecast the state quantities; for the central heating network, the temperature of water in the pipeline is used as a state quantity, and a partial differential equation is adopted for forecasting. The power system and the natural gas network are coupled through a gas turbine, and the power system and the centralized heating network are coupled through a cogeneration unit. The constructed dynamic state estimation method not only can forecast the dynamic process of the comprehensive energy system, but also can estimate the state quantity, thereby obtaining accurate dynamic information.
The invention specifically adopts the following technical scheme to solve the technical problems:
a dynamic state estimation method of an integrated energy system comprises the following steps:
defining a specific form of a state quantity x, determining a measurement quantity z of the comprehensive energy system, calculating a step length delta t, and forming a forecast error variance matrix Q and a measurement error variance matrix R;
forming a system state transition matrix F according to the topological structure and element parameters of the comprehensive energy system and the power balance relationship of the gas turbine unit and the cogeneration unittAnd a measurement matrix Ht；
According to the state transition matrix FtForecasting the state quantity to obtain a forecast value of the state quantity x at the t +1 moment
Figure BDA0002594174370000021
And forecast error moment of varianceArray Pt+1|t；
Measuring phasor z at time t +1t+1Calculating a filter gain Kt+1Calculating the state quantity estimated value at the t +1 moment
Figure BDA0002594174370000022
And estimate the error variance matrix Pt+1|t+1；
Judging whether the calculation is finished or not, if so, outputting an estimation result, otherwise, enabling t to be t +1 to form a system state transition matrix Ft+1And a measurement matrix Ht+1The next time estimation is performed.
Further, as a preferred technical solution of the present invention, the state quantity x defined in the method is in a specific form:
Figure BDA0002594174370000023
wherein xEFor real part e of voltage of power grid bus ppAnd imaginary part fpConstructed power system state variables, i.e. xE＝[e1f1e2f2… epfp…]T，
Figure BDA0002594174370000024
nBThe number of the buses is; x is the number ofGFor the gas density rho of a node in a natural gas pipelinekAnd the gas flow from the node i end to the node j end in the pipeline ij
Figure BDA0002594174370000025
Gas flow from node j end to node i end in pipeline ji
Figure BDA0002594174370000026
Constructed natural gas network state variables, i.e.
Figure BDA0002594174370000028
nNAnd nLRespectively being natural gas network nodesNumber of points and number of pipes; x is the number ofHFor central heating network node temperature ThConstructed heat network state quantity, i.e. xH＝[T1,T2,…，Th，…]T，
Figure BDA0002594174370000029
nHThe number of nodes of the heat supply network;
Figure BDA00025941743700000210
for the time dimension, superscript represents the real number spatial dimension;
the specific form of the measurement quantity z determined in the method is as follows:
Figure BDA0002594174370000031
wherein
Figure BDA0002594174370000032
Figure BDA0002594174370000033
zH＝[Tz1,Tz2,…,Tzl,…]T；zVAnd zIMeasuring vectors formed by real parts and imaginary parts of the power grid bus voltage, branch current and real parts and imaginary parts of node injection current respectively; z is a radical ofPAnd zMRespectively measuring vectors formed by measuring the pressure intensity of a natural gas network node and measuring the gas flow at two ends of a pipeline; z is a radical ofHMeasuring vector formed by measuring the temperature of the central heating network node; t iszlMeasuring the temperature of the node l of the centralized heating network.
Further, as a preferred technical solution of the present invention, the system state transition matrix FtThe system is determined by a dynamic equation, a topological relation and element parameters which are met by the state quantities of a power system, a natural gas pipe network and a central heating network;
for the power system, the state quantity is forecasted by the following exponential smoothing algorithm:
Lt＝αSxEt+(1-αS)(Lt-1+Tt-1)
Tt＝βS(Lt-Lt-1)+(1-βS)Tt-1
xEt+1＝Lt+Tt
in the formula, LtAnd TtRespectively a smooth horizontal component and a trend component; l ist-1And Tt-1Respectively a smooth horizontal component and a trend component at the time of t-1; x is the number ofEtAnd xEt+1Respectively obtaining exponential smoothing prediction values of state quantities at the time t and the time t + 1; alpha is alphaSAnd betaSRespectively are smoothing parameters;
and (3) establishing a natural gas network model: taking the flow at two ends of each pipeline as state quantity, wherein the direction is from the small node number to the large node number as positive; each pipeline dynamic equation is determined by the following partial differential equation:
Figure BDA0002594174370000034
Figure BDA0002594174370000035
in the formula: tau and zeta are the spatial distance in time and pipeline direction respectively;
Figure BDA0002594174370000036
is the flow rate of gas in the pipeline; d and a are the diameter and the sectional area of the pipeline respectively; c. C2Is the speed of sound; gamma is the coefficient of friction;
Figure BDA0002594174370000037
is the average flow velocity of the gas in the pipeline; in addition, the gas in the natural gas pipeline network meets the mass balance principle at the node, so the state quantity of the partial differential equation also needs to meet the following boundary conditions:
Figure BDA0002594174370000038
in the formula:
Figure BDA0002594174370000039
and
Figure BDA00025941743700000310
the gas flow of the i end in the pipeline ik at the time t and the gas flow of the outflow node i are respectively;
Figure BDA00025941743700000311
the gas flow from the node i to the node j in the pipeline ij at the moment t; k, j belongs to i and represents a node k and j are connected with the node i;
discretizing the partial differential equation by a finite element method:
Figure BDA0002594174370000041
Figure BDA0002594174370000042
Figure BDA0002594174370000043
in the formula, LPijIs the length of conduit ij; rhoj,t+1And ρj,tThe node j gas densities at the t +1 moment and the t moment respectively; rhoi,t+1And ρi,tThe node i gas densities at the t +1 moment and the t moment respectively;
Figure BDA0002594174370000044
and
Figure BDA0002594174370000045
the gas flow from the node j to the node i in the pipeline ji at the time t +1 and the time t respectively;
Figure BDA0002594174370000046
the gas flow from the node i to the node j in the pipeline ij at the moment t + 1; dijAnd aijThe diameter and the cross-sectional area of the pipeline ij are respectively;
for a centralized heating network, the flow rate of hot water in a pipeline is usually controlled to be constant, and only the temperature dynamics is considered; the temperatures of the water at the various nodes in the heating network are determined by the following partial differential equations:
Figure BDA0002594174370000047
in the formula: t isHIs the temperature of water in the heat supply network; zetaHThe space distance in the direction of the heat supply pipeline; v. ofHThe flow rate of hot water in the pipeline;
discretizing partial differential equations of the natural gas network and the centralized heat supply network by a finite element method to obtain a system state transition matrix Ft：
Figure BDA0002594174370000048
In the formula: i is an identity matrix; fGAnd FHRespectively obtaining state transition matrixes after discretization of partial differential equations of a natural gas network and a centralized heat supply network;
the power balance relationship of the gas turbine set is as follows:
Figure BDA0002594174370000049
Figure BDA00025941743700000410
in the formula: e.g. of the typept、eqtAnd fpt、fqtThe real part and the imaginary part of the bus voltage at the time t are respectively, and subscripts p and q are bus numbers;
Figure BDA00025941743700000411
the natural gas flow out of node s at time t; pGp,tThe output power of the gas turbine generator connected to the bus p at the moment t; gYpqAnd BYpqRespectively are the p row and q column elements of the node admittance matrix;
Figure BDA0002594174370000051
represents a bus p connecting the bus q; etapThe energy conversion efficiency of the gas turbine generator connected to the bus p;
Figure BDA0002594174370000052
a gas turbine generator representing a bus p is connected to the s-th node of the natural gas pipeline network.
The power balance relation of the cogeneration unit is as follows:
Figure BDA0002594174370000053
wherein: pHThe heat power is adopted; c is the specific heat of water; delta TcIs the temperature difference of the circulating water of the cogeneration unit.
Further, as a preferred embodiment of the present invention, in the method, the measurement matrix H istThe concrete form is as follows:
Figure BDA0002594174370000054
in the formula: hE、HGAnd HHThe measurement matrixes are respectively a power system, a natural gas network and a centralized heating network.
Further, as a preferred embodiment of the present invention, in the method, a predicted value of the state quantity x is obtained
Figure BDA0002594174370000055
And forecast error variance matrix Pt+1|tThe specific calculation method comprises the following steps:
Figure BDA0002594174370000056
Pt+1|t＝FtPt|tFT+Q。
in the formula:
Figure BDA0002594174370000057
is the estimated value of the state quantity at the time t; pt|tIs the estimated error variance matrix at time t.
Further, as a preferred technical solution of the present invention, in the method, a filter gain K is obtainedt+1State quantity estimated value at time t +1
Figure BDA0002594174370000058
And estimate the error variance matrix Pt+1|t+1The specific calculation method comprises the following steps:
Figure BDA0002594174370000059
Figure BDA00025941743700000510
Pt+1|t+1＝Pt+1|t-Kk+1HtPt+1|t。
by adopting the technical scheme, the invention can produce the following technical effects:
1. the comprehensive energy dynamic state estimation method provided by the invention can accurately forecast and estimate the dynamic process of the system and obtain accurate system dynamic information.
2. The comprehensive energy dynamic state estimation method provided by the invention can couple the power system, the natural gas network and the centralized heat supply network through the gas turbine generator set and the cogeneration unit, thereby further improving the estimation precision.
3. The dynamic state estimation method of the comprehensive energy system selects a proper state transition equation to predict the state quantity according to the time scale characteristic of the dynamic process, and ensures the close coupling relation of the dynamic process among different systems.
4. The system process equation and the measurement equation of the comprehensive energy dynamic state estimation method are linear equations, and the calculation efficiency is higher compared with a nonlinear method.
Drawings
Fig. 1 is a schematic diagram of a natural gas pipeline model provided by the invention.
Detailed Description
The following describes embodiments of the present invention with reference to the drawings.
The invention provides a dynamic state estimation method of an integrated energy system, which specifically comprises the following steps:
step 1, defining a specific form of a state quantity x, determining a measurement quantity z of a comprehensive energy system, calculating a step length delta t, and forming a forecast error variance matrix Q and a measurement error variance matrix R;
and 2, initializing variables. Setting the iteration time t to be 0; and setting a forecast error variance and a measurement error variance of the initial value of the state quantity.
Step 3, forming a system state transition matrix F according to the topological structure and the element parameters of the comprehensive energy system and the power balance relation of the gas turbine unit and the cogeneration unittAnd a measurement matrix Ht；
Step 4, according to the state transition matrix FtForecasting the state quantity to obtain a forecast value of the state quantity x at the t +1 moment
Figure BDA0002594174370000061
And forecast error variance matrix Pt+1|t；
Step 5. forming t +1 moment to measure phasor zt+1Calculating a filter gain Kt+1Calculating the state quantity estimated value at the t +1 moment
Figure BDA0002594174370000062
And estimate the error variance matrix Pt+1|t+1。
Step 6, judging whether the calculation is finished or not, if so, executing the step 7, otherwise, turning to the step 3 to form a system state transition matrix Ft+1And a measurement matrix Ht+1Estimating the next moment;
and 7, finishing calculation and outputting an estimation result.
In step 1, the state quantity x defined by the method is in a specific form:
Figure BDA0002594174370000063
wherein x isEFor real part e of voltage of power grid bus ppAnd imaginary part fpConstructed power system state variables, i.e. xE＝[e1f1e2f2…epfp…]T，
Figure BDA0002594174370000064
nBThe number of the buses is; x is the number ofGFor the gas density rho of a node in a natural gas pipelinekAnd the gas flow from the node i end to the node j end in the pipeline ij
Figure BDA0002594174370000065
Gas flow from node j end to node i end in pipeline ji
Figure BDA0002594174370000066
Constructed natural gas network state variables, i.e.
Figure BDA0002594174370000067
nNAnd nLThe number of the natural gas network nodes and the number of the pipelines are respectively; x is the number ofHFor central heating network node temperature ThConstructed heat network state quantity, i.e. xH＝[T1,T2,…，Th，…]T，
Figure BDA0002594174370000071
nHThe number of nodes of the heat supply network;
Figure BDA0002594174370000072
for the time dimension, superscript denotes the real number spatial dimension.
In step 1, the specific form of the measurement quantity z determined by the method is
Figure BDA0002594174370000073
Wherein the content of the first and second substances,
Figure BDA0002594174370000074
Figure BDA0002594174370000075
zH＝[Tz1,Tz2,…,Tzl,…]T；zVand zIMeasuring vectors formed by real parts and imaginary parts of the power grid bus voltage, branch current and real parts and imaginary parts of node injection current respectively; z is a radical ofPAnd zMRespectively measuring vectors formed by measuring the pressure intensity of a natural gas network node and measuring the gas flow at two ends of a pipeline; z is a radical ofHMeasuring vector formed by measuring the temperature of the central heating network node; t iszlMeasuring the temperature of the node l of the centralized heating network. The calculation step length delta t is 600s, and Q and R are diagonal matrixes formed by taking the forecast error variance and the measurement error variance as diagonal elements respectively.
Further, in the step 3, the system state transition matrix FtThe system is determined by a dynamic equation, a topological relation and element parameters which are satisfied by the state quantities of the power system, the natural gas pipe network and the central heating network.
For the power system, the state quantity is forecasted by the following exponential smoothing algorithm:
Lt＝αSxEt+(1-αS)(Lt-1+Tt-1)
Tt＝βS(Lt-Lt-1)+(1-βS)Tt-1
xEt+1＝Lt+Tt
in the formula, LtAnd TtRespectively a smooth horizontal component and a trend component; l ist-1And Tt-1Respectively a smooth horizontal component and a trend component at the time of t-1; x is the number ofEtAnd xEt+1Respectively obtaining exponential smoothing prediction values of state quantities at the time t and the time t + 1; alpha is alphaSAnd betaSRespectively are smoothing parameters; the subscript t denotes the time of day.
The natural gas network model was established as shown in figure 1. The flow at the two ends of each pipeline is used as a state quantity, and the direction is from the small node number to the large node number. Each pipeline dynamic equation is determined by the following partial differential equation:
Figure BDA0002594174370000076
Figure BDA0002594174370000077
in the formula: tau and zeta are the spatial distance in time and pipeline direction respectively;
Figure BDA0002594174370000078
is the flow rate of gas in the pipeline; d and a are the diameter and the sectional area of the pipeline respectively; c. C2Is the speed of sound; gamma is the coefficient of friction;
Figure BDA0002594174370000079
is the average flow rate of gas in the pipe. In addition, the gas in the natural gas pipeline network should meet the mass balance principle at the node, so the state quantity of the partial differential equation also needs to meet the following boundary conditions:
Figure BDA00025941743700000710
in the formula:
Figure BDA0002594174370000081
and
Figure BDA0002594174370000082
the gas flow of the i end in the pipeline ik at the time t and the gas flow of the outflow node i are respectively;
Figure BDA0002594174370000083
the gas flow from the node i to the node j in the pipeline ij at the moment t; k, j ∈ i denotes that nodes k and j connect node i.
Discretizing the partial differential equation by a finite element method:
Figure BDA0002594174370000084
Figure BDA0002594174370000085
Figure BDA0002594174370000086
in the formula, LPijIs the length of the conduit ij. Rhoj,t+1And ρj,tThe node j gas densities at the t +1 moment and the t moment respectively; rhoi,t+1And ρi,tThe node i gas densities at the t +1 moment and the t moment respectively;
Figure BDA0002594174370000087
and
Figure BDA0002594174370000088
the gas flow from the node j to the node i in the pipeline ji at the time t +1 and the time t respectively;
Figure BDA0002594174370000089
the gas flow from the node i to the node j in the pipeline ij at the moment t + 1; dijAnd aijRespectively the diameter and the cross-sectional area of the pipe ij.
For central heating networks, the hot water flow rate in the pipe is usually controlled to be constant, taking into account only the temperature dynamics. Neglecting heat losses, the temperature of the water at the various nodes in the heating network is determined by the following partial differential equation:
Figure BDA00025941743700000810
in the formula: t isHIs the temperature of water in the heat supply network; zetaHThe space distance in the direction of the heat supply pipeline; v. ofHIs the hot water flow rate in the pipe.
Discretizing partial differential equations of the natural gas network and the centralized heat supply network by a finite element method to obtain a system state transition matrix Ft：
Figure BDA00025941743700000811
In the formula: i is an identity matrix; fGAnd FHRespectively are state transition matrixes obtained after discretization of partial differential equations of a natural gas network and a centralized heat supply network.
The power balance relationship of the gas turbine set is as follows:
Figure BDA00025941743700000812
Figure BDA00025941743700000813
in the formula: e.g. of the typept、eqtAnd fpt、fqtThe real part and the imaginary part of the bus node voltage at the time t are respectively, and subscripts p and q are bus numbers;
Figure BDA0002594174370000091
the natural gas flow out of node s at time t; pGp,tThe output power connected to the bus p at time t; gYpqAnd BYpqRespectively are the p row and q column elements of the node admittance matrix;
Figure BDA0002594174370000092
represents a bus p connecting the bus q; etapThe energy conversion efficiency of the gas turbine generator connected to the bus p;
Figure BDA0002594174370000093
showing the pth gas turbine connected at the s-th node of the natural gas pipeline network.
The power balance relation of the cogeneration unit is as follows:
Figure BDA0002594174370000094
wherein: pHThe heat power is adopted; c is the specific heat of water; delta TcFor temperature difference of circulating water of cogeneration unit。
Further, the measurement matrix H formed in the step 3tThe concrete form is as follows:
Figure BDA0002594174370000095
in the formula: hE、HGAnd HHThe measurement matrixes are respectively a power system, a natural gas network and a centralized heating network.
Further, in the step 4, the predicted value
Figure BDA0002594174370000096
And forecast error variance matrix Pt+1|tThe specific calculation method comprises the following steps:
Figure BDA0002594174370000097
Pt+1|t＝FtPt|tFT+Q
in the formula:
Figure BDA0002594174370000098
is the estimated value of the state quantity at the time t; pt|tIs the estimated error variance matrix at time t.
Further, in the step 5, the filter gain Kt+1State quantity estimated value at time t +1
Figure BDA0002594174370000099
And estimate the error variance matrix Pt+1|t+1The specific calculation method comprises the following steps:
Figure BDA00025941743700000910
Figure BDA00025941743700000911
Pt+1|t+1＝Pt+1|t-Kk+1HtPt+1|t
in summary, the method of the present invention adopts the kalman filter algorithm to perform dynamic state estimation on the integrated energy system. For the electric power system, the node voltage phasor is used as a state quantity, and an exponential smoothing model is adopted for forecasting; for a natural gas network, gas pressure and flow are used as state quantities, and algebraic equations and partial differential equations meeting mass balance and momentum conservation constraints are adopted to forecast the state quantities; for the central heating network, the temperature of water in the pipeline is used as a state quantity, and a partial differential equation is adopted for forecasting. The power system and the natural gas network are coupled through a gas turbine, and the power system and the centralized heating network are coupled through a cogeneration unit. The constructed dynamic state estimation method not only can forecast the dynamic process of the comprehensive energy system, but also can estimate the state quantity, thereby obtaining accurate dynamic information.
The embodiments of the present invention have been described in detail with reference to the drawings, but the present invention is not limited to the above embodiments, and various changes can be made within the knowledge of those skilled in the art without departing from the gist of the present invention.

Claims (6)

1. A method for estimating the dynamic state of an integrated energy system is characterized by comprising the following steps:
defining a specific form of a state quantity x, determining a measurement quantity z of the comprehensive energy system, calculating a step length delta t, and forming a forecast error variance matrix Q and a measurement error variance matrix R;
forming a system state transition matrix F according to the topological structure and element parameters of the comprehensive energy system and the power balance relationship of the gas turbine unit and the cogeneration unittAnd a measurement matrix Ht；
According to the state transition matrix FtForecasting the state quantity to obtain a forecast value of the state quantity x at the t +1 moment
Figure FDA0002594174360000011
And forecast error variance matrix Pt+1|t；
Measuring phasor z at time t +1t+1Calculating a filter gain Kt+1Calculating the state quantity estimated value at the t +1 moment
Figure FDA0002594174360000012
And estimate the error variance matrix Pt+1|t+1；
Judging whether the calculation is finished or not, if so, outputting an estimation result, otherwise, enabling t to be t +1 to form a system state transition matrix Ft+1And a measurement matrix Ht+1The next time estimation is performed.
2. The method according to claim 1, wherein the method further comprises: the state quantity x defined in the method is in a specific form:
Figure FDA0002594174360000013
wherein xEFor real part e of voltage of power grid bus ppAnd imaginary part fpConstructed power system state variables, i.e. xE＝[e1f1e2f2… epfp…]T，
Figure FDA0002594174360000014
nBThe number of the buses is; x is the number ofGFor the gas density rho of a node in a natural gas pipelinekAnd the gas flow from the node i end to the node j end in the pipeline ij
Figure FDA0002594174360000015
Gas flow from node j end to node i end in pipeline ji
Figure FDA0002594174360000016
Constructed natural gas network state variables, i.e.
Figure FDA0002594174360000017
Figure FDA0002594174360000018
nNAnd nLThe number of the natural gas network nodes and the number of the pipelines are respectively; x is the number ofHFor central heating network node temperature ThConstructed heat network state quantity, i.e. xH＝[T1,T2,…，Th，…]T，
Figure FDA0002594174360000019
nHThe number of nodes of the heat supply network;
Figure FDA00025941743600000110
for the time dimension, superscript represents the real number spatial dimension;
the specific form of the measurement quantity z determined in the method is as follows:
Figure FDA00025941743600000111
wherein
Figure FDA00025941743600000112
zG＝
Figure FDA00025941743600000113
Figure FDA00025941743600000114
zVAnd zIMeasuring vectors formed by real parts and imaginary parts of the power grid bus voltage, branch current and real parts and imaginary parts of node injection current respectively; z is a radical ofPAnd zMRespectively measuring vectors formed by measuring the pressure intensity of a natural gas network node and measuring the gas flow at two ends of a pipeline; z is a radical ofHMeasuring vector formed by measuring the temperature of the central heating network node; t iszlMeasuring the temperature of the node l of the centralized heating network.
3. The method according to claim 1, wherein the method further comprises: the system state transition matrix FtIs generated by electricityThe system, the natural gas pipe network and the centralized heating network meet the state quantity, and the topological relation and the element parameters are jointly determined;
for the power system, the state quantity is forecasted by the following exponential smoothing algorithm:
Lt＝αSxEt+(1-αS)(Lt-1+Tt-1)
Tt＝βS(Lt-Lt-1)+(1-βS)Tt-1
xEt+1＝Lt+Tt
in the formula, LtAnd TtRespectively a smooth horizontal component and a trend component; l ist-1And Tt-1Respectively a smooth horizontal component and a trend component at the time of t-1; x is the number ofEtAnd xEt+1Respectively obtaining exponential smoothing prediction values of state quantities at the time t and the time t + 1; alpha is alphaSAnd betaSRespectively are smoothing parameters;
and (3) establishing a natural gas network model: taking the flow at two ends of each pipeline as state quantity, wherein the direction is from the small node number to the large node number as positive; each pipeline dynamic equation is determined by the following partial differential equation:
Figure FDA0002594174360000021
Figure FDA0002594174360000022
in the formula: tau and zeta are the spatial distance in time and pipeline direction respectively;
Figure FDA0002594174360000023
is the flow rate of gas in the pipeline; d and a are the diameter and the sectional area of the pipeline respectively; c. C2Is the speed of sound; gamma is the coefficient of friction;
Figure FDA0002594174360000024
is the average flow velocity of gas in the pipeline(ii) a In addition, the gas in the natural gas pipeline network meets the mass balance principle at the node, so the state quantity of the partial differential equation also needs to meet the following boundary conditions:
Figure FDA0002594174360000025
in the formula:
Figure FDA0002594174360000026
and
Figure FDA0002594174360000027
the gas flow of the i end in the pipeline ik at the time t and the gas flow of the outflow node i are respectively;
Figure FDA0002594174360000028
the gas flow from the node i to the node j in the pipeline ij at the moment t; k, j belongs to i and represents a node k and j are connected with the node i;
discretizing the partial differential equation by a finite element method:
Figure FDA0002594174360000029
Figure FDA00025941743600000210
in the formula, LPijIs the length of conduit ij; rhoj,t+1And ρj,tThe node j gas densities at the t +1 moment and the t moment respectively; rhoi,t+1And ρi,tThe node i gas densities at the t +1 moment and the t moment respectively;
Figure FDA00025941743600000211
and
Figure FDA00025941743600000212
from node j to node i in pipeline ji at time t +1 and time t, respectivelyThe gas flow rate;
Figure FDA0002594174360000031
the gas flow from the node i to the node j in the pipeline ij at the moment t + 1; dijAnd aijThe diameter and the cross-sectional area of the pipeline ij are respectively;
for a centralized heating network, the flow rate of hot water in a pipeline is usually controlled to be constant, and only the temperature dynamics is considered; the temperatures of the water at the various nodes in the heating network are determined by the following partial differential equations:
Figure FDA0002594174360000032
in the formula: t isHIs the temperature of water in the heat supply network; zetaHThe space distance in the direction of the heat supply pipeline; v. ofHThe flow rate of hot water in the pipeline;
discretizing partial differential equations of the natural gas network and the centralized heat supply network by a finite element method to obtain a system state transition matrix Ft：
Figure FDA0002594174360000033
In the formula: i is an identity matrix; fGAnd FHRespectively obtaining state transition matrixes after discretization of partial differential equations of a natural gas network and a centralized heat supply network;
the power balance relationship of the gas turbine set is as follows:
Figure FDA0002594174360000034
Figure FDA0002594174360000035
in the formula: e.g. of the typept、eqtAnd fpt、fqtThe real part and the imaginary part of the bus voltage at the time t are respectively, and subscripts p and q are bus numbers;
Figure FDA0002594174360000036
the natural gas flow out of node s at time t; pGp,tThe output power of the gas turbine generator connected to the bus p at the moment t; gYpqAnd BYpqRespectively are the p row and q column elements of the node admittance matrix;
Figure FDA0002594174360000037
represents a bus p connecting the bus q; etapThe energy conversion efficiency of the gas turbine generator connected to the bus p; p < s represents that a gas turbine generator of a bus p is connected to the s-th node of the natural gas pipe network;
the power balance relation of the cogeneration unit is as follows:
Figure FDA0002594174360000038
wherein: pHThe heat power is adopted; c is the specific heat of water; delta TcIs the temperature difference of the circulating water of the cogeneration unit.
4. The method according to claim 1, wherein the method further comprises: in which method a measurement matrix H is usedtThe concrete form is as follows:
Figure FDA0002594174360000041
in the formula: hE、HGAnd HHThe measurement matrixes are respectively a power system, a natural gas network and a centralized heating network.
5. The method according to claim 1, wherein the method further comprises: predicted values of the state quantities x in the method
Figure FDA0002594174360000042
And forecast error variance matrix Pt+1|tThe specific calculation method comprises the following steps:
Figure FDA0002594174360000043
Pt+1|t＝FtPt|tFt T+Q
in the formula:
Figure FDA0002594174360000044
is the estimated value of the state quantity at the time t; pt|tIs the estimated error variance matrix at time t.
6. The method according to claim 1, wherein the method further comprises: in the method, the filter gain Kt+1State quantity estimated value at time t +1
Figure FDA0002594174360000045
And estimate the error variance matrix Pt+1|t+1The specific calculation method comprises the following steps:
Kt+1＝Pt+1|tHt T(HtPt+1|tHt T+R)-1
Figure FDA0002594174360000046
Pt+1|t+1＝Pt+1|t-Kk+1HtPt+1|t。
CN202010704522.5A 2020-07-21 2020-07-21 Dynamic state estimation method for comprehensive energy system Pending CN111898816A (en)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN202010704522.5A CN111898816A (en) 2020-07-21 2020-07-21 Dynamic state estimation method for comprehensive energy system

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN202010704522.5A CN111898816A (en) 2020-07-21 2020-07-21 Dynamic state estimation method for comprehensive energy system

Publications (1)

Publication Number Publication Date
CN111898816A true CN111898816A (en) 2020-11-06

Family

ID=73190768

Family Applications (1)

Application Number Title Priority Date Filing Date
CN202010704522.5A Pending CN111898816A (en) 2020-07-21 2020-07-21 Dynamic state estimation method for comprehensive energy system

Country Status (1)

Country Link
CN (1) CN111898816A (en)

Cited By (4)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN112906317A (en) * 2021-02-09 2021-06-04 南京信息工程大学 Robust dynamic state estimation method for natural gas pipe network
CN113486532A (en) * 2021-07-22 2021-10-08 东南大学 Dynamic safety control method for electric heating comprehensive energy system
CN117195778A (en) * 2023-11-08 2023-12-08 天津市津安热电有限公司 Parameter identification correction method for hydraulic simulation model of heating pipe network
CN117291445A (en) * 2023-11-27 2023-12-26 国网安徽省电力有限公司电力科学研究院 Multi-target prediction method based on state transition under comprehensive energy system

Citations (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN108448568A (en) * 2018-03-08 2018-08-24 国网山东省电力公司潍坊供电公司 Power distribution network admixture method of estimation based on a variety of time cycle measurement data
CN110866688A (en) * 2019-11-12 2020-03-06 广州供电局有限公司 Real-time risk assessment method for multi-energy coupling system
CN111428351A (en) * 2020-03-11 2020-07-17 国网辽宁省电力有限公司大连供电公司 Electric-thermal comprehensive energy system load flow calculation method based on forward-backward substitution method

Patent Citations (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN108448568A (en) * 2018-03-08 2018-08-24 国网山东省电力公司潍坊供电公司 Power distribution network admixture method of estimation based on a variety of time cycle measurement data
CN110866688A (en) * 2019-11-12 2020-03-06 广州供电局有限公司 Real-time risk assessment method for multi-energy coupling system
CN111428351A (en) * 2020-03-11 2020-07-17 国网辽宁省电力有限公司大连供电公司 Electric-thermal comprehensive energy system load flow calculation method based on forward-backward substitution method

Non-Patent Citations (4)

* Cited by examiner, † Cited by third party
Title
LIANG CHEN 等: "Trajectory tracking method of natural gas, district heating and power systems", 《ENERGY CONVERSION AND MANAGEMENT》, vol. 259, pages 1 - 12 *
李波: "混合量测下电力系统动态状态估计研究", 《中国优秀硕士学位论文全文数据库 工程科技Ⅱ辑》, no. 07, pages 042 - 750 *
艾小猛 等: "一种考虑天然气系统动态过程的气电联合系统优化运行模型", 《电网技术》, vol. 42, no. 2, pages 409 - 412 *
钟俊杰: "综合能源系统的能量流优化与故障辨识", 《中国优秀硕士学位论文全文数据库 工程科技Ⅱ辑》, no. 06, pages 039 - 43 *

Cited By (8)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN112906317A (en) * 2021-02-09 2021-06-04 南京信息工程大学 Robust dynamic state estimation method for natural gas pipe network
CN112906317B (en) * 2021-02-09 2023-08-22 南京信息工程大学 Robust dynamic state estimation method for natural gas pipe network
CN113486532A (en) * 2021-07-22 2021-10-08 东南大学 Dynamic safety control method for electric heating comprehensive energy system
CN113486532B (en) * 2021-07-22 2023-12-19 东南大学 Dynamic safety control method for electric heating comprehensive energy system
CN117195778A (en) * 2023-11-08 2023-12-08 天津市津安热电有限公司 Parameter identification correction method for hydraulic simulation model of heating pipe network
CN117195778B (en) * 2023-11-08 2024-02-20 天津市津安热电有限公司 Parameter identification correction method for hydraulic simulation model of heating pipe network
CN117291445A (en) * 2023-11-27 2023-12-26 国网安徽省电力有限公司电力科学研究院 Multi-target prediction method based on state transition under comprehensive energy system
CN117291445B (en) * 2023-11-27 2024-02-13 国网安徽省电力有限公司电力科学研究院 Multi-target prediction method based on state transition under comprehensive energy system

Similar Documents

Publication Publication Date Title
CN111898816A (en) Dynamic state estimation method for comprehensive energy system
WO2020093296A1 (en) Interval power flow calculation method for power-heat integrated energy system
CN108921727B (en) Regional comprehensive energy system reliability assessment method considering thermal load dynamic characteristics
CN107842908B (en) Real-time heat supply load control method based on environmental parameter compensation
WO2019071763A1 (en) Method for estimating state of electro-thermal coupling system in combination with dynamic characteristics of pipeline
CN110518583B (en) Comprehensive energy system reliability assessment method considering dynamic characteristics
CN109377008B (en) Electric-thermal coupling comprehensive energy system risk assessment method
CN111928335B (en) Secondary network hydraulic balance method based on intelligent valve
CN108920866B (en) Heat supply network dynamic regulation operating parameter estimation method based on moving horizon estimation theory
CN111780933B (en) Method and system for diagnosing leakage fault of high-pressure heater based on neural network and thermodynamic modeling
CN112290536A (en) Online scheduling method of electricity-heat comprehensive energy system based on near-end strategy optimization
CN108021742B (en) A kind of steam heating pipeline steady-state operating condition estimation method considering hydrophobic model
CN110752605B (en) Optimal power flow calculation method for electric-thermal coupling comprehensive energy system
KR20130017368A (en) Heating supply determination method for district heating network
CN111046594A (en) Hot water heating network dynamic simulation method based on cross iteration principle
Taler et al. Control of the temperature in the hot liquid tank by using a digital PID controller considering the random errors of the thermometer indications
CN111523210A (en) Prediction analysis method and system for temperature rise and drop process of urban central heating system
Tian et al. Modeling and simulation for multi energy flow coupled network computing
CN113065249B (en) Method and device for predicting supply and return water temperature of heating system
CN116611706A (en) Dynamic carbon emission factor measuring and calculating method based on multi-energy main body
CN107783941B (en) Solar energy compensation type electric boiler heat supply control method based on valley electricity price change
CN111783309A (en) Dynamic simulation method of steam heating network based on internal conservation
CN110705066B (en) Dynamic simulation method for gas-electricity coupling park comprehensive energy system based on projection integration
Qin et al. Hybrid physics and data-driven method for modeling and analysis of electricity–heat integrated energy systems
CN115081193A (en) Power grid-heat grid model construction method, fusion simulation method and system

Legal Events

Date Code Title Description
PB01 Publication
PB01 Publication
SE01 Entry into force of request for substantive examination
SE01 Entry into force of request for substantive examination