CN112036003A - 一种考虑不完全量测的质调节热力系统静态状态估计方法 - Google Patents

一种考虑不完全量测的质调节热力系统静态状态估计方法 Download PDF

Info

Publication number
CN112036003A
CN112036003A CN202010642001.1A CN202010642001A CN112036003A CN 112036003 A CN112036003 A CN 112036003A CN 202010642001 A CN202010642001 A CN 202010642001A CN 112036003 A CN112036003 A CN 112036003A
Authority
CN
China
Prior art keywords
node
measurement
formula
temperature
nodes
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.)
Granted
Application number
CN202010642001.1A
Other languages
English (en)
Other versions
CN112036003B (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.)
Southeast University
Original Assignee
Southeast University
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 Southeast University filed Critical Southeast University
Priority to CN202010642001.1A priority Critical patent/CN112036003B/zh
Publication of CN112036003A publication Critical patent/CN112036003A/zh
Application granted granted Critical
Publication of CN112036003B publication Critical patent/CN112036003B/zh
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Classifications

    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06FELECTRIC DIGITAL DATA PROCESSING
    • G06F30/00Computer-aided design [CAD]
    • G06F30/20Design optimisation, verification or simulation
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06FELECTRIC DIGITAL DATA PROCESSING
    • G06F17/00Digital computing or data processing equipment or methods, specially adapted for specific functions
    • G06F17/10Complex mathematical operations
    • G06F17/11Complex mathematical operations for solving equations, e.g. nonlinear equations, general mathematical optimization problems
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06FELECTRIC DIGITAL DATA PROCESSING
    • G06F17/00Digital computing or data processing equipment or methods, specially adapted for specific functions
    • G06F17/10Complex mathematical operations
    • G06F17/16Matrix or vector computation, e.g. matrix-matrix or matrix-vector multiplication, matrix factorization
    • HELECTRICITY
    • H02GENERATION; CONVERSION OR DISTRIBUTION OF ELECTRIC POWER
    • H02JCIRCUIT ARRANGEMENTS OR SYSTEMS FOR SUPPLYING OR DISTRIBUTING ELECTRIC POWER; SYSTEMS FOR STORING ELECTRIC ENERGY
    • H02J3/00Circuit arrangements for ac mains or ac distribution networks
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06FELECTRIC DIGITAL DATA PROCESSING
    • G06F2119/00Details relating to the type or aim of the analysis or the optimisation
    • G06F2119/08Thermal analysis or thermal optimisation
    • HELECTRICITY
    • H02GENERATION; CONVERSION OR DISTRIBUTION OF ELECTRIC POWER
    • H02JCIRCUIT ARRANGEMENTS OR SYSTEMS FOR SUPPLYING OR DISTRIBUTING ELECTRIC POWER; SYSTEMS FOR STORING ELECTRIC ENERGY
    • H02J2203/00Indexing scheme relating to details of circuit arrangements for AC mains or AC distribution networks
    • H02J2203/20Simulating, e g planning, reliability check, modelling or computer assisted design [CAD]

Landscapes

  • Engineering & Computer Science (AREA)
  • Physics & Mathematics (AREA)
  • General Physics & Mathematics (AREA)
  • Mathematical Physics (AREA)
  • Theoretical Computer Science (AREA)
  • Data Mining & Analysis (AREA)
  • Pure & Applied Mathematics (AREA)
  • Computational Mathematics (AREA)
  • Mathematical Analysis (AREA)
  • Mathematical Optimization (AREA)
  • General Engineering & Computer Science (AREA)
  • Algebra (AREA)
  • Databases & Information Systems (AREA)
  • Software Systems (AREA)
  • Operations Research (AREA)
  • Evolutionary Computation (AREA)
  • Geometry (AREA)
  • Computer Hardware Design (AREA)
  • Computing Systems (AREA)
  • Power Engineering (AREA)
  • Control Of Temperature (AREA)
  • Management, Administration, Business Operations System, And Electronic Commerce (AREA)

Abstract

本发明公开了一种考虑不完全量测的质调节热力系统静态状态估计方法,所述方法包括以下步骤:步骤10)根据典型热网模型,建立热力统一方程;步骤20)对热力模型进行变换,给出不完全量测模型下的量测和状态量,建立不完全量测模型;步骤30)确定系统内的量测配置,建立质调节热力系统中的信息矩阵,利用静态最小二乘法直接计算最优估计状态向量。该方法考虑实际系统中量测配置的技术性与经济性约束,所建立的不完全量测模型适用于各种量测配置情况,基于信息矩阵直接给出状态量的最优估计,十分具有工程意义。

Description

