CN111199012A

CN111199012A - 一种能源互联网热力网络潮流计算方法

Info

Publication number: CN111199012A
Application number: CN202010052603.1A
Authority: CN
Inventors: 齐文斌; 平原; 谭志海; 乔真; 袁洁; 高卓
Original assignee: Beijing Sifang Automation Co Ltd; Beijing Sifang Project Co ltd
Current assignee: Beijing Sifang Ruihe Technology Co ltd; Beijing Sifang Automation Co Ltd; Beijing Sifang Engineering Co Ltd
Priority date: 2020-01-17
Filing date: 2020-01-17
Publication date: 2020-05-26
Anticipated expiration: 2040-01-17
Also published as: CN111199012B

Abstract

本发明提出一种能源互联网潮流计算的热力网络节点模型及其计算方法，通过建立热力网络潮流节点模型，用热力节点的供水温度、回水温度、供水压力、回水压力作为网络状态量，用节点净注入流量为零、节点的注入能量和为零作为基础模型，实现热力潮流计算。与现有技术相比，减少了计算最小热力网路回路环节、不采用流量作为状态量，方便与电网的节点模型组合，采用牛顿法等常规一体化潮流算法，收敛性好，编码容易，提高了潮流计算效率。

Description

一种能源互联网热力网络潮流计算方法

技术领域

本发明涉及能源互联网技术领域，更具体地，涉及一种能源互联网的热力网络潮流计算方法。

背景技术

能源互联网(IOE)是以电能为核心，集成热、冷、燃气等能源，综合利用互联网等技术，深度融合能源系统与信息通信系统，协调多种能源的生产、传输、分配、转换、消费及交易，具有高效、清洁、低碳、安全特征的开放式能源互联网络。

能源互联网是在多能源生产、传输和消费环节相关信息采集的基础上，通过能效管理技术和需求侧管理技术，实现全局能源的综合协调优化，同时利用大数据技术，为能源生产运营、能源消费提供服务，为能源交易提供平台和技术支撑。

潮流计算是能源互联网优化控制的基础，已有的能源互联网潮流模型采用的是电力网络节点模型与热力网络回路模型，热力网络用回水供水温度、管线流量作为状态量进行潮流计算的。其指定其中一个热源为平衡节点，给定该节点处的水头压力值和供水热温度，供热功率和回热温度未知；其它非平衡节点热源的供热功率和供热温度已知，回热温度未知；热负荷节点的热功率和输出温度已知，供热、回热温度未知，各管道的流量也作为未知量。热力潮流模型采用回路模型实时计算。但上述潮流计算方法在计算过程中，需要在形成网络拓扑节点之后，搜索供热网络回路，形成多个独立的最小热力回路，在最小回路模型的基础上，做潮流计算，造成计算工作量庞大，不利于快速高效的算法实现。

发明内容

为了克服现有技术中存在的缺陷，本发明提出一种热力网络潮流计算方法，用热力节点的供水温度、回水温度、供水压力、回水压力作为热力网络状态量，用节点净注入流量为零、节点的注入能量作为基础模型，建立热力网络潮流节点潮流模型和热力潮流方程，所述方程以节点为计算基本单元，形成热力网络潮流计算模型，并进行潮流计算。

进一步地，所述潮流计算方法包括如下步骤：

S1：选取状态变量，建立状态方程，

其中状态变量：x＝[h_R1,h_S1,T_R1,T_S1,...,h_Ri,h_Si,T_Ri,T_Si,...,h_RN,h_SN,T_RN,T_SN]^T

状态方程：f(x)＝[⊿m_R1,⊿m_s1,⊿mT_R1,⊿mT_s1,....,⊿m_Ri,⊿m_si,⊿mT_Ri,⊿mTs_i,...,⊿m_RN,⊿m_sN,⊿mT_RN,⊿mTs_N]^T＝0

式中i＝1,2,...,N,N为网络节点数，h_Ri为节点i的回水压力，h_Si为节点i的供水压力

