CN101344761A - 热力站节能改造控制对象建模方法 - Google Patents

Abstract

热力站节能改造控制对象建模方法，它涉及一种控制对象建模方法。本发明解决了现有的建模方法存在模型精度低、无法辨识热过程的滞后时间、需要完整的阶跃响应、不适应节能改造应用等问题。本发明的主要内容为：取得单位脉冲信号、单位阶跃信号和单位斜坡信号的响应曲线；由典型传递函数取得辨识方程；利用最小二乘法辨识二阶加滞后对象模型的未知参数；经离散化后取得供热对象的差分方程。本发明所述建模方法可代替常规的飞升曲线法，可获高精度模型。模型精度相对常规的飞升曲线法建模精度由7％提高到3％以内，而且可直接辨识对象滞后时间，进而代替常规方法中由实验测取。

Description

热力站节能改造控制对象建模方法

技术领域

本发明涉及供热节能改造控制对象建模方法。

背景技术

现有建模法不适应建筑供热节能改造应用：我国三北地区建筑供热能耗占全社会能耗27％，“十一五”期间供热节能是重点，目前供热系统运行和管理多数停留在人工和半自动状态。从节能目的出发，完善城市供热系统，进行供热节能改造是异常重要的。在实施节能改造节能减排领域中，设计或改造供热控制系统，首先需要确定控制系统动态特性的数学模型。目前动态特性建模法主要有阶跃响应法、脉冲响应法、频率响应法、相关分析法、谱分析法等，但这些建模法往往对数据的长度比较敏感，需要大量的数据，另外建模的精度也比较低，因此不能满足节能改造应用的要求。

现有建模法中与本发明相似的飞升曲线法建模精度低：当阶跃响应曲线比较规则的时候，常用的飞升曲线法有近似法、半对数法、切线法和两点法。当阶跃响应曲线不规则时，上述的方法就会失效，此时可以使用面积法辨识模型参数。但面积法对阶跃响应曲线的完整性要求比较高，必须获得完全进入稳态的过程输出曲线。当数据不满足完整性要求时，面积法会失效。上述方法普遍存在一个问题，在有噪声的情况下，精度会有所下降，甚至失效。

现有建模法不能辨识过程的滞后时间：供热过程具有大惯性、大时滞的特点，因此过程的滞后时间是供热模型辨识的重要参数之一。在许多工业生产过程中，大部分纯滞后时间是由物料传输所致，对此可以通过近似估计得到，对离散时间模型，纯滞后时间τ₀一般取采样时间T的整数倍，如τ₀＝k₀T，k₀＝1，2，3，…。但在实际供热过程中，有些过程变量无法测量，这导致滞后时间也无法进行测量。

现有建模法需要完整的阶跃响应，但有时工艺不允许：阶跃响应曲线虽是一种测定过程动态特性的简单易行方法，但当过程长时间处于较大干扰信号作用时，被控量输出的变化范围可能会超过实际工艺所允许的范围，这使得建模阶跃响应曲线是一条未进入稳态的不完整曲线，因此这种情况不能使用常规飞升曲线法建模。

发明内容

基于现行建模方法存在的模型精度低、无法辨识热过程的滞后时间、需要完整的阶跃响应、不适应节能改造应用等问题，提出一种热力站节能改造控制对象建模方法。

本发明解决上述技术问题采用的技术方案是：

本发明所述建模方法是按照以下步骤实现的：

步骤一、取得单位脉冲信号、单位阶跃信号和单位斜坡信号的响应模型：

供热过程采用二阶加滞后来描述，二阶加滞后对象的模型结构为 $G\left(s\right)=\frac{K}{(T_{1}s+1)(T_{2}s+1)}{e}^{-\mathrm{\τs}},$ 其被辨识的参数为K，T₁，T₂和τ；其中s为拉普拉斯算子，单位为秒^-1；K为比例系数；T₁，T₂均为时间常数；τ为滞后时间常数；