一种考虑不完全量测的质调节热力系统静态状态估计方法
技术领域
本发明属于能源系统运行优化技术领域,具体来说,涉及一种考虑不完全量测的质调节热力系统静态状态估计方法。
背景技术
日益增长的能源消耗与环境压力引起了社会对于高效率、低排放的能源利用方式的需求。综合能源系统通过将多种异质能流耦合,为不同能源系统之间的运行优化提供了更大的灵活性空间,实现了不同能量子系统之间的协同规划和联合经济调度。以热电联产机组、热泵和电锅炉等典型的能源转换设备为核心的电热型综合能源系统主要应用于园区级系统,在世界范围内得到了广泛发展。随着经济发展与技术进步,大量的热电联产设备投入使得电力系统与热力系统之间的耦合越发紧密,针对电热综合能源系统在工程中的应用亟需深入研究。
状态分析是一种用于评估系统运行状态的统计方法,为系统的在线运行控制提供了可靠的量测信息,是电热联合系统安全分析、经济运行以及优化调度的分析基础。电力系统已针对状态估计方法展开了大量研究,但是相关研究在热力系统中相对不够成熟。现有方法一般假设热力系统处于全量测状态,以此为基础进行状态估计,该方法虽然估计精度较高,然而,考虑到经济性与技术性约束,实际工程中往往不存在全量测的热力系统。考虑系统中不完全量测的情况,当前的部分研究通过智能算法对缺失的信息进行预测进而不全量测量,该方法本质上仍是一种全量测模型,且估计结果依赖于预测值而缺乏稳定性。因此,有必要从机理建模的角度出发,研究一种通用的考虑不完全量测的热力系统状态估计模型,使其更符合实际工程中的应用条件。
发明内容
本发明所要解决的技术问题是:提出一种考虑不完全量测的质调节热力系统静态状态估计方法,该方法首先根据现有热力模型中,管道温降方程与节点温度混合方程相互耦合的特点,建立一种热力联合方程,将两种热力方程统一。在此基础上,建立考虑不完全量测的质调节热力系统量测模型,并给出了不同的量测配置的建模方法。基于量测方程推导量测矩阵进而建立热力系统信息矩阵,从而直接计算最优状态估计值。该模型考虑了实际工程中的经济性与技术性约束,可应用于各种量测配置系统,十分具有应用价值。
为解决上述技术问题,本技术方案采用一种考虑不完全量测的质调节热力系统静态状态估计方法,该方法包括以下步骤:
步骤10)根据典型热网模型,建立热力联合方程;
步骤20)对热力模型进行变换,给出不完全量测模型下的量测和状态量,建立不完全量测模型;
步骤30)确定系统内的量测配置,建立质调节热力系统中的信息矩阵,计算最优估计状态向量。
作为本发明的进一步介绍,所述的步骤10)具体包括:
步骤101)质调节热网中的水力状态为常数,因此,主要针对热力模型建模,包含管道温降方程和节点温度混合方程,可分别表示为:
Figure BDA0002571834880000021
(∑mout)Tout=(∑minTin) (2)
式中,Tpe和Tps分别表示管道的末端温度和首端温度;l为管道长度;λ为管道导热系数;Ta为环境温度;Cρ为工质比热容;mout为流入节点的管道流量,min为节点的流出流量;Tout为节点流出温度,Tin为节点流入管道的末端温度。
步骤102)根据不同的节点类型,节点流出流量包含由节点流出至管道的流量Ms,out以及节点注出流量Md,out,分别表示为
Figure BDA0002571834880000022
式中,As,out为节点和流出管道的关联矩阵,as,out,ij=1表示从节点i流出流量至管道j,as,out,ij=0表示由节点i流出的流量与管道j无关,dout为节点注出流量,仅在负荷节点处不为0。节点流入流量包含由管道流入节点的流量Ms,in和节点注入流量Md,in,分别表示为
Figure BDA0002571834880000023
式中,As,in为节点和流入管道的关联矩阵,as,in,ij=1表示从管道j流出流量至节点i,as,out,ij=0表示由节点i流入的流量与管道j无关,din为节点注入流量,仅在热源节点处不为0。
将式(1),式(3)和式(4)代入式(2),可得到式(5)所示方程,其中Ms为热力系统各节点的绝对流出流量矩阵。
Figure BDA0002571834880000031
由于各管道的起始温度等于节点的流出温度,其关系可表示为:
Tps=AsTs (6)
式中,as,ij=1表示管道i的起始温度等于节点j的流出温度,as,ij=0表示管道i的起始温度与节点j的流出温度无关。将式(6)待入式(5),可得到
Figure BDA0002571834880000032
步骤103)式(7)适用于热力系统内所有节点,将其展开,可表示为:
Figure BDA0002571834880000033
式中,Ms,sr,Ms,int和Ms,ld分别表示供水网络中源节点、联络节点和负荷节点的绝对流出流量;Ts,sr,Ts,int和Ts,ld分别表示源节点、联络节点和负荷节点的供水温度;Jbs,sr,Jbs,int和Jbs,ld分别表示源节点、联络节点和负荷节点对应的常数向量;Jsij表示由节点经管道流入节点i的流量之和;以上向量均可根据式(7)对应元素获取。根据热力系统的拓扑关系,Js11,Js12,Js13,Js23和Js33均为零矩阵,其余元素可以非零。
步骤104)类比于供水网,推导适用于回水网式(8)适用于热力系统内所有节点,将其展开,可表示为:
Figure BDA0002571834880000034
式中,Mr和Jr分别表示回水网中的节点的绝对流量流出矩阵,Jr为系统内的节点流出流量矩阵,Mr,in为节点的流入流量矩阵,Ar回水网的管道起始温度-节点流量关联温度,Tr为节点的回热温度。
类比于式(8),将式(9)展开,可表示为:
Figure BDA0002571834880000041
式中,Mr,sr,Mr,int和Mr,ld分别表示回水网中源节点、联络节点和负荷节点的绝对流出流量;Tr,sr,Tr,int和Tr,ld分别表示源节点、联络节点和负荷节点的供水温度;Jbr,sr,Jbr,int和Jbr,ld分别表示源节点、联络节点和负荷节点在回水网中对应的常数向量;Jrij表示回水网由节点经管道流入节点i的流量之和;以上向量均可根据式(9)对应元素获取。根据热力系统的拓扑关系,Jr11,Jr21,Jr31,Jr32和Jr33均为零矩阵,其余元素可以非零。
作为本发明的进一步介绍,所述步骤20)包括:
步骤201)质调节热网中的量测量包括源节点的供水温度(Ts,ld),负荷节点的回水温度(Tr,sr),联络节点的热功率(φint)以及节点的注入流量;状态量包括非源节点的供水温度(Ts,ns),非负荷节点的回水温度(Tr,nl)以及负荷、热源节点的热功率(φldsr),冗余量测包含部分非源节点的供水温度(Ts,ns)和非负荷节点的回水温度(Tr,nl)。状态量可表示为:
xh={Ts,nl,Tr,nssrld} (11)
量测量可表示为:
Zh={Ts,ld,Tr,nsint,Ts,nl,Tr,nl} (12)
步骤202)对热力模型式(9)和式(10)进行变换,建立不完全量测模型,可表示为:
0=Js21Ts,sr+(Js22-Ms,int)Ts,int+Jbs,int (13)
0=(Jr22-Mr,int)Tr,int+Jr23Tr,ld+Jbr,int (14)
Figure BDA0002571834880000042
Figure BDA0002571834880000043
CρdoutTs,ld=φld+CρdoutTr,ld (17)
CρdinTr,sr=-φsr+CρdinTs,sr (18)
Ts,nl=Ts,nl (19)
Tr,ns=Tr,ns (20)
式(13)-式(18)表示通过间接量测信息对状态量进行估计,式(19)-式(20)表示通过直接量测信息对状态量进行估计。式(13)-式(20)中,等式均由量测量构成,等式右侧均由状态量构成,其中,式(13)和式(14)左侧分别表示联络节点在供热网和回热网中的联络节点处的广义能量不平衡量,式(17)和式(18)左侧分别表示负荷节点的供热能量和热源节点的回热能量。
步骤203)以确保可观测性为前提,质调节热力系统中的不同量测配置对应的量测方程可表示为:
表1.质调节热力系统中不同量测配置下的量测模型
Figure BDA0002571834880000051
其中,ths和thr分别表示联络节点在供热网和回热网中的联络节点处的广义能量不平衡量,即式(13)和式(14)左侧量;φld,s和φsr,r分别表示负荷节点的供热能量和热源节点的回热能量;nh和nh3分别表示热力系统内的节点数目与联络节点数目。表1的类型1中,量测量维度等于状态量,即不存在冗余,对应方程为能流计算方程,类型3中所有所有状态量均可以通过直接量测进行估计,是一种全量测配置,实际情况中,任意类型的量测配置均可以表示为类型1和类型2的组合。
作为本发明的进一步介绍,所述的步骤30)包括:
步骤301)确定系统内的量测配置,根据表1选择量测方程,建立不完全量测下的质调节热力系统静态状态估计方法,可表示为:
Figure BDA0002571834880000052
式中,rh表示状态估计的残差,Wh为量测量的权重矩阵,f(xh)为量测方程。根据式(13)-式(19),量测量与状态量呈线性关系,可表示为:
f(xh)=Hhxh+bh (22)
其中,Hh为质调节热力模型的量测矩阵,可表示为:
Figure BDA0002571834880000061
常数向量bh可表示为
Figure BDA0002571834880000062
步骤302)根据量测矩阵Hh建立质调节热力模型中的信息矩阵,直接计算最优状态估计值。其中,信息矩阵Ch可表示为:
Ch=(Hh)TWhHh (25)
根据信息矩阵,质调节热力模型中的最优状态估计值可表示为:
Figure BDA0002571834880000063
其中,Xh表示状态量的最优估计值。
与现有技术相比,本发明具有以下有益效果:该方法考虑工程中的经济性与技术性约束,提出了一种考虑不完全量测的状态估计模型,使之更符合工程实际应用;同时,该模型可适用于不同类型的量测配置情况,具有普适性。在此基础上,类比于电力系统推导了热力系统中的信息矩阵与最优状态估计值,计算简单,避免了复杂的优化问题求解,而热力系统状态矩阵直观反映了系统中不同变量对于状态精度的增益情况,衡量了量测冗余对于估计精度的改善效果,为热力系统中量测的优化配置提供理论指导。
附图说明
图1为本发明实施例的具体方法;
图2为本发明实施例中考虑不完全量测的质调节热力系统静态状态估计方法流程图;
图3为本发明实施例中采用的热力系统结构图;
图4为本发明实施例中状态估计结果与真实值对比图;
图5为本发明实施例中随着量测维度增加,热力状态估计的平均相对误差与滤波效果演变图。
具体实施方式
下面结合实例和附图,对本发明实施例的技术方案做进一步的说明。
实施例1:一种考虑不完全量测的质调节热力系统静态状态估计方法,该方法包括以下步骤:
步骤10)根据典型热网模型,建立热力联合方程;
步骤20)对热力模型进行变换,给出不完全量测模型下的量测和状态量,建立不完全量测模型;
步骤30)确定系统内的量测配置,建立质调节热力系统中的信息矩阵,计算最优估计状态向量。
所述的步骤10)具体包括:
步骤101)质调节热网中的水力状态为常数,因此,主要针对热力模型建模,包含管道温降方程和节点温度混合方程,可分别表示为:
Figure BDA0002571834880000071
(∑mout)Tout=(∑minTin) (2)
式中,Tpe和Tps分别表示管道的末端温度和首端温度;l为管道长度;λ为管道导热系数;Ta为环境温度;Cρ为工质比热容;mout为流入节点的管道流量,min为节点的流出流量;Tout为节点流出温度,Tin为节点流入管道的末端温度。
步骤102)根据不同的节点类型,节点流出流量包含由节点流出至管道的流量Ms,out以及节点注出流量Md,out,分别表示为
Figure BDA0002571834880000072
式中,As,out为节点和流出管道的关联矩阵,as,out,ij=1表示从节点i流出流量至管道j,as,out,ij=0表示由节点i流出的流量与管道j无关,dout为节点注出流量,仅在负荷节点处不为0。节点流入流量包含由管道流入节点的流量Ms,in和节点注入流量Md,in,分别表示为
Figure BDA0002571834880000081
式中,As,in为节点和流入管道的关联矩阵,as,in,ij=1表示从管道j流出流量至节点i,as,out,ij=0表示由节点i流入的流量与管道j无关,din为节点注入流量,仅在热源节点处不为0。
将式(1),式(3)和式(4)代入式(2),可得到式(5)所示方程。其中Ms为热力系统各节点的绝对流出流量矩阵。
Figure BDA0002571834880000082
由于各管道的起始温度等于节点的流出温度,其关系可表示为:
Tps=AsTs (6)
式中,as,ij=1表示管道i的起始温度等于节点j的流出温度,as,ij=0表示管道i的起始温度与节点j的流出温度无关。将式(6)待入式(5),可得到
Figure BDA0002571834880000083
步骤103)式(7)适用于热力系统内所有节点,将其展开,可表示为:
Figure BDA0002571834880000084
式中,Ms,sr,Ms,int和Ms,ld分别表示供水网络中源节点、联络节点和负荷节点的绝对流出流量;Ts,sr,Ts,int和Ts,ld分别表示源节点、联络节点和负荷节点的供水温度;Jbs,sr,Jbs,int和Jbs,ld分别表示源节点、联络节点和负荷节点对应的常数向量;Jsij表示由节点经管道流入节点i的流量之和;以上向量均可根据式(7)对应元素获取。根据热力系统的拓扑关系,Js11,Js12,Js13,Js23和Js33均为零矩阵,其余元素可以非零。
步骤104)类比于供水网,推导适用于回水网式(8)适用于热力系统内所有节点,将其展开,可表示为:
Figure BDA0002571834880000091
式中,Mr和Jr分别表示回水网中的节点的绝对流量流出矩阵,Jr为系统内的节点流出流量矩阵,Mr,in为节点的流入流量矩阵,Ar回水网的管道起始温度-节点流量关联温度,Tr为节点的回热温度。
类比于式(8),将式(9)展开,可表示为:
Figure BDA0002571834880000092
式中,Mr,sr,Mr,int和Mr,ld分别表示回水网中源节点、联络节点和负荷节点的绝对流出流量;Tr,sr,Tr,int和Tr,ld分别表示源节点、联络节点和负荷节点的供水温度;Jbr,sr,Jbr,int和Jbr,ld分别表示源节点、联络节点和负荷节点在回水网中对应的常数向量;Jrij表示回水网由节点经管道流入节点i的流量之和;以上向量均可根据式(9)对应元素获取。根据热力系统的拓扑关系,Jr11,Jr21,Jr31,Jr32和Jr33均为零矩阵,其余元素可以非零。
所述的步骤20)包括:
步骤201)质调节热网中的量测量包括源节点的供水温度(Ts,ld),负荷节点的回水温度(Tr,sr),联络节点的热功率(φint)以及节点的注入流量;状态量包括非源节点的供水温度(Ts,ns),非负荷节点的回水温度(Tr,nl)以及负荷、热源节点的热功率(φldsr),冗余量测包含部分非源节点的供水温度(Ts,ns)和非负荷节点的回水温度(Tr,nl)。状态量可表示为:
xh={Ts,nl,Tr,nssrld} (11)
量测量可表示为:
Zh={Ts,ld,Tr,nsint,Ts,nl,Tr,nl} (12)
步骤202)对热力模型式(9)和式(10)进行变换,建立不完全量测模型,可表示为:
0=Js21Ts,sr+(Js22-Ms,int)Ts,int+Jbs,int (13)
0=(Jr22-Mr,int)Tr,int+Jr23Tr,ld+Jbr,int (14)
Figure BDA0002571834880000093
Figure BDA0002571834880000101
CρdoutTs,ld=φld+CρdoutTr,ld (17)
CρdinTr,sr=-φsr+CρdinTs,sr (18)
Ts,nl=Ts,nl (19)
Tr,ns=Tr,ns (20)
式(13)-式(18)表示通过间接量测信息对状态量进行估计,式(19)-式(20)表示通过直接量测信息对状态量进行估计。式(13)-式(20)中,等式均由量测量构成,等式右侧均由状态量构成,其中,式(13)和式(14)左侧分别表示联络节点在供热网和回热网中的联络节点处的广义能量不平衡量,式(17)和式(18)左侧分别表示负荷节点的供热能量和热源节点的回热能量。
步骤203)以确保可观测性为前提,质调节热力系统中的不同量测配置对应的量测方程可表示为:
表1.质调节热力系统中不同量测配置下的量测模型
Figure BDA0002571834880000102
其中,ths和thr分别表示联络节点在供热网和回热网中的联络节点处的广义能量不平衡量,即式(13)和式(14)左侧量;φld,s和φsr,r分别表示负荷节点的供热能量和热源节点的回热能量;nh和nh3分别表示热力系统内的节点数目与联络节点数目。表1的类型1中,量测量维度等于状态量,即不存在冗余,对应方程为能流计算方程,类型3中所有所有状态量均可以通过直接量测进行估计,是一种全量测配置,实际情况中,任意类型的量测配置均可以表示为类型1和类型2的组合。
所述步骤30)包括:
步骤301)确定系统内的量测配置,根据表1选择量测方程,建立不完全量测下的质调节热力系统静态状态估计模型,可表示为:
Figure BDA0002571834880000103
式中,rh表示状态估计的残差,Wh为量测量的权重矩阵,f(xh)为量测方程。根据式(13)-式(19),量测量与状态量呈线性关系,可表示为:
f(xh)=Hhxh+bh (22)
其中,Hh为质调节热力模型的量测矩阵,可表示为:
Figure BDA0002571834880000111
常数向量bh可表示为
Figure BDA0002571834880000112
步骤302)根据量测矩阵Hh建立质调节热力模型中的信息矩阵,直接计算最优状态估计值。其中,信息矩阵Ch可表示为:
Ch=(Hh)TWhHh (25)
根据信息矩阵,质调节热力模型中的最优状态估计值可表示为:
Figure BDA0002571834880000113
其中,Xh表示状态量的最优估计值。
应用实施例:以巴厘岛网状热力系统为算例进行说明。为便于对比估计精度,在节点2和节点5处添加冗余的供水温度量测,在节点13和节点22处添加冗余的回水温度量测,系统结构图如图3所示,
如图1所示,本发明实施例提供一种考虑不完全量测的质调节热力系统静态状态估计方法,包括以下步骤:
步骤10)根据典型热网模型,建立热力联合方程;
步骤20)对热力模型进行变换,给出不完全量测模型下的量测和状态量,建立不完全量测模型;
步骤30)确定系统内的量测配置,建立质调节热力系统中的信息矩阵,计算最优估计状态向量。
在上述实施例中,所述的步骤10)具体包括:
步骤101)质调节热网中的水力状态为常数,因此,主要针对热力模型建模,包含管道温降方程和节点温度混合方程,可分别表示为:
Figure BDA0002571834880000121
(∑mout)Tout=(∑minTin) (2)
式中,Tpe和Tps分别表示管道的末端温度和首端温度;l为管道长度;λ为管道导热系数;Ta为环境温度;Cρ为工质比热容;mout为流入节点的管道流量,min为节点的流出流量;Tout为节点流出温度,Tin为节点流入管道的末端温度。
步骤102)根据不同的节点类型,节点流出流量包含由节点流出至管道的流量Ms,out以及节点注出流量Md,out,分别表示为
Figure BDA0002571834880000122
式中,As,out为节点和流出管道的关联矩阵,as,out,ij=1表示从节点i流出流量至管道j,as,out,ij=0表示由节点i流出的流量与管道j无关,dout为节点注出流量,仅在负荷节点处不为0。节点流入流量包含由管道流入节点的流量Ms,in和节点注入流量Md,in,分别表示为
Figure BDA0002571834880000123
式中,As,in为节点和流入管道的关联矩阵,as,in,ij=1表示从管道j流出流量至节点i,as,out,ij=0表示由节点i流入的流量与管道j无关,din为节点注入流量,仅在热源节点处不为0。
将式(1),式(3)和式(4)代入式(2),可得到式(5)所示方程。其中Ms为热力系统各节点的绝对流出流量矩阵。
Figure BDA0002571834880000124
由于各管道的起始温度等于节点的流出温度,其关系可表示为:
Tps=AsTs (6)
式中,as,ij=1表示管道i的起始温度等于节点j的流出温度,as,ij=0表示管道i的起始温度与节点j的流出温度无关。将式(6)待入式(5),可得到
Figure BDA0002571834880000131
步骤103)式(7)适用于热力系统内所有节点,将其展开,可表示为:
Figure BDA0002571834880000132
式中,Ms,sr,Ms,int和Ms,ld分别表示供水网络中源节点、联络节点和负荷节点的绝对流出流量;Ts,sr,Ts,int和Ts,ld分别表示源节点、联络节点和负荷节点的供水温度;Jbs,sr,Jbs,int和Jbs,ld分别表示源节点、联络节点和负荷节点对应的常数向量;Jsij表示由节点经管道流入节点i的流量之和;以上向量均可根据式(7)对应元素获取。根据热力系统的拓扑关系,Js11,Js12,Js13,Js23和Js33均为零矩阵,其余元素可以非零。
步骤104)类比于供水网,推导适用于回水网式(8)适用于热力系统内所有节点,将其展开,可表示为:
Figure BDA0002571834880000133
式中,Mr和Jr分别表示回水网中的节点的绝对流量流出矩阵,Jr为系统内的节点流出流量矩阵,Mr,in为节点的流入流量矩阵,Ar回水网的管道起始温度-节点流量关联温度,Tr为节点的回热温度。
类比于式(8),将式(9)展开,可表示为:
Figure BDA0002571834880000134
式中,Mr,sr,Mr,int和Mr,ld分别表示回水网中源节点、联络节点和负荷节点的绝对流出流量;Tr,sr,Tr,int和Tr,ld分别表示源节点、联络节点和负荷节点的供水温度;Jbr,sr,Jbr,int和Jbr,ld分别表示源节点、联络节点和负荷节点在回水网中对应的常数向量;Jrij表示回水网由节点经管道流入节点i的流量之和;以上向量均可根据式(9)对应元素获取。根据热力系统的拓扑关系,Jr11,Jr21,Jr31,Jr32和Jr33均为零矩阵,其余元素可以非零。
在上述实施例中,所述的步骤20)具体包括:
步骤201)质调节热网中的量测量包括源节点的供水温度(Ts,ld),负荷节点的回水温度(Tr,sr),联络节点的热功率(φint)以及节点的注入流量;状态量包括非源节点的供水温度(Ts,ns),非负荷节点的回水温度(Tr,nl)以及负荷、热源节点的热功率(φldsr),冗余量测包含部分非源节点的供水温度(Ts,ns)和非负荷节点的回水温度(Tr,nl)。状态量可表示为:
xh={Ts,nl,Tr,nssrld} (11)
量测量可表示为:
Zh={Ts,ld,Tr,nsint,Ts,nl,Tr,nl} (12)
步骤202)对热力模型式(9)和式(10)进行变换,建立不完全量测模型,可表示为:
0=Js21Ts,sr+(Js22-Ms,int)Ts,int+Jbs,int (13)
0=(Jr22-Mr,int)Tr,int+Jr23Tr,ld+Jbr,int (14)
Figure BDA0002571834880000141
Figure BDA0002571834880000142
CρdoutTs,ld=φld+CρdoutTr,ld (17)
CρdinTr,sr=-φsr+CρdinTs,sr (18)
Ts,nl=Ts,nl (19)
Tr,ns=Tr,ns (20)
式(13)-式(18)表示通过间接量测信息对状态量进行估计,式(19)-式(20)表示通过直接量测信息对状态量进行估计。式(13)-式(20)中,等式均由量测量构成,等式右侧均由状态量构成,其中,式(13)和式(14)左侧分别表示联络节点在供热网和回热网中的联络节点处的广义能量不平衡量,式(17)和式(18)左侧分别表示负荷节点的供热能量和热源节点的回热能量。
步骤203)以确保可观测性为前提,质调节热力系统中的不同量测配置对应的量测方程可表示为:
表1.质调节热力系统中不同量测配置下的量测模型
Figure BDA0002571834880000151
其中,ths和thr分别表示联络节点在供热网和回热网中的联络节点处的广义能量不平衡量,即式(13)和式(14)左侧量;φld,s和φsr,r分别表示负荷节点的供热能量和热源节点的回热能量;nh和nh3分别表示热力系统内的节点数目与联络节点数目。表1的类型1中,量测量维度等于状态量,即不存在冗余,对应方程为能流计算方程,类型3中所有所有状态量均可以通过直接量测进行估计,是一种全量测配置,实际情况中,任意类型的量测配置均可以表示为类型1和类型2的组合。本实施例中,在节点2和节点5处添加冗余的供水温度量测,在节点13和节点22处添加冗余的回水温度量测
在上述实施例中,所述的步骤30)具体包括:
步骤301)确定系统内的量测配置,根据表1选择量测方程,建立不完全量测下的质调节热力系统静态状态估计方法,可表示为:
Figure BDA0002571834880000152
式中,rh表示状态估计的残差,Wh为量测量的权重矩阵,f(xh)为量测方程。根据式(13)-式(19),量测量与状态量呈线性关系,可表示为:
f(xh)=Hhxh+bh (22)
其中,Hh为质调节热力模型的量测矩阵,可表示为:
Figure BDA0002571834880000161
常数向量bh可表示为
Figure BDA0002571834880000162
步骤302)根据量测矩阵Hh建立质调节热力模型中的信息矩阵,直接计算最优状态估计值。其中,信息矩阵Ch可表示为:
Ch=(Hh)TWhHh (25)
根据信息矩阵,质调节热力模型中的最优状态估计值可表示为:
Figure BDA0002571834880000163
其中,Xh表示状态量的最优估计值。本实例中的状态估计结果如图4所示。逐步增加量测冗余度,热力系统的滤波效果与预估误差如图5所示。
本发明实施例的一种考虑不完全量测的质调节热力系统静态状态估计方法,考虑实际热力系统中的经济性与技术性约束,建立热力联合方程从而将系统中的量测量与状态量完全分离。在此基础上,建立不完全量测模型与量测矩阵,从而推到热力系统中的信息矩阵,并直接估计最优状态值。所建立模型适用于不同量测配置下的热力系统,计算简单,适用广泛,十分具有工程应用价值。

