CN114114958B - 一种考虑源荷温度控制特性的供热系统状态离散仿真方法 - Google Patents

Abstract

一种考虑源荷温度控制特性的供热系统状态离散仿真方法：针对待仿真的供热系统，输入拓扑连接关系、元件参数、控制参数、仿真计算参数，进行仿真初始化；建立供热系统状态空间模型，包括热力环节微分方程组t控制环节离散代数方程组；针对热力环节微分方程组，计算热力环节状态变量在当前仿真时刻的1阶、2阶、3阶导数；计算积分步长；检测是否有控制器离散动作事件发生；使用量子化积分公式将热力环节状态变量积分至仿真时刻，置控制环节离散代数变量的历史量，计算仿真时刻的控制环节离散代数变量；当仿真结束，输出仿真结果。本发明为考虑源荷温度控制特性的供热系统的高效率仿真计算提供了必要工具。

Description

一种考虑源荷温度控制特性的供热系统状态离散仿真方法

技术领域

本发明涉及一种供热系统状态离散仿真方法。特别是涉及一种考虑源荷温度控制特性的供热系统状态离散仿真方法。

背景技术

综合能源系统通过能量的梯级利用提高一次能源利用率，在节能减排方面表现出较大的优势。随着源荷侧控制器、智能代理等装置的快速发展，具有可控热负荷和灵活响应能力的供热系统成为未来智慧综合能源系统的重要组成部分。为了开展能源系统规划设计与运行控制，有必要研究发展高效的建模仿真分析方法。供热系统源荷侧的用户行为和控制器被建模为离散变量，并和状态变量频繁交互，极大地增加了能源系统仿真问题的复杂性，需要一种建模仿真方法高效处理供热系统中的连续-离散混合特性。

传统的仿真方法包括龙格库塔法等时域积分方法，但是难以仿真包含大量离散控制环节的系统。量子化状态系统(quantized state system，QSS)方法是用于求解以离散变量主导的系统的积分算法，QSS方法通过状态变量的量子化来代替经典数值积分的时间离散化，系统的状态变量以“量子”为单位变化，依次计算每次状态变量变化所需要的时间，从而推进积分。在包含大量紧密耦合的连续变量和离散变量的供热系统中，QSS方法的仿真效率受限。

发明内容

本发明所要解决的技术问题是，提供一种能够准确、高效地刻画连续-离散变量的耦合特性的考虑源荷温度控制特性的供热系统状态离散仿真方法。

本发明所采用的技术方案是：一种考虑源荷温度控制特性的供热系统状态离散仿真方法，包括如下步骤：

1)针对待仿真的供热系统，输入拓扑连接关系、元件参数、控制参数、仿真计算参数，进行仿真初始化；

2)建立供热系统状态空间模型，其中包括管道、热源、用户建筑供热负荷在系统拓扑约束下构成的热力环节微分方程组，以及热源温度控制器和用户室内温度控制器构成的控制环节离散代数方程组；

3)针对供热系统状态空间模型中的热力环节微分方程组，计算热力环节状态变量xd在当前仿真时刻t的1阶、2阶、3阶导数

4)根据步骤3)中得到的热力环节状态变量x_d在当前仿真时刻t的1阶、2阶、3阶导数，计算积分步长h；

5)检测是否有控制器离散动作事件发生在区间(t，t+h]，若有则进行事件定位并更新积分步长h，否则直接进入步骤6)；

6)使用量子化积分公式将热力环节状态变量x_d积分至仿真时刻t+h，置控制环节离散代数变量的历史量将x_d(t+h)和/>代入控制环节离散代数方程组，计算仿真时刻t+h的控制环节离散代数变量u_d(t+h)，令当前仿真时刻t＝t+h；

7)判断仿真时刻t是否达到仿真终止时间T；若t≥T，则仿真结束，输出仿真结果；否则返回步骤3)。

本发明的一种考虑源荷温度控制特性的供热系统状态离散仿真方法，适用于考虑源荷温度控制特性的供热系统的动态仿真，能够考虑离散变量和连续变量的耦合特性，结合状态变量的量子化和时域积分的时间离散，为连续-离散混合模型提供了一种高效的仿真方法。即本发明将QSS方法的状态变量量子化和时域离散相结合，准确、高效地刻画了连续-离散变量的耦合特性，适用于求解包含大量离散控制器的供热系统。因此，本发明的方法为考虑源荷温度控制特性的供热系统的高效率仿真计算提供了必要工具。

附图说明