以下分别计算二阶系统的单位脉冲信号、单位阶跃信号、单位斜坡信号的响应模型，并分别令其为A(t)，B(t)，C(t)；当单位脉冲信号的响应模型A(t)、单位阶跃信号的响应模型B(t)、单位斜坡信号的响应模型C(t)三者之一由实验获取，则另两个响应模型可通过式(1)求得：

$\left\{\begin{array}{c}B\left(t\right)={\∫}_{t_0}^{t}A\left(\mathrm{\ξ}\right)\mathrm{d\ξ}\\ C\left(t\right)={\∫}_{t_0}^{t}B\left(\mathrm{\ξ}\right)\mathrm{d\ξ}\end{array}\right.---\left(1\right)$

热力站供热对象的模型为 $G\left(s\right)=\frac{K}{(T_{1}s+1)(T_{2}s+1)}{e}^{-\mathrm{\τs}}$ 时，其中令T₂＝βT₁，β为比例系数，则单位脉冲信号的响应模型为：

$A\left(t\right)=\frac{K}{T_1(1-\mathrm{\β})}[\exp (-\frac{t-\mathrm{\τ}}{T_1})-\exp (-\frac{t-\mathrm{\τ}}{\mathrm{\β}{T}_1})]+\mathrm{\ω}\left(t\right)---\left(2\right)$

上式中ω(t)为噪声项，t为时间，以下各式同；

单位阶跃信号的响应模型为：

$B\left(t\right)={\∫}_{t_0}^{t}A\left(\mathrm{\ξ}\right)\mathrm{d\ξ}=K+\frac{K}{1-\mathrm{\β}}[\mathrm{\βexp}(-\frac{t-\mathrm{\τ}}{\mathrm{\β}{T}_1})-\exp (-\frac{t-\mathrm{\τ}}{T_1})]+{\∫}_{t_0}^t\mathrm{\ω}\left(\mathrm{\ξ}\right)\mathrm{d\ξ}---\left(3\right)$

单位斜坡信号的响应模型为：

$C\left(t\right)={\∫}_{t_0}^{t}B\left(\mathrm{\ξ}\right)\mathrm{d\ξ}$

$=K(t-\mathrm{\τ})-(\mathrm{\β}+1){\mathrm{KT}}_1+\frac{{\mathrm{KT}}_1}{1-\mathrm{\β}}[\exp (-\frac{t-\mathrm{\τ}}{T_1})-{\mathrm{\β}}^2\exp (-\frac{t-\mathrm{\τ}}{\mathrm{\β}{T}_1})]+\∫\∫\mathrm{\ω}\left(\mathrm{\ξ}\right)\mathrm{d\ξd\σ}---\left(4\right)$

步骤二、由典型传递函数取得辨识方程

由步骤一中所述的式(2)、式(3)和式(4)

可求得供热对象的二阶加滞后最小二乘辨识模型方程：

$C\left(t\right)=\left[\begin{array}{cccc}t& -1& -A\left(t\right)& -B\left(t\right)\end{array}\right]\left[\begin{array}{c}K\\ \mathrm{K\τ}\\ T_1{T}_2\\ T_1+T_2\end{array}\right]+\mathrm{\Ω}\left(t\right)---\left(5\right)$

其中Ω(t)为噪声项，是对ω(t)的滤波；

取不同时刻t₁，t₂，…，t_m，可得线性差分方程组

C＝Uθ+Ω                            (6)

其中， $C=\left[\begin{array}{c}C\left(1\right)\\ C\left(2\right)\\ \·\\ \·\\ \·\\ C\left(m\right)\end{array}\right],$ $U=\left[\begin{array}{cccc}1& -1& -A\left(1\right)& -B\left(1\right)\\ 2& -1& -A\left(2\right)& -B\left(2\right)\\ \·& \·& \·& \·\\ \·& \·& \·& \·\\ \·& \·& \·& \·\\ m& -1& -A\left(m\right)& -B\left(m\right)\end{array}\right],$ $\mathrm{\θ}=\left[\begin{array}{c}K\\ \mathrm{K\τ}\\ T_1{T}_2\\ T_1+T_2\end{array}\right],$ $\mathrm{\Ω}=\left[\begin{array}{c}{\mathrm{\Ω}}_1\\ {\mathrm{\Ω}}_2\\ \·\\ \·\\ \·\\ {\mathrm{\Ω}}_m\end{array}\right]$

差分方程中θ未知，通过测量C和U来估计未知参数θ，即被辨识的参数K，T₁，T₂和τ，为了估计θ的n个参数，必须有m≥n；

步骤三、利用最小二乘法辨识二阶加滞后最小二乘辨识模型方程中的未知参数

由于存在模型误差和观测噪声，未知参数θ的估计要选择最小误差平方最小性能指标，定义误差向量平方J＝Ω^TΩ，由式(6)可得：

J＝(C-Uθ)^T(C-Uθ)

估计

使误差平方J最小，即使 $\frac{\∂J}{\∂\mathrm{\θ}}=0,$ 得最小二乘参数估计值：

$\widehat{\mathrm{\θ}}={\left(U^{T}U\right)}^{-1}{U}^{T}C---\left(7\right)$

当未知参数估计值

取得后，可得对象的传递函数模型(二阶加滞后对象的模型结构)：

$G\left(s\right)=\frac{K}{(T_{1}s+1)(T_{2}s+1)}{e}^{-\mathrm{\τs}}---\left(8\right)$

步骤四、经离散化后取得供热对象的差分方程

将步骤三中得到的传递函数模型G(s)进行z变换，离散化后得到脉冲传递函数：

$G\left(z\right)=\frac{Y\left(z\right)}{R\left(z\right)}=\frac{b_{m^{\′}}{z}^{m^{\′}}+b_{m^{\′}-1}{z}^{m^{\′}-1}+\·\·\·+b_0}{a_{n^{\′}}{z}^{n^{\′}}+a_{n^{\′}-1}{z}^{n^{\′}-1}+\·\·\·+a_0}{z}^{-k_0}---\left(9\right)$

其中z为z变换算子；Y(z)和R(z)分别为对象的输出与输入z变换；m’，n’分别为分子、分母的阶数，且m′≥0，n′≥0；a₀，a₁，a₂，…，a_n′与b₀，b₁，b₂，…，b_m′分别是分母、分子多项式系数；k₀为输出纯滞后步数，且k₀≥0；

将所述脉冲传递函数的分子、分母同除以zⁿ′并进行数学变换，再将所得的式子中的Y(z)替换为y(k)、z^-iY(z)替代为y(k-i)、z^-jR(z)替换为r(k-j)，最终得到对象的差分方程：

$y\left(k\right)=[-\underset{i=0}{\overset{n^{\′}-1}{\mathrm{\Σ}}}{a}_{i}y(k-n^{\′}+i)+\underset{j=0}{\overset{m^{\′}}{\mathrm{\Σ}}}{b}_{j}r(k-n^{\′}+j-k_0)]/a_{n^{\′}}---\left(10\right)$

其中k为第k个时间点。

本发明的有益效果：

本发明所述建模方法可代替常规的飞升曲线法，可获高精度模型。模型精度相对常规的飞升曲线法建模精度由7％提高到3％以内，而且可直接辨识对象滞后时间，进而代替常规方法中由实验测取。本方法亦可在过程控制中，整定控制器参数时应用。本方法对实验条件及现场环境要求低，只须取部分实验数据而不必观测全过程，减少由建模带来的负面效应。本方法自身具有滤波作用，减少干扰，提高建模可靠性。本发明所述建模方法直接为供热节能改造分析和综合控制系统应用。

附图说明

图1是本发明的脉冲信号响应曲线图(横坐标k表示时间点，单位是秒；纵坐标A(t)单位是℃/k)；图2是本发明的阶跃信号响应曲线图(横坐标k表示时间点，单位是秒；纵坐标B(t)单位是℃)；图3是本发明的斜坡信号响应曲线图(横坐标k表示时间点，单位是秒；纵坐标C(t)单位是℃·k)；图4是数据长度不同时的二阶阶跃响应曲线(1-m＝30，2-m＝50，3-m＝70)；图5是数据长度m＝70时，一般飞升曲线法与最小二乘法阶跃建模的曲线比较图(1-实测曲线，2-阶跃响应法曲线，3-最小二乘法曲线)。

具体实施方式

具体实施方式一：本实施方式所述的热力站节能改造控制对象建模方法是按照以下步骤实现的：

步骤一、取得单位脉冲信号、单位阶跃信号和单位斜坡信号的响应模型：

供热过程采用二阶加滞后来描述，二阶加滞后对象的模型结构为 $G\left(s\right)=\frac{K}{(T_{1}s+1)(T_{2}s+1)}{e}^{-\mathrm{\τs}},$ 其被辨识的参数为K，T₁，T₂和τ；其中s为拉普拉斯算子，单位为秒^-1；K为比例系数；T₁，T₂均为时间常数；τ为滞后时间常数；

以下分别计算二阶系统的单位脉冲信号、单位阶跃信号、单位斜坡信号的响应模型，并分别令其为A(t)，B(t)，C(t)；当单位脉冲信号的响应模型A(t)、单位阶跃信号的响应模型B(t)、单位斜坡信号的响应模型C(t)三者之一由实验获取，则另两个响应模型可通过式(1)求得：

$\left\{\begin{array}{c}B\left(t\right)={\∫}_{t_0}^{t}A\left(\mathrm{\ξ}\right)\mathrm{d\ξ}\\ C\left(t\right)={\∫}_{t_0}^{t}B\left(\mathrm{\ξ}\right)\mathrm{d\ξ}\end{array}\right.---\left(1\right)$

热力站供热对象的模型为 $G\left(s\right)=\frac{K}{(T_{1}s+1)(T_{2}s+1)}{e}^{-\mathrm{\τs}}$ 时，其中令T₂＝βT₁，β为比例系数，则单位脉冲信号的响应模型为：

$A\left(t\right)=\frac{K}{T_1(1-\mathrm{\β})}[\exp (-\frac{t-\mathrm{\τ}}{T_1})-\exp (-\frac{t-\mathrm{\τ}}{\mathrm{\β}{T}_1})]+\mathrm{\ω}\left(t\right)---\left(2\right)$

上式中ω(t)为噪声项，t为时间，以下各式同；

单位阶跃信号的响应模型为：

$B\left(t\right)={\∫}_{t_0}^{t}A\left(\mathrm{\ξ}\right)\mathrm{d\ξ}=K+\frac{K}{1-\mathrm{\β}}[\mathrm{\βexp}(-\frac{t-\mathrm{\τ}}{\mathrm{\β}{T}_1})-\exp (-\frac{t-\mathrm{\τ}}{T_1})]+{\∫}_{t_0}^t\mathrm{\ω}\left(\mathrm{\ξ}\right)\mathrm{d\ξ}---\left(3\right)$

单位斜坡信号的响应模型为：

$C\left(t\right)={\∫}_{t_0}^{t}B\left(\mathrm{\ξ}\right)\mathrm{d\ξ}$

$=K(t-\mathrm{\τ})-(\mathrm{\β}+1){\mathrm{KT}}_1+\frac{{\mathrm{KT}}_1}{1-\mathrm{\β}}[\exp (-\frac{t-\mathrm{\τ}}{T_1})-{\mathrm{\β}}^2\exp (-\frac{t-\mathrm{\τ}}{\mathrm{\β}{T}_1})]+\∫\∫\mathrm{\ω}\left(\mathrm{\ξ}\right)\mathrm{d\ξd\σ}---\left(4\right)$

步骤二、由典型传递函数取得辨识方程

由步骤一中所述的式(2)、式(3)和式(4)

可求得供热对象的二阶加滞后最小二乘辨识模型方程：

$C\left(t\right)=\left[\begin{array}{cccc}t& -1& -A\left(t\right)& -B\left(t\right)\end{array}\right]\left[\begin{array}{c}K\\ \mathrm{K\τ}\\ T_1{T}_2\\ T_1+T_2\end{array}\right]+\mathrm{\Ω}\left(t\right)---\left(5\right)$

其中Ω(t)为噪声项，是对ω(t)的滤波，有抑制噪声的作用；

取不同时刻t₁，t₂，…，t_m，可得线性差分方程组

C＝Uθ+Ω                            (6)

其中， $C=\left[\begin{array}{c}C\left(1\right)\\ C\left(2\right)\\ \·\\ \·\\ \·\\ C\left(m\right)\end{array}\right],$ $U=\left[\begin{array}{cccc}1& -1& -A\left(1\right)& -B\left(1\right)\\ 2& -1& -A\left(2\right)& -B\left(2\right)\\ \·& \·& \·& \·\\ \·& \·& \·& \·\\ \·& \·& \·& \·\\ m& -1& -A\left(m\right)& -B\left(m\right)\end{array}\right],$ $\mathrm{\θ}=\left[\begin{array}{c}K\\ \mathrm{K\τ}\\ T_1{T}_2\\ T_1+T_2\end{array}\right],$ $\mathrm{\Ω}=\left[\begin{array}{c}{\mathrm{\Ω}}_1\\ {\mathrm{\Ω}}_2\\ \·\\ \·\\ \·\\ {\mathrm{\Ω}}_m\end{array}\right]$

差分方程中θ未知，通过测量C和U来估计未知参数θ，即被辨识的参数K，T₁，T₂和τ，为了估计θ的n个参数，必须有m≥n；

步骤三、利用最小二乘法辨识二阶加滞后最小二乘辨识模型方程中的未知参数

由于存在模型误差和观测噪声，未知参数θ的估计要选择最小误差平方最小性能指标，定义误差向量平方J＝Ω^TΩ，由式(6)可得：

J＝(C-Uθ)^T(C-Uθ)

估计

使误差平方J最小，即使 $\frac{\∂J}{\∂\mathrm{\θ}}=0,$ 得最小二乘参数估计值：

$\widehat{\mathrm{\θ}}={\left(U^{T}U\right)}^{-1}{U}^{T}C---\left(7\right)$

当未知参数估计值

取得后，可得对象的传递函数模型(二阶加滞后对象的模型结构)：

$G\left(s\right)=\frac{K}{(T_{1}s+1)(T_{2}s+1)}{e}^{-\mathrm{\τs}}---\left(8\right)$

步骤四、经离散化后取得供热对象的差分方程

将步骤三中得到的传递函数模型G(s)进行z变换，离散化后得到脉冲传递函数：

$G\left(z\right)=\frac{Y\left(z\right)}{R\left(z\right)}=\frac{b_{m^{\′}}{z}^{m^{\′}}+b_{m^{\′}-1}{z}^{m^{\′}-1}+\·\·\·+b_0}{a_{n^{\′}}{z}^{n^{\′}}+a_{n^{\′}-1}{z}^{n^{\′}-1}+\·\·\·+a_0}{z}^{-k_0}---\left(9\right)$

其中z为z变换算子；Y(z)和R(z)分别为对象的输出与输入z变换；m’，n’分别为分子、分母的阶数，且m′≥0，n′≥0；a₀，a₁，a₂，…，a_n′与b₀，b₁，b₂，…，b_m′分别是分母、分子多项式系数；k₀为输出纯滞后步数，且k₀≥0；

将所述脉冲传递函数的分子、分母同除以zⁿ′并进行数学变换，再将所得的式子中的Y(z)替换为y(k)、z^-iY(z)替代为y(k-i)、z^-jR(z)替换为r(k-j)，最终得到对象的差分方程：

$y\left(k\right)=[-\underset{i=0}{\overset{n^{\′}-1}{\mathrm{\Σ}}}{a}_{i}y(k-n^{\′}+i)+\underset{j=0}{\overset{m^{\′}}{\mathrm{\Σ}}}{b}_{j}r(k-n^{\′}+j-k_0)]/a_{n^{\′}}---\left(10\right)$

其中k为第k个时间点。

具体建模应用实例：

如图1～5所示，以供热站二级网回水温度为输出量，一级网阀门开度为输入量，建立此对象模型为例，说明本发明的实施方案：

1.获取斜坡、阶跃、脉冲响应曲线：

测取阀门开度阶跃变化回水温度响应曲线，由此得斜坡、阶跃、脉冲响应，参见图1～3，表1为典型响应数值，

表1典型响应数值

2.选择样本个数m，列写差分方程：

在表1中取A，B，C，t构成差分方程：

$U=\left[\begin{array}{cccc}t\left(1\right)& -1& -A\left(1\right)& -B\left(1\right)\\ t\left(2\right)& -1& -A\left(2\right)& -B\left(2\right)\\ \·& \·& \·& \·\\ \·& \·& \·& \·\\ \·& \·& \·& \·\\ t\left(m\right)& -1& -A\left(m\right)& -B\left(m\right)\end{array}\right],$ $C=\left[\begin{array}{c}C\left(1\right)\\ C\left(2\right)\\ \·\\ \·\\ \·\\ C\left(m\right)\end{array}\right]$

3.利用最小二乘法求模型参数K，T₁，T₂和τ：

根据最小二乘法可以求得 $\widehat{\mathrm{\θ}}={\left(U^{T}U\right)}^{-1}{U}^{T}C,$

其中， $\widehat{\mathrm{\θ}}=\left[\begin{array}{c}K\\ \mathrm{K\τ}\\ T_1{T}_2\\ T_1+T_2\end{array}\right],$ $U=\left[\begin{array}{cccc}t\left(1\right)& -1& -A\left(1\right)& -B\left(1\right)\\ t\left(2\right)& -1& -A\left(2\right)& -B\left(2\right)\\ \·& \·& \·& \·\\ \·& \·& \·& \·\\ \·& \·& \·& \·\\ t\left(m\right)& -1& -A\left(m\right)& -B\left(m\right)\end{array}\right],$ $C=\left[\begin{array}{c}C\left(1\right)\\ C\left(2\right)\\ \·\\ \·\\ \·\\ C\left(m\right)\end{array}\right]$

从而可以得到参数K，T₁，T₂和τ。

4.取得回水温度变化的传递函数模型：

对m分别取30，50，70进行辨识，得模型为

$G_1\left(s\right)=\frac{9.97}{(12.64s+1)(34.49s+1)}{e}^{-9.56s}$

$G_2\left(s\right)=\frac{9.94}{(13.14s+1)(33.45s+1)}{e}^{-9.6s}$

$G_3\left(s\right)=\frac{9.86}{(14.43s+1)(31.25s+1)}{e}^{-9.57s}$

由图4可知此辨识法对数据长度不敏感，可不必选择完全进入稳态的数据。

由传递函数得差分方程模型：

由第4步可知采样个数m＝70时的供热过程动态模型传递函数：

$G\left(s\right)=\frac{9.86}{(14.43s+1)(31.25s+1)}{e}^{-9.57s},$

则其脉冲传递函数为：

$G\left(z\right)=z^{-9}\×\frac{0.001992z^2+0.01553z+0.003272}{z^3-1.902z^2+0.9037}$

将分子、分母同除以z³：

$G\left(z\right)=z^{-9}\×\frac{0.001992z^{-1}+0.01553z^{-2}+0.003272z^{-3}}{1-1.902z^{-1}+0.9037z^{-3}}$

根据 $G\left(z\right)=\frac{Y\left(z\right)}{R\left(z\right)}$ 可得：

(1-1.902z^-1+0.9037z^-3)Y(z)＝(0.001992z^-10+0.01553z^-11+0.003272z^-12)R(z)替换得到模型的差分方程：

y(k)＝1.902y(k-1)-0.9037y(k-3)+0.001992r(k-10)+

      0.01553r(k-11)+0.003272r(k-12)

5.与一般飞升曲线建模法对比(参见图5和表2)：

利用较高精度的阶跃响应曲线法给出二阶模型参数，K，τ直接由阶跃响应曲线得到，而T₁，T₂由式(11)解出，

$\left\{\begin{array}{c}{t}_1-\frac{T_1{T}_2}{T_1-T_2}\ln \frac{T_1}{T_2}=T_1+T_2\\ t_1-t_0=T_1{\left(\frac{T_1}{T_2}\right)}^{\frac{T_2}{T_1-T_2}}\end{array}\right.---\left(11\right)$

取得阶跃响应辨识模型为：

$G\left(s\right)=\frac{10.06}{(10.57s+1)(36.87s+1)}{e}^{-6.7s}$

表2两种方法精度比较

表2给出了本发明法与一般飞升曲线的误差比较。经比较说明，本发明法的精度要高于常规方法的辨识精度。其中精度用 $\frac{1}{n}\underset{k=m}{\overset{m+n}{\mathrm{\Σ}}}|y\left(k\right)-\hat{y}\left(k\right)|/y\left(k\right)\×100\%$ 表示。

Claims (1)

1、一种热力站节能改造控制对象建模方法，其特征在于所述建模方法是按照以下步骤实现的：

步骤一、取得单位脉冲信号、单位阶跃信号和单位斜坡信号的响应模型：

供热过程采用二阶加滞后来描述，二阶加滞后对象的模型结构为 $G\left(s\right)=\frac{K}{(T_{1}s+1)(T_{2}s+1)}{e}^{-\mathrm{\τs}},$ 其被辨识的参数为K，T₁，T₂和τ；其中s为拉普拉斯算子，单位为秒^-1；K为比例系数；T₁，T₂均为时间常数；τ为滞后时间常数；

以下分别计算二阶系统的单位脉冲信号、单位阶跃信号、单位斜坡信号的响应模型，并分别令其为A(t)，B(t)，C(t)；当单位脉冲信号的响应模型A(t)、单位阶跃信号的响应模型B(t)、单位斜坡信号的响应模型C(t)三者之一由实验获取，则另两个响应模型可通过式(1)求得：

$\left\{\begin{array}{c}B\left(t\right)={\∫}_{t_0}^{t}A\left(\mathrm{\ξ}\right)\mathrm{d\ξ}\\ C\left(t\right)={\∫}_{t_0}^{t}B\left(\mathrm{\ξ}\right)\mathrm{d\ξ}\end{array}\right.---\left(1\right)$

热力站供热对象的模型为 $G\left(s\right)=\frac{K}{(T_{1}s+1)(T_{2}s+1)}{e}^{-\mathrm{\τs}}$ 时，其中令T₂＝βT₁，β为比例系数，则单位脉冲信号的响应模型为：

$A\left(t\right)=\frac{K}{T_1(1-\mathrm{\β})}[\exp (-\frac{t-\mathrm{\τ}}{T_1})-\exp (-\frac{t-\mathrm{\τ}}{\mathrm{\β}{T}_1})]+\mathrm{\ω}\left(t\right)---\left(2\right)$

上式中ω(t)为噪声项，t为时间，以下各式同；

单位阶跃信号的响应模型为：

$B\left(t\right)={\∫}_{t_0}^{t}A\left(\mathrm{\ξ}\right)\mathrm{d\ξ}=K+\frac{K}{1-\mathrm{\β}}[\mathrm{\βexp}(-\frac{t-\mathrm{\τ}}{\mathrm{\β}{T}_1})-\exp (-\frac{t-\mathrm{\τ}}{T_1})]+{\∫}_{t_0}^t\mathrm{\ω}\left(\mathrm{\ξ}\right)\mathrm{d\ξ}---\left(3\right)$

单位斜坡信号的响应模型为：

$C\left(t\right)={\∫}_{t_0}^{t}B\left(\mathrm{\ξ}\right)\mathrm{d\ξ}$

$=K(t-\mathrm{\τ})-(\mathrm{\β}+1)K{T}_1+\frac{{\mathrm{KT}}_1}{1-\mathrm{\β}}[\exp (-\frac{t-\mathrm{\τ}}{T_1})-{\mathrm{\β}}^2\exp (-\frac{t-\mathrm{\τ}}{\mathrm{\β}{T}_1})]+\∫\∫\mathrm{\ω}\left(\mathrm{\ξ}\right)\mathrm{d\ξd\σ}---\left(4\right)$

步骤二、由典型传递函数取得辨识方程

由步骤一中所述的式(2)、式(3)和式(4)

可求得供热对象的二阶加滞后最小二乘辨识模型方程：

$C\left(t\right)=\left[\begin{array}{cccc}t& -1& -A\left(t\right)& -B\left(t\right)\end{array}\right]\left[\begin{array}{c}K\\ \mathrm{K\τ}\\ T_1{T}_2\\ T_1+T_2\end{array}\right]+\mathrm{\Ω}\left(t\right)---\left(5\right)$

其中Ω(t)为噪声项，是对ω(t)的滤波；

取不同时刻t₁，t₂，…，t_m，可得线性差分方程组

C＝Uθ+Ω                (6)

其中， $C=\left[\begin{array}{c}C\left(1\right)\\ C\left(2\right)\\ \·\\ \·\\ \·\\ C\left(m\right)\end{array}\right],$ $U=\left[\begin{array}{cccc}1& -1& -A\left(1\right)& -B\left(1\right)\\ 2& -1& -A\left(2\right)& -B\left(2\right)\\ \·& \·& \·& \·\\ \·& \·& \·& \·\\ \·& \·& \·& \·\\ m& -1& -A\left(m\right)& -B\left(m\right)\end{array}\right],$ $\mathrm{\θ}=\left[\begin{array}{c}K\\ \mathrm{K\τ}\\ T_1{T}_2\\ T_1+T_2\end{array}\right],$ $\mathrm{\Ω}=\left[\begin{array}{c}{\mathrm{\Ω}}_1\\ {\mathrm{\Ω}}_2\\ \·\\ \·\\ \·\\ {\mathrm{\Ω}}_m\end{array}\right]$

差分方程中θ未知，通过测量C和U来估计未知参数θ，即被辨识的参数K，T₁，T₂和τ，为了估计θ的n个参数，必须有m≥n；

步骤三、利用最小二乘法辨识二阶加滞后最小二乘辨识模型方程中的未知参数

由于存在模型误差和观测噪声，未知参数θ的估计要选择最小误差平方最小性能指标，定义误差向量平方J＝Ω^TΩ，由式(6)可得：

J＝(C-Uθ)^T(C-Uθ)

估计

，使误差平方J最小，即使 $\frac{\∂J}{\∂\mathrm{\θ}}=0,$ 得最小二乘参数估计值：

$\widehat{\mathrm{\θ}}={\left(U^{T}U\right)}^{-1}{U}^{T}C---\left(7\right)$

当未知参数估计值

取得后，可得对象的传递函数模型：

$G\left(s\right)=\frac{K}{(T_{1}s+1)(T_{2}s+1)}{e}^{-\mathrm{\τs}}---\left(8\right)$

步骤四、经离散化后取得供热对象的差分方程

将步骤三中得到的传递函数模型G(s)进行z变换，离散化后得到脉冲传递函数：

$G\left(z\right)=\frac{Y\left(z\right)}{R\left(z\right)}=\frac{b_{m^{\′}}{z}^{m^{\′}}+b_{m^{\′}-1}{z}^{m^{\′}-1}+\·\·\·+b_0}{a_{n^{\′}}{z}^{n^{\′}}+a_{n^{\′}-1}{z}^{n^{\′}-1}+\·\·\·+a_0}{z}^{-k_0}---\left(9\right)$

其中z为z变换算子；Y(z)和R(z)分别为对象的输出与输入z变换；m’，n’分别为分子、分母的阶数，且m′≥0，n′≥0；a₀，a₁，a₂，…，a_n′与b₀，b₁，b₂，…，b_m′分别是分母、分子多项式系数；k₀为输出纯滞后步数，且k₀≥0；

将所述脉冲传递函数的分子、分母同除以zⁿ′并进行数学变换，再将所得的式子中的Y(z)替换为y(k)、z^-iY(z)替代为y(k-i)、z^-jR(z)替换为r(k-j)，最终得到对象的差分方程：

$y\left(k\right)=[-\underset{i=0}{\overset{n^{\′}-1}{\mathrm{\Σ}}}{a}_{i}y(k-n^{\′}+i)+\underset{j=0}{\overset{m^{\′}}{\mathrm{\Σ}}}{b}_{j}r(k-n^{\′}+j-k_0)]/a_{n^{\′}}---\left(10\right)$

其中k为第k个时间点。

Images