T_Ri为节点i的回水温度，T_Si为节点i的供水温度；⊿m_Ri和⊿m_Si分别为节点i的回水、供水净流量，⊿mT_Ri表示为节点i的回水回路注入能量和，⊿mT_si为节点i的供水回路注入能量和；S2：初始化状态变量h_Ri、h_Si、T_Ri、T_Si；

S3：计算热力网络各节点压力和温度残差△m_Ri、△m_Si、△mT_Ri、△mT_Si；

S4：测试步骤S3中计算得到的△m_Ri、△m_Si、△mT_Ri、△mT_Si标幺值的绝对值是否达到收敛标准ε，若是，则判为收敛，输出热力潮流状态变量值，若否，则转步骤S5；

S5：计算计算状态方程f(x)对状态变量x的雅可比矩阵J，其中

S6：求解雅克比矩阵方程，计算修正向量，并进一步对状态变量进行迭代修正，直到结果收敛，从而得到热力潮流状态变量值。

进一步地，所述状态方程包括节点流量平衡方程、回水节点混合温度平衡方程和供水节点混合温度平衡方程。

进一步地，所述节点流量平衡方程包括：

式中:Δm_Ri和Δm_Si分别为节点i的回水、供水净流量，h_Rj为节点j的回水压力，h_Ri为节点i的回水压力，h_Si为节点i的供水压力，h_Sj为节点j的供水压力；s_ij为流向标志参数：当h_Si>h_Sj则s_ij>0,反之，s_ij<0；K_ij为热力网络管线的阻力系数；Φ_Gi为机组的热功率，Φ_Li为负荷的功率；T_Si0为机组供水温度，T_oi0为热负荷出口温度；T_Ri为节点i的回水温度，T_Si为节点i的供水温度，C_p为水的比热容。

进一步地，所述回水节点温度平衡方程为：

其中ΔmT_Ri表示节点i的回水回路注入能量和，h_Rj为节点j的回水压力，h_Ri为节点i的回水压力；K_ij为热力网络管线的阻力系数；Φ_Gi为机组的热功率，Φ_Li为负荷的功率；T_Si0为机组供水温度，T_oi0为热负荷出口温度；T_Ri为节点i的回水温度，T_Rj为节点j的回水温度，T_Si为节点i的供水温度，C_p为水的比热容；λ为管道的热传导系数,L为管道的长度,T_a为环境温度。

进一步地，所述供水节点混合温度平衡方程为：

式中ΔmT_si为节点i的供水回路注入能量和；h_Si为节点i的供水压力，h_Sj为节点j的供水压力；K_ij为热力网络管线的阻力系数；Φ_Gi为机组的热功率，Φ_Li为负荷的功率；T_Si0为机组供水温度，T_oi0为热负荷出口温度；T_Ri为节点i的回水温度，T_Si为节点i的供水温度，T_Sj为节点j的供水温度，C_p为水的比热容；λ为管道的热传导系数,L为管道的长度,T_a为环境温度。

进一步地，所述步骤S2中，节点i状态量状态变量为：h_Ri、h_Si、T_Ri、T_Si，初值选取如下：

h_Si＝h_ni(1-rand(1.0))

h_Ri＝0.8*h_ni(1+rand(1.0))

T_Ri、T_Ri取管道的基准温度；

式中，i＝1,2,...,N,N为网络节点数，h_Ri为节点i的回水压力，h_Si为节点i的供水压力，T_Ri为节点i的回水温度，T_si为节点i的供水温度，h_ni为节点i的基准压力，rand(1.0)为小于一的随机数。

进一步地，步骤S3中，各节点的压力和温度残差计算方法如下：