图1是本发明一种考虑源荷温度控制特性的供热系统状态离散仿真方法的流程图；

图2是热力系统实例拓扑图；

图3本发明实例中1号热源出口温度和控制器输出状态随时间变化图；

图4是本发明实例中热源热功率时间变化图；

图5是本发明实例中用户室内温度平均温度随时间变化图；

图6是本发明实例中室内温度和控制器输出状态随时间变化图；

图7是本发明实例中连续状态变量平均相对误差随时间变化图。

具体实施方式

下面结合实施例和附图对本发明的一种考虑源荷温度控制特性的供热系统状态离散仿真方法做出详细说明。

如图1所示，本发明的一种考虑源荷温度控制特性的供热系统状态离散仿真方法，包括如下步骤：

1)针对待仿真的供热系统，输入拓扑连接关系、元件参数、控制参数、仿真计算参数，进行仿真初始化；

所述的元件参数包括热网各管道截面积、长度、热阻、微元段长度、热工质流量，各热源的加热单元功率、数量、热工质质量和流量，各用户散热器的散热系数、散热面积、热工质质量和流量，各用户建筑的散热系数和体积，热工质密度、热工质比热容、空气密度、空气比热容、环境温度；所述的控制参数包括热源出口温度控制区间、用户室内温度控制区间；所述的仿真计算参数包括仿真终止时间、量子化积分阈值、各仿真变量及其历史量初值；置当前仿真时刻t＝0。

2)建立供热系统状态空间模型，其中包括管道、热源、用户建筑供热负荷在系统拓扑约束下构成的热力环节微分方程组，以及热源温度控制器和用户室内温度控制器构成的控制环节离散代数方程组；其中，

(1)所述的热力环节微分方程组包括：

(1.1)系统中管道热力方程

根据管道热力偏微分方程：

式中，为管道的温度变量，ρ_p为热工质密度，c_p为热工质比热容，/>为管道截面积，/>为管道热阻，T^a为环境温度，/>为热工质流量，α为管道编号；

将管道支路划分为长度为Δx微元段，根据单个微元段长度Δx和管道长度l_α得到微元段数量 符号表示向上取整；管道支路α的第n个微元段的热力微分方程形式为：

边界条件为：

式中，为管道支路α的第n个微元段的温度；/>为管道支路α的第n-1个微元段的温度，当n取1时，为/>为管道支路α的入口温度，/>为管道支路α的出口温度；

(1.2)热源热力微分方程：

式中，为热源中热工质的质量，/>和/>为热源的入口和出口温度，/>为热源热工质流量，/>为一个加热单元的功率，/>为接入热源的加热单元数量，β为热源编号；

(1.3)用户建筑供热负荷热力微分方程：

式中，为散热器中热工质的质量，/>和/>为散热器的入口和出口温度，/>为散热器散热量，/>为散热器的开关变量，/>为散热器的散热系数，/>为散热器的散热面积，/>为散热器热工质流量，/>为建筑的室内温度，ρ_a为空气的密度，c_a为空气的比热容，和/>为建筑的散热系数和体积，y为散热器编号，ε为建筑的编号；

(1.4)热力拓扑约束方程：

每条支路两端各连接着一个节点，根据输入的拓扑连接关系，得到节点i的热力拓扑约束方程为：

式中，为节点i的温度，g_k为支路k的流量，/>为支路k的出口温度，R_i为注入节点i的支路编号集合，/>为支路j的入口温度，S_i为流出节点i的支路编号集合；

(2)所述的控制环节离散代数方程组包括：

(2.1)源、荷温度控制器离散代数方程：

式中，为接入热源的加热单元数量，/>为散热器的开关变量，/>为热源的出口温度，/>为用户室内温度，/>和/>为热源出口温度的上限和下限值，/>和/>为用户室内温度的上限和下限值；

(3)所述的供热系统状态空间模型，是通过联立式(2)-(10)得到如下形式：

式中，x_d为热力环节状态变量，包含管道温度变量、热源出口温度变量、用户散热器出口温度变量、用户室内温度变量；f_d为热力环节微分方程组；u_d为控制环节离散代数变量，包含各热源控制器输出的加热单元数量、各用户温度控制器输出的散热器开关状态，为控制环节离散代数变量的历史量，z_d为控制环节离散代数方程组。

3)针对供热系统状态空间模型中的热力环节微分方程组，计算热力环节状态变量x_d在当前仿真时刻t的1阶、2阶、3阶导数