Claims (4)

1.一种考虑不完全量测的质调节热力系统静态状态估计方法,其特征在于,该方法包括以下步骤:
步骤10)根据典型热网模型,建立热力联合方程;
步骤20)对热力模型进行变换,给出不完全量测模型下的量测和状态量,建立不完全量测模型;
步骤30)确定系统内的量测配置,建立质调节热力系统中的信息矩阵,计算最优估计状态向量。
2.根据权利要求1所述的一种考虑不完全量测的质调节热力系统静态状态估计方法,其特征在于,所述步骤10)根据典型热网模型,建立热力联合方程,具体如下:
步骤101)质调节热网中的水力状态为常数,因此,主要针对热力模型建模,包含管道温降方程和节点温度混合方程,可分别表示为:
Figure FDA0002571834870000011
(∑mout)Tout=(∑minTin) (2)
式中,Tpe和Tps分别表示管道的末端温度和首端温度;l为管道长度;λ为管道导热系数;Ta为环境温度;Cρ为工质比热容;mout为流入节点的管道流量,min为节点的流出流量;Tout为节点流出温度,Tin为节点流入管道的末端温度;
步骤102)根据不同的节点类型,节点流出流量包含由节点流出至管道的流量Ms,out以及节点注出流量Md,out,分别表示为
Figure FDA0002571834870000012
式中,As,out为节点和流出管道的关联矩阵,as,out,ij=1表示从节点i流出流量至管道j,as,out,ij=0表示由节点i流出的流量与管道j无关,dout为节点注出流量,仅在负荷节点处不为0。节点流入流量包含由管道流入节点的流量Ms,in和节点注入流量Md,in,分别表示为
Figure FDA0002571834870000013
式中,As,in为节点和流入管道的关联矩阵,as,in,ij=1表示从管道j流出流量至节点i,as,out,ij=0表示由节点i流入的流量与管道j无关,din为节点注入流量,仅在热源节点处不为0;
将式(1),式(3)和式(4)代入式(2),可得到式(5)所示方程,其中Ms为热力系统各节点的绝对流出流量矩阵;
Figure FDA0002571834870000021
由于各管道的起始温度等于节点的流出温度,其关系可表示为:
Tps=AsTs (6)
式中,as,ij=1表示管道i的起始温度等于节点j的流出温度,as,ij=0表示管道i的起始温度与节点j的流出温度无关。将式(6)待入式(5),可得到
Figure FDA0002571834870000022
步骤103)式(7)适用于热力系统内所有节点,将其展开,可表示为:
Figure FDA0002571834870000023
式中,Ms,sr,Ms,int和Ms,ld分别表示供水网络中源节点、联络节点和负荷节点的绝对流出流量;Ts,sr,Ts,int和Ts,ld分别表示源节点、联络节点和负荷节点的供水温度;Jbs,sr,Jbs,int和Jbs,ld分别表示源节点、联络节点和负荷节点对应的常数向量;Jsij表示由节点经管道流入节点i的流量之和;以上向量均可根据式(7)对应元素获取,根据热力系统的拓扑关系,Js11,Js12,Js13,Js23和Js33均为零矩阵,其余元素可以非零;
步骤104)类比于供水网,推导适用于回水网式(8)适用于热力系统内所有节点,将其展开,可表示为:
Figure FDA0002571834870000024
式中,Mr和Jr分别表示回水网中的节点的绝对流量流出矩阵,Jr为系统内的节点流出流量矩阵,Mr,in为节点的流入流量矩阵,Ar回水网的管道起始温度-节点流量关联温度,Tr为节点的回热温度;
类比于式(8),将式(9)展开,可表示为:
Figure FDA0002571834870000031
式中,Mr,sr,Mr,int和Mr,ld分别表示回水网中源节点、联络节点和负荷节点的绝对流出流量;Tr,sr,Tr,int和Tr,ld分别表示源节点、联络节点和负荷节点的供水温度;Jbr,sr,Jbr,int和Jbr,ld分别表示源节点、联络节点和负荷节点在回水网中对应的常数向量;Jrij表示回水网由节点经管道流入节点i的流量之和;以上向量均可根据式(9)对应元素获取;根据热力系统的拓扑关系,Jr11,Jr21,Jr31,Jr32和Jr33均为零矩阵,其余元素可以非零。
3.根据权利要求2所述的一种考虑不完全量测的质调节热力系统静态状态估计方法,其特征在于,所述步骤20)确定系统内冗余量测配置,建立考虑不完全量测的质调节热力系统静态状态估计方法,并针对系统内不同的量测配置,给出建模方法,包括:
步骤201)质调节热网中的量测量包括源节点的供水温度(Ts,ld),负荷节点的回水温度(Tr,sr),联络节点的热功率(φint)以及节点的注入流量;状态量包括非源节点的供水温度(Ts,ns),非负荷节点的回水温度(Tr,nl)以及负荷、热源节点的热功率(φldsr),冗余量测包含部分非源节点的供水温度(Ts,ns)和非负荷节点的回水温度(Tr,nl)。状态量可表示为:
xh={Ts,nl,Tr,nssrld} (11)
量测量可表示为:
Zh={Ts,ld,Tr,nsint,Ts,nl,Tr,nl} (12)
步骤202)对热力模型式(9)和式(10)进行变换,建立不完全量测模型,可表示为:
0=Js21Ts,sr+(Js22-Ms,int)Ts,int+Jbs,int (13)
0=(Jr22-Mr,int)Tr,int+Jr23Tr,ld+Jbr,int (14)
Figure FDA0002571834870000032
Figure FDA0002571834870000033
CρdoutTs,ld=φld+CρdoutTr,ld (17)
CρdinTr,sr=-φsr+CρdinTs,sr (18)
Ts,nl=Ts,nl (19)
Tr,ns=Tr,ns (20)
式(13)-式(18)表示通过间接量测信息对状态量进行估计,式(19)-式(20)表示通过直接量测信息对状态量进行估计,式(13)-式(20)中,等式均由量测量构成,等式右侧均由状态量构成,其中,式(13)和式(14)左侧分别表示联络节点在供热网和回热网中的联络节点处的广义能量不平衡量,式(17)和式(18)左侧分别表示负荷节点的供热能量和热源节点的回热能量;
步骤203)以确保可观测性为前提,质调节热力系统中的不同量测配置对应的量测方程可表示为:
表1.质调节热力系统中不同量测配置下的量测模型
Figure FDA0002571834870000041
其中,ths和thr分别表示联络节点在供热网和回热网中的联络节点处的广义能量不平衡量,即式(13)和式(14)左侧量;φld,s和φsr,r分别表示负荷节点的供热能量和热源节点的回热能量;nh和nh3分别表示热力系统内的节点数目与联络节点数目,表1的类型1中,量测量维度等于状态量,即不存在冗余,对应方程为能流计算方程,类型3中所有所有状态量均可以通过直接量测进行估计,是一种全量测配置,实际情况中,任意类型的量测配置均可以表示为类型1和类型2的组合。
4.根据权利要求3所述的一种考虑不完全量测的质调节热力系统静态状态估计方法,其特征在于,所述步骤30)确定系统内的量测配置,建立质调节热力系统中的信息矩阵,计算最优估计状态向量,包括:
步骤301)确定系统内的量测配置,根据表1选择量测方程,建立不完全量测下的质调节热力系统静态状态估计方法,可表示为:
Figure FDA0002571834870000051
式中,rh表示状态估计的残差,Wh为量测量的权重矩阵,f(xh)为量测方程;根据式(13)-式(19),量测量与状态量呈线性关系,可表示为:
f(xh)=Hhxh+bh (22)
其中,Hh为质调节热力模型的量测矩阵,可表示为:
Figure FDA0002571834870000052
常数向量bh可表示为
Figure FDA0002571834870000053
步骤302)根据量测矩阵Hh建立质调节热力模型中的信息矩阵,直接计算最优状态估计值。其中,信息矩阵Ch可表示为:
Ch=(Hh)TWhHh (25)
根据信息矩阵,质调节热力模型中的最优状态估计值可表示为:
Figure FDA0002571834870000054
其中,Xh表示状态量的最优估计值。
CN202010642001.1A 2020-07-06 2020-07-06 一种考虑不完全量测的质调节热力系统静态状态估计方法 Active CN112036003B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN202010642001.1A CN112036003B (zh) 2020-07-06 2020-07-06 一种考虑不完全量测的质调节热力系统静态状态估计方法

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN202010642001.1A CN112036003B (zh) 2020-07-06 2020-07-06 一种考虑不完全量测的质调节热力系统静态状态估计方法