式中，h_Rj为节点j的回水压力，h_Ri为节点i的回水压力，h_Si为节点i的供水压力，h_Sj为节点j的供水压力；s_ij为流向标志参数：当h_Si>h_Sj则s_ij>0,反之，s_ij<0；K_ij为热力网络管线的阻力系数；Φ_Gi为机组的热功率，Φ_Li为负荷的功率；T_Si0为机组供水温度，T_oi0为热负荷出口温度，为已知量测；T_Ri为节点i的回水温度，T_Si为节点i的供水温度，T_Sj为节点j的供水温度，C_p为水的比热容；λ为管道的热传导系数,L为管道的长度,T_a为环境温度。

进一步地，所述步骤S4中，雅克比矩阵J表示为：

进一步地，所述步骤S6包括：

S61：求解雅克比矩阵方程：

S62：计算修正向量：

S63：对状态变量进行一步迭代修正：

h_Si ^(k+1)＝h_Si ^(k)+△h_Si ^(k)

h_Ri ^(k+1)＝h_Ri ^(k)+△h_Ri ^(k)

T_Si ^(k+1)＝T_Si ^(k)+△T_Si ^(k)

T_Ri ^(k+1)＝T_Ri ^(k)+△T_Ri ^(k)

式中(k)表示迭代次数。

本发明有益的技术效果为：

通过建立热力网络潮流节点模型，进行潮流计算，用热力节点的供水温度、回水温度、供水压力、回水压力作为热力网络状态量，用节点净注入流量为零、节点的注入能量和为零作为基础模型，实现热力潮流计算，与现有的热力回路模型算法比较，减少了计算最小热力网路回路环节、不采用流量作为状态量，方便与电网的节点模型组合，采用牛顿法等常规一体化潮流算法，收敛性好，编码容易，提高了算法效率。

附图说明

图1示出了本发明热力网络计算模型图。

图2示出了本发明所述潮流计算方法的流程图。

具体实施方式

为使本发明实施例的目的、技术方案和优点更加清楚，下面将结合本发明实施例的附图，对本发明实施例的技术方案进行清楚、完整地描述。显然，所描述的实施例是本发明的一部分实施例，而不是全部的实施例。基于所描述的本发明的实施例，本领域普通技术人员所获得的所有其他实施例，都属于本发明保护的范围。

本发明以附图1所示的热力网络来说明节点发热力模型与计算实现方法。

如图1所示，其中：(1)、(2)、(3)为节点号；1、2、3为支路号，可见，图1中的热力网络包含3条热力管网1,2,3、2个热力负荷load1，load2、1个热力机组的热力网络source1。以下将通过计算图1所述热力网络的热力潮流说明本发明的模型及计算实现过程。

图1中热力网络模型参数如表1所示：

表1

支路序号	首端节点ID	尾端节点ID	长度(km)	直径(m)	参数K<sub>ij</sub>
						1	(1)	(3)	400	0.15	0.0179
2	(1)	(2)	400	0.15	0.0179
						3	(2)	(3)	600	0.15	0.0269

其中已知量包括：

功率量测：负荷功率Φ_L2＝Φ_L1＝0.3MW.

热力机组功率Φ_G3＝0.3MW.

温度量测：

机组供水温度T_S30＝100℃

负荷供水温度To₁₀＝T_O20＝50℃

环境温度T_a＝10℃

1.热力网络潮流节点模型

本发明通过节点的流量方程，热能流方程，用节点的压力h、温度T作为状态变量，建立热力节点潮流模型，热力潮流方程，用节点净注入流量(⊿m_Rj、⊿m_Sj)为零、节点的注入能量和(⊿mT_Ri、⊿mTs_i)为零作为状态方程，建立基础模型模型，方程以节点为计算基本单元，形成潮流计算模型，潮流方程如下：

1)节点流量平衡

Δm_Ri和Δm_Si分别为节点的回水、供水净流量。

2)回水节点混合温度

3)供水节点混合温度