所述的计算热力环节状态变量x_d在当前仿真时刻t的1阶、2阶、3阶导数 具体计算公式为：

式中，x_d(t)、u_d(t)分别为当前仿真时刻t的热力环节状态变量和控制环节离散代数变量，分别为热力环节微分方程组的1阶、2阶导数表达式。

4)根据步骤3)中得到的热力环节状态变量x_d在当前仿真时刻t的1阶、2阶、3阶导数，计算积分步长h；所述的计算积分步长h，具体计算过程如下：

量子化积分公式为：

量子化状态变量的积分公式为：

式中，q_d(t)、和/>分别为量子化状态变量在当前仿真时刻t的值，以及量子化状态变量的1阶、2阶导数。

量子化状态变量的更新公式为：

量子化约束方程为：

h＝min{h||x_d，u(t+h)-q_d，u(t+h)|＝ΔQ，u＝1，2，...N_q} (16)

式中，h为积分步长，x_d，u(t+h)、q_d，u(t+h)分别为x_d(t+h)和q_d(t+h)中的第u个元素，ΔQ为输入的量子化积分阈值，N_q为热力环节状态变量的向量维度。

将式(13)、(14)、(15)代入式(16)，得到如下步长公式：

式中，ΔQ为输入的量子化积分阈值，为对步骤3)所得向量/>中的每个元素求绝对值后，绝对值最大的元素。

5)检测是否有控制器离散动作事件发生在区间(t，t+h]，若有则进行事件定位并更新积分步长h，否则直接进入步骤6)；具体包括：

(5.1)根据当前仿真时刻t热力环节状态变量x_d(t)及1阶、2阶、3阶导数 进行上限和下限条件判断：

上限条件：且

下限条件：且

式中，x_d，p(t)为x_d(t)中的第i个元素，(x_d，p(t))^(m)表示x_d，p(t)的m阶导数，Ω_c为x_d中热源出口温度和用户室内温度变量序号集合，分别为热源出口温度或用户室内温度控制区间的上限值和下限值，h为积分步长；

(5.2)若上限和下限条件均不满足，则没有控制器离散动作事件发生在区间(t，t+h]，直接进入步骤6)；

(5.3)若上限或下限条件满足，则求解对应的上限或下限事件的定位方程，计算控制器离散动作事件的准确发生时刻并根据积分步长调整方程，更新h：

上限事件定位方程：

下限事件定位方程：

积分步长调整方程：

式中，t为当前仿真时刻，Ω_s为x_d(t)中触发上限或者下限条件的状态变量序号集合。

6)使用量子化积分公式将热力环节状态变量x_d积分至仿真时刻t+h，置控制环节离散代数变量的历史量将x_d(t+h)和/>代入控制环节离散代数方程组，计算仿真时刻t+h的控制环节离散代数变量u_d(t+h)，令当前仿真时刻t＝t+h；

所述的量子化积分公式：

式中，x_d(t+h)为仿真时刻t+h的热力环节状态变量，x_d(t)为仿真时刻t的热力环节状态变量，h为积分步长，分别为热力环节状态变量x_d在仿真时刻t的1阶、2阶、3阶导数。

7)判断仿真时刻t是否达到仿真终止时间T；若t≥T，则仿真结束，输出仿真结果；否则返回步骤3)。

下面给出具体实例：

本实例以MATLAB编程语言环境为基础，实现了考虑源荷温度控制特性的供热系统状态离散仿真方法，通过供热系统算例对本发明的方法进行验证与分析。仿真测试的硬件平台为Intel Core(TM)i7-8700 CPU@3.20GHz，8GB RAM的4核PC机；软件环境为64位Windows 10操作系统。

如图2所示，供热系统分为供水网和回水网两个部分，供水网和回水网完全对称，对应的管道支路参数相同，热源支路、热用户支路连接供水网和回水网，形成闭合的热工质流动回路，供水网和回水网由管道支路和节点组成，供水网节点编号为1～32，回水网中节点编号为33～64，编号顺序与供水网一致；供水网管道支路编号为1～32，回水网管道支路编号为33～64，编号顺序与供水网一致；回水网的节点和支路编号顺序与供水网一致，不在图中列出，可以比较容易的推断出来。热源和热用户支路编号为65～88。系统拓扑如图1所示。65号支路为1号热源所在位置；87号支路为2号热源所在位置；88号支路为3号热源所在位置。1～21号热负荷分别位于66～86号支路上，66～86号支路中每一条支路都为10栋建筑供热。