Publications (2)

Publication Number Publication Date
CN112036003A true CN112036003A (zh) 2020-12-04
CN112036003B CN112036003B (zh) 2024-01-12

Family

ID=73579105

Family Applications (1)

Application Number Title Priority Date Filing Date
CN202010642001.1A Active CN112036003B (zh) 2020-07-06 2020-07-06 一种考虑不完全量测的质调节热力系统静态状态估计方法

Country Status (1)

Country Link
CN (1) CN112036003B (zh)

Citations (6)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN106998079A (zh) * 2017-04-28 2017-08-01 东南大学 一种热电联合优化调度模型的建模方法
CN107817681A (zh) * 2017-10-16 2018-03-20 清华大学 一种基于双侧等效模型的热网稳态运行状态估计方法
CN108494021A (zh) * 2018-04-20 2018-09-04 东北大学 电-热-气综合能源系统的稳定评估与静态控制方法
CN111082417A (zh) * 2019-12-01 2020-04-28 国网辽宁省电力有限公司经济技术研究院 一种基于综合能源系统电气热联合网络的状态估计方法
CN111191182A (zh) * 2019-12-17 2020-05-22 东南大学 基于线性化热力模型的静态热电联合潮流混合求解方法
CN111310310A (zh) * 2020-01-20 2020-06-19 东南大学 一种用于量调节的热力系统静态潮流快速解耦计算方法