上式中，in为流入支路的节点集合，out为流入支路的节点集合，对于供水回路，热力机组的T_Si0已知。T_oi0为热负荷出口温度，为已知量测。h_Rj为节点j的回水压力，h_Ri为节点i的回水压力，h_Si为节点i的供水压力，h_Sj为节点j的供水压力；s_ij为流向标志参数：当h_Si>h_Sj则s_ij>0,反之，s_ij<0；K_ij为热力网络管线的阻力系数；Φ_Gi为机组的热功率，Φ_Li为负荷的功率；T_Si0为机组供水温度，T_oi0为热负荷出口温度，为已知量测；T_Ri为节点i的回水温度，T_Si为节点i的供水温度，T_Sj为节点j的供水温度，C_p为水的比热容；λ为管道的热传导系数,L为管道的长度,T_a为环境温度。

图1所示的热力网络中有三个拓扑节点，其中(3)号节点为平衡节点，为参考节点量，潮流计算的未知节点为(1)、(2)。由此，图2的热力网络中，

状态方程为：f(x)＝[m_R1,m_S1,mT_R1,mT_S1,,m_R2,m_S2,mT_R2,mTs₂]^T＝0

状态变量为：x＝[h_R1,h_S1,T_R1,T_S1,h_R2,h_S2,T_R2,T_S2]^T

2.图2示出了本发明的热力网络潮流节点法计算，如图所示，其包含以下步骤：

步骤1：节点压力、温度初始化；

[h_R1,h_S1,T_R1,T_S1,...,h_Ri,h_Si,T_Ri,T_Si,...,h_Ri,h_Si,T_Ri,T_Si]^T

节点i状态量状态变量为：h_Ri、h_Si、T_rR、T_Si,赋初值，初值选取如下：

h_Si＝h_ni(1-rand(1.0))

h_Ri＝0.8*h_ni(1+rand(1.0))

T_Ri、T_Ri取管道的基准温度；

i＝1,2,...,N,N为网络节点数

其中h_Ri为节点i的回水压力；h_Si为节点i的供水压力；T_Ri为节点i的回水温度；T_Si为节点i的供水温度；h_ni为节点i的基准压力；rand(1.0)为小于一的随机数；供水节点，用节点的基准压力减去基准压力乘以小于1的随机数，作为初始压力；回水节点，用节点的0.8倍的基准压力加上0.8倍基准压力乘以小于1的随机数，作为初始压力。

本实施例中，取h_ni＝2.0

得到初始值x＝[h_R1,h_S1,T_R1,T_S1,h_R2,h_S2,T_R2,T_S2]^T为：

x＝[1.68,1.92,50,100,1.64,1.96,50,100]^T

步骤2：计算各节点残差△mRi、△mSi、△mTRi、△mTSi；

利用下式所述的热力网络潮流节点模型来计算节点残差

首次迭代：

f(x)＝[⊿m_R1,⊿m_S1,⊿mT_R1,⊿mT_S1,⊿m_R2,⊿m_S2,⊿mT_R2,⊿mTs₂]^T为：

f(x)＝[1.6211,-1.6211,81.8187,-160.3976,4.3266,-4.3266,216.3330,-428.3807]^T

步骤3：测试△m_Ri、△m_Si、△mT_Ri、△mT_Si标幺值的绝对值是否达到收敛标准ε，达到判为收敛，否则转步骤4；

第一次迭代：△m_Ri、△m_Si、△mT_Ri、△mT_Si标幺值的绝对值大于ε，ε取0.001,

基准值：m取100MW，mT取10000MWT

步骤4：计算f(x)对x的雅可比矩阵

各元素导数计算公式如下:

计算雅可比矩阵元素，导数元素依次如下:

[22.6221,0.0,0.0,0.0286,-16.7131,0.0,0.0,0.0；0.0,22.6221,0.0,-0.02869,0.0000,-16.7131,0.0,0.0；1131.0634,0.0,4.6985,1.4347,-835.6137,0.0,-1.6522,0.0；0.0,2262.2094,0.00,0.2365,0.00,-1671.3145,0.0,0.0；-16.7131,0.0,0.0,0.0,21.2576,0.0,0.0,0.0286；0.0,-16.7131,0.0,0.0,0.0,21.2576,0.0,-0.0286；-835.6572,0.0,0.0,0.0,1062.8830,0.0,5.7326,1.4347；0.0,-1671.2167,0.0,-1.6522,0.0,2125.6584,0.0,-1.4347；]