设置仿真初始时间t＝0，仿真终止时间T＝3600s，量子化积分阈值ΔQ＝1×10^-6，状态变量初值x_d(0)，离散变量初值u_d(0)，离散变量历史值其余参数如表1和表2所示，设置仿真场景如下：

用户室内温度控制区间为24～25摄氏度，热源出口温度控制区间为80～82摄氏度。采用三种方法对比仿真精度，基准为向后微分公式(backward differentiationformula，BDF)积分方法的高精度仿真结果，对比方法为3阶量子化状态系统(third-orderquantized state system，QSS3)方法，待验证的方法为本专利提出的考虑源荷温度控制特性的状态离散仿真方法。

图3是1号热源出口温度随时间变化图，当温度到达设定上限时，1号热源关闭一个加热单元，温度下降；图4为1号热源热功率，加热单元关闭后，1号热源的加热功率由0.8MW下降至0.6MW；5为室内温度及其平均值，所有建筑的室内温度设定为24～25摄氏度，控制散热器的开关保持室内温度在区间内；图6为一个建筑的室内温度和控制器输出状态，当温度到达设定上限时，关闭散热器，当温度到达设定下限时，开启散热器，使得温度精确保持在给定区间；图7为与基准对比的所有变量的平均相对误差随时间变化图。仿真效率如表3所示，由图7和表3可知，本专利提出的考虑源荷温度控制特性的供热系统状态离散仿真方法相比于QSS3算法有更高的效率，能够高效仿真源荷侧包含离散控制器的供热系统，在离散-连续耦合系统建模仿真中具有优势。

表1支路参数

注：与供水管道对应的回水管道(编号为33-64)具有相同参数，不再列出。

表2其他参数

表3仿真效率

方法	仿真耗时(s)
		QSS3	26.8
本发明方法	5.2

Claims (7)

1.一种考虑源荷温度控制特性的供热系统状态离散仿真方法，其特征在于，包括如下步骤：

1)针对待仿真的供热系统，输入拓扑连接关系、元件参数、控制参数、仿真计算参数，进行仿真初始化；

2)建立供热系统状态空间模型，其中包括管道、热源、用户建筑供热负荷在系统拓扑约束下构成的热力环节微分方程组，以及热源温度控制器和用户室内温度控制器构成的控制环节离散代数方程组；

3)针对供热系统状态空间模型中的热力环节微分方程组，计算热力环节状态变量x_d在当前仿真时刻t的1阶、2阶、3阶导数

4)根据步骤3)中得到的热力环节状态变量x_d在当前仿真时刻t的1阶、2阶、3阶导数，计算积分步长h；

5)检测是否有控制器离散动作事件发生在区间(t,t+h]，若有则进行事件定位并更新积分步长h，否则直接进入步骤6)；

6)使用量子化积分公式将热力环节状态变量x_d积分至仿真时刻t+h，置控制环节离散代数变量的历史量将x_d(t+h)和/>代入控制环节离散代数方程组，计算仿真时刻t+h的控制环节离散代数变量u_d(t+h)，令当前仿真时刻t＝t+h；

7)判断仿真时刻t是否达到仿真终止时间T；若t≥T，则仿真结束，输出仿真结果；否则返回步骤3)。

2.根据权利要求1所述的一种考虑源荷温度控制特性的供热系统状态离散仿真方法，其特征在于，步骤1)所述的元件参数包括热网各管道截面积、长度、热阻、微元段长度、热工质流量，各热源的加热单元功率、数量、热工质质量和流量，各用户散热器的散热系数、散热面积、热工质质量和流量，各用户建筑的散热系数和体积，热工质密度、热工质比热容、空气密度、空气比热容、环境温度；所述的控制参数包括热源出口温度控制区间、用户室内温度控制区间；所述的仿真计算参数包括仿真终止时间、量子化积分阈值、各仿真变量及其历史量初值；置当前仿真时刻t＝0。

3.根据权利要求1所述的一种考虑源荷温度控制特性的供热系统状态离散仿真方法，其特征在于，步骤2)中

(1)所述的热力环节微分方程组包括：

(1.1)系统中管道热力方程

将管道支路划分为长度为Δx微元段，根据单个微元段长度Δx和管道长度l_α得到微元段数量 符号表示向上取整；管道支路α的第n个微元段的热力微分方程形式为：

边界条件为：