Patent Citations (6)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN106998079A (zh) * 2017-04-28 2017-08-01 东南大学 一种热电联合优化调度模型的建模方法
CN107817681A (zh) * 2017-10-16 2018-03-20 清华大学 一种基于双侧等效模型的热网稳态运行状态估计方法
CN108494021A (zh) * 2018-04-20 2018-09-04 东北大学 电-热-气综合能源系统的稳定评估与静态控制方法
CN111082417A (zh) * 2019-12-01 2020-04-28 国网辽宁省电力有限公司经济技术研究院 一种基于综合能源系统电气热联合网络的状态估计方法
CN111191182A (zh) * 2019-12-17 2020-05-22 东南大学 基于线性化热力模型的静态热电联合潮流混合求解方法
CN111310310A (zh) * 2020-01-20 2020-06-19 东南大学 一种用于量调节的热力系统静态潮流快速解耦计算方法

Non-Patent Citations (1)

* Cited by examiner, † Cited by third party
Title
HAIXIANG ZANG等: "A Robust State Estimator for Integrated Electrical and Heating Networks", IEEE ACCESS, vol. 7, pages 109990 - 110001, XP011740794, DOI: 10.1109/ACCESS.2019.2933525 *

Also Published As

Publication number Publication date
CN112036003B (zh) 2024-01-12