步骤5：解方程

计算修正向量△h_Ri、△h_Si、△T_Ri、△T_Si

进行一步迭代修正：

h_Si ^(k+1)＝h_Si ^(k)+△h_Si ^(k)

h_Ri ^(k+1)＝h_Ri ^(k)+△h_Ri ^(k)

T_Si ^(k+1)＝T_Si ^(k)+△T_Si ^(k)

T_Ri ^(k+1)＝T_Ri ^(k)+△T_Ri ^(k)

式中(k)表示迭代次数。转步骤2继续迭代。

首次迭代，得到结果如下：

上述计算迭代9次,最大delt＝0.00846252590,收敛。

最后结果,状态量x＝[h_R1,h_S1,T_R1,T_S1,h_R2,h_S2,T_R2,T_S2]^T为：

x＝[1.6553,1.9446,50,97.8845,1.6415,1.9584,48.9457,98.1448]^T

至此，完成了热力网络潮流计算全部过程，得到该热力潮流状态变量值。

可见，本发明能够快速准确地完成能源互联网热力网络潮流计算，并输出潮流状态变量值。

申请人结合说明书附图对本发明的实施例做了详细的说明与描述，但是本领域技术人员应该理解，以上实施例仅为本发明的优选实施方案，详尽的说明只是为了帮助读者更好地理解本发明精神，而并非对本发明保护范围的限制，相反，任何基于本发明的发明精神所作的任何改进或修饰都应当落在本发明的保护范围之内。

Claims

1.一种热力网络潮流计算方法，其特征在于，用热力节点的供水温度、回水温度、供水压力、回水压力作为热力网络状态量，用节点净注入流量为零、节点的注入能量为零作为基础模型，建立热力网络潮流节点潮流模型和热力潮流方程，所述方程以节点为计算基本单元，形成热力网络潮流计算模型，并进行潮流计算。

2.如权利要求2所述的方法，其特征在于，所述潮流计算方法包括如下步骤：

S1：选取热力网络状态变量，建立状态方程，

其中状态变量为：x＝[h_R1，h_S1，T_R1，T_S1，...，h_Ri，h_Si，T_Ri，T_Si，...，h_RN，h_SN，T_RN，T_SN]^T

状态方程为：f(x)＝[Δm_R1，Δm_S1，ΔmT_R1，ΔmT_S1，....，Δm_Ri，Δm_Si，ΔmT_Ri，ΔmT_Si，...，Δm_RN，Δm_SN，ΔmT_RN，ΔmT_SN]^T＝0

式中i＝1，2，...，N，N为网络节点数，h_Ri为节点i的回水压力，h_Si为节点i的供水压力，T_Ri为节点i的回水温度，T_Si为节点i的供水温度；Δm_Ri和Δm_Si分别为节点的回水、供水净流量，ΔmT_Ri表示为节点j的回水回路注入能量和，ΔmT_Ri表示为节点j的回水回路注入能量和；

S2：初始化状态变量h_Ri、h_Si、T_Ri、T_Si；

S3：计算热力网络各节点压力和温度残差Δm_Ri、Δm_Si、ΔmT_Ri、ΔmT_Si；

S4：测试步骤S3中计算得到的Δm_Ri、Δm_Si、ΔmT_Ri、ΔmT_Si标幺值的绝对值是否达到收敛标准ε，若是，则判为收敛，输出热力潮流状态变量值，若否，则转步骤S5；

S5：计算计算状态方程f例对状态变量x的雅可比矩阵J，其中

3.如权利要求2所述的方法，其特征在于，所述状态方程包括节点流量平衡方程、回水节点混合温度平衡方程和供水节点混合温度平衡方程。

4.如权利要求3所述的方法，其特征在于，所述节点流量平衡方程包括：