式中，ρ_p为热工质密度，c_p为热工质比热容，为管道截面积，/>为管道热阻，T^a为环境温度，/>为热工质流量，α为管道编号；/>为管道支路α的第n个微元段的温度；/>为管道支路α的第n-1个微元段的温度，当n取1时，为/>为管道支路α的入口温度，/>为管道支路α的出口温度；

(1.2)热源热力微分方程：

式中，为热源中热工质的质量，/>和/>为热源的入口和出口温度，/>为热源热工质流量，/>为一个加热单元的功率，/>为接入热源的加热单元数量，β为热源编号；

(1.3)用户建筑供热负荷热力微分方程：

式中，为散热器中热工质的质量，/>和/>为散热器的入口和出口温度，/>为散热器散热量，/>为散热器的开关变量，/>为散热器的散热系数，/>为散热器的散热面积，/>为散热器热工质流量，/>为建筑的室内温度，ρ_a为空气的密度，c_a为空气的比热容，/>和/>为建筑的散热系数和体积，γ为散热器编号，ε为建筑的编号；

(1.4)热力拓扑约束方程：

每条支路两端各连接着一个节点，根据输入的拓扑连接关系，得到节点i的热力拓扑约束方程为：

式中，为节点i的温度，g_k为支路k的流量，/>为支路k的出口温度，R_i为注入节点i的支路编号集合，/>为支路j的入口温度，S_i为流出节点i的支路编号集合；

(2)所述的控制环节离散代数方程组包括：

(2.1)源、荷温度控制器离散代数方程：

式中，为接入热源的加热单元数量，/>为散热器的开关变量，/>为热源的出口温度，/>为用户室内温度，/>和/>为热源出口温度的上限和下限值，/>和/>为用户室内温度的上限和下限值；

(3)所述的供热系统状态空间模型，是通过联立式(1)-(9)得到如下形式：

式中，x_d为热力环节状态变量，包含管道温度变量、热源出口温度变量、用户散热器出口温度变量、用户室内温度变量；f_d为热力环节微分方程组；u_d为控制环节离散代数变量，包含各热源控制器输出的加热单元数量、各用户温度控制器输出的散热器开关状态，为控制环节离散代数变量的历史量，z_d为控制环节离散代数方程组。

4.根据权利要求1所述的一种考虑源荷温度控制特性的供热系统状态离散仿真方法，其特征在于，步骤3)中所述的计算热力环节状态变量x_d在当前仿真时刻t的1阶、2阶、3阶导数具体计算公式为：

式中，x_d(t)、u_d(t)分别为当前仿真时刻t的热力环节状态变量和控制环节离散代数变量，分别为热力环节微分方程组的1阶、2阶导数表达式。

5.根据权利要求1所述的一种考虑源荷温度控制特性的供热系统状态离散仿真方法，其特征在于，步骤4)所述的计算积分步长h，具体公式如下：

式中，ΔQ为输入的量子化积分阈值，为对步骤3)所得向量/>中的每个元素求绝对值后，绝对值最大的元素。

6.根据权利要求1所述的一种考虑源荷温度控制特性的供热系统状态离散仿真方法，其特征在于，步骤5)具体步骤包括：

(5.1)根据当前仿真时刻t热力环节状态变量x_d(t)及1阶、2阶、3阶导数 进行上限和下限条件判断：

上限条件：且/>

下限条件：且/>

式中，x_d,p(t)为x_d(t)中的第p个元素，(x_d,p(t))^(m)表示x_d,p(t)的m阶导数，Ω_c为x_d中热源出口温度和用户室内温度变量序号集合，分别为热源出口温度或用户室内温度控制区间的上限值和下限值，h为积分步长；

(5.2)若上限和下限条件均不满足，则没有控制器离散动作事件发生在区间(t,t+h]，直接进入步骤6)；

(5.3)若上限或下限条件满足，则求解对应的上限或下限事件的定位方程，计算控制器离散动作事件的准确发生时刻并根据积分步长调整方程，更新h：

上限事件定位方程：

下限事件定位方程：

积分步长调整方程：

式中，t为当前仿真时刻，Ω_s为x_d(t)中触发上限或者下限条件的状态变量序号集合。

7.根据权利要求1所述的一种考虑源荷温度控制特性的供热系统状态离散仿真方法，其特征在于，步骤6)中所述的量子化积分公式：

式中，x_d(t+h)为仿真时刻t+h的热力环节状态变量，x_d(t)为仿真时刻t的热力环节状态变量，h为积分步长，分别为热力环节状态变量x_d在仿真时刻t的1阶、2阶、3阶导数。