Similar Documents

Publication Publication Date Title
US11435265B2 (en) Method for estimating state of combined heat and power system
CN107817681B (zh) 一种基于双侧等效模型的热网稳态运行状态估计方法
CN109492325B (zh) 一种基于扩展能源集线器的多能耦合系统潮流分析方法
CN111082417A (zh) 一种基于综合能源系统电气热联合网络的状态估计方法
CN111222213A (zh) 一种热力网络动态仿真方法及装置
CN110661266B (zh) 用于热电联产系统动态最优能流计算的差分步长优选方法
CN112257279B (zh) 一种电热综合能源系统可行域构建方法
CN112330127A (zh) 一种多时间尺度电热综合能源系统静态安全控制方法
CN111898816A (zh) 一种综合能源系统动态状态估计方法
CN111898224A (zh) 基于分布式能源系统管网能量损失模型的优化控制装置
CN112036003A (zh) 一种考虑不完全量测的质调节热力系统静态状态估计方法
CN111310343B (zh) 一种用于综合能源系统调度的供热网络热路建模方法
CN112883662A (zh) 一种蒸汽供热网络动态运行水力状态估计方法及系统
CN111998423A (zh) 电蓄热循环风量-水温控制系统及其预测控制方法
CN116611706A (zh) 基于多能源主体的动态碳排放因子测算方法
CN113111555B (zh) 一种基于叠加解耦法的质调节热力系统能流快速计算方法
CN113690891B (zh) 一种基于解析法的电-热互联综合能源系统概率潮流确定方法
CN113283190A (zh) 一种基于沸腾换热模型的缸盖温度场的仿真方法及设备
KR101515374B1 (ko) 열병합발전 시스템으로 구축된 개별 열/전기 에너지 및 중앙 열/전기 에너지 공급 연계 시스템, 및 이를 운용하는 방법
CN112906220A (zh) 综合能源微网园区系统状态的估计方法
CN115906587B (zh) 一种综合能源系统供热管网动态耦合模型建立方法
CN113836739B (zh) 一种计及供热管网多损耗的电热联合系统潮流计算方法
Qin et al. Heating network quasi-dynamic model of multi-energy flow system based on forward method
CN112861448B (zh) 一种电-气耦合系统区间线性能流模型的求解方法及装置
CN113313369B (zh) 基于改进雅可比矩阵的热网潮流计算方法

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
GR01 Patent grant
GR01 Patent grant