式中：Δm_Ri和Δm_Si分别为节点i的回水、供水净流量，h_Rj为节点j的回水压力，h_Ri为节点i的回水压力，h_Si为节点i的供水压力，h_Sj为节点j的供水压力；s_ij为流向标志参数：当h_Si＞h_Sj则s_ij＞0，反之，s_ij＜0；K_ij为热力网络管线的阻力系数；Φ_Gi为机组的热功率，Φ_Li为负荷的功率；T_Si0为机组供水温度，T_oi0为热负荷出口温度；T_Ri为节点i的回水温度，T_Si为节点i的供水温度，C_p为水的比热容。

5.如权利要求3所述的方法，其特征在于，所述回水节点温度平衡方程为：

其中ΔmT_Ri表示节点i的回水回路注入能量和，h_Rj为节点j的回水压力，h_Ri为节点i的回水压力；K_ij为热力网络管线的阻力系数；Φ_Gi为机组的热功率，Φ_Li为负荷的功率；T_Si0为机组供水温度，T_oi0为热负荷出口温度；T_Ri为节点i的回水温度，T_Rj为节点j的回水温度，T_Si为节点i的供水温度，C_p为水的比热容；λ为管道的热传导系数，L为管道的长度，T_a为环境温度。

6.如权利要求3所述的方法，其特征在于，所述供水节点混合温度平衡方程为：

式中ΔmT_si为节点i的供水回路注入能量和；h_Si为节点i的供水压力，h_Sj为节点j的供水压力；K_ij为热力网络管线的阻力系数；Φ_Gi为机组的热功率，Φ_Li为负荷的功率；T_Si0为机组供水温度，T_oi0为热负荷出口温度；T_Ri为节点i的回水温度，T_Si为节点i的供水温度，T_Sj为节点j的供水温度，C_p为水的比热容；λ为管道的热传导系数，L为管道的长度，T_a为环境温度。

7.如权利要求2所述的方法，其特征在于，所述步骤S2中，节点i状态量状态变量为：h_Ri、h_Si、T_rR、T_Si，初值选取如下：

h_Si＝h_ni(1-rand(1.0))

h_Ri＝0.8*h_ni(1+rand(1.0))

T_Ri、T_Ri取管道的基准温度；

式中，i＝1，2，...，N，N为网络节点数，h_Ri为节点i的回水压力，h_Si为节点i的供水压力，T_Ri为节点i的回水温度，T_si为节点i的供水温度，h_ni为节点i的基准压力，rand(1.0)为小于一的随机数。

8.如权利要求7所述的方法，其特征在于，步骤S3中，各节点的压力和温度残差计算方法如下：

式中，h_Rj为节点j的回水压力，h_Ri为节点i的回水压力，h_Si为节点i的供水压力，h_Sj为节点j的供水压力；s_ij为流向标志参数：当h_Si＞h_Sj则s_ij＞0，反之，s_ij＜0；K_ij为热力网络管线的阻力系数；Φ_Gi为机组的热功率，Φ_Li为负荷的功率；T_Si0为机组供水温度，T_oi0为热负荷出口温度，为已知量测；T_Ri为节点i的回水温度，T_Si为节点i的供水温度，T_Sj为节点j的供水温度，C_p为水的比热容；λ为管道的热传导系数，L为管道的长度，T_a为环境温度。

9.如权利要求8所述的方法，其特征在于，所述步骤S5中，雅克比矩阵J表示为：

10.如权利要求9所述的方法，其特征在于，所述步骤S6包括：

S61：求解雅克比矩阵方程：

S62：计算修正向量：

S63：对状态变量进行一步迭代修正：

h_Si ^(k+1)＝h_Si ^(k)+△h_Si ^(k)

h_Ri ^(k+1)＝h_Ri ^(k)+Δh_Ri ^(k)

T_Si ^(k+1)＝T_Si ^(k)+ΔT_Si ^(k)

T_Ri ^(k+1)＝T_Ri ^(k)+ΔT_Ri ^(k)

式中(k)表示迭代次数。