CN104850008A - 一种层燃锅炉的多变量控制与实时优化的方法 - Google Patents

Abstract

本发明提供了一种层燃锅炉的多变量控制与实时优化的方法，包括如下步骤：1)将层燃锅炉的控制转化为3变量输入和4变量输出的多变量控制问题；2)将所述多变量控制问题拆分为两个带前馈反馈功能的多变量控制器相互耦合的控制问题；对所述两个带前反馈功能的多变量控制问题分别进行状态空间升阶转换，将当前时刻输入信号的增量作为新的输入信号进行控制；3)对经过状态空间转换后的两个带前馈反馈功能的多变量控制问题基于当前时刻拆分成过去和将来作用两部分，并转化为求带有约束条件的二次规划最优解问题；4)通过锅炉排速度区间，划分不同的层燃锅炉模型组，通过模型组自动切换解决层燃锅炉的强非线性问题。

Description

一种层燃锅炉的多变量控制与实时优化的方法

技术领域

本发明涉及一种层燃锅炉的多变量控制与实时优化的方法。

背景技术

由于历史原因，层燃锅炉是目前我国民间和工业界应用最为广泛的锅炉类型之一。然而，该类锅炉不仅燃烧效率低，而且对环境影响特别恶劣。近些年来，随着原材价格不断高升，以及我国在低碳环保方面的要求也越来越严苛，如何提高层燃锅炉的燃烧效率，已成为社会各界广为关注的问题之一。目前我国层燃锅炉的燃烧效率远低于国际平均水平。其主要原因是，目前层燃锅炉的燃烧熄灭界限还主要依靠人工通过裸眼观测手动完成。然而熄灭界限的误差将直接关系到层燃锅炉的燃烧效率。一来，如果熄灭点离炉膛前端太近，会造成燃料燃烧不充分，而且大大降低了锅炉的有效燃烧面积，使得炉膛温度很难达到预定设定值。为了保证足够的热源，锅炉不得不延长运行时间，结果造成进一步的能源浪费；另一方面，如果熄灭界限离炉膛后端太近，不仅燃料无法充分燃烧；而且，温度过高的燃料层在离开炉排末端时，不仅会严重损坏炉排末端设备，而且过高温度的炉渣在进入渣槽后遇水会产生大量的酸性气体。挥发后的酸性气体会直接腐蚀炉膛内部机构，从而大大降低了锅炉的使用寿命。因此，本专利提出了一种锅炉控制系统，它可以使炉膛燃烧熄灭界线始终稳定工作在最佳区域，便于充分利用炉膛内燃料，减少酸性气体的产生，保护炉膛内部设备。同时，在确保供暖品质的同时，使得燃料处于最合理、最节省的供给工况。

层燃锅炉系统是一个存在大时延、非线性、强耦合的系统，经线性化处理后，系统易体现出高阶、非最小相位特性。为了提高层燃锅炉的实时控制效果，目前业界比较推崇的控制方法为：基于多变量MPC(Model Predictive Control，模型预测控制)算法的APC(Advanced Process Control，先进过程控制)技术。然而，目前MPC在实际工程应用中存在诸多问题，比如锅炉的大时延问题、强非线性问题、输入/输出强耦合问题、关键参数难确定(如熄灭界限)等。该发明专利提供了一套高效、智能、应用广泛的层燃锅炉多变量控制与实时优化策略，以解决目前存在的诸多问题。

发明内容

本发明所要解决主要技术问题是提供一种层燃锅炉的多变量控制与实时优化的方法，控制过程迅速、稳定、抗干扰性强。

为了解决上述的技术问题，本发明提供了一种层燃锅炉的多变量控制与实时优化的方法，其特征在于包括如下步骤：

1)将层燃锅炉的控制转化为3变量输入和4变量输出的多变量控制问题；所述3变量输入包括引风转速u₁、鼓风转速u₂、炉排速度u₃；所述4变量输出包括炉膛温度y₁、炉膛负压y₂、含氧量y₃、火焰熄灭位置y₄；

2)将所述多变量控制问题拆分为两个带前馈反馈功能的多变量控制器相互耦合的控制问题；其中，第一控制器的反馈控制变量在第二控制器中以前馈预测变量的形式出现；第二控制器的反馈控制变量在第一控制器中以前馈预测变量的形式出现；对所述两个带前反馈功能的多变量控制问题分别进行状态空间升阶转换，将当前时刻输入信号的增量作为新的输入信号进行控制；

3)对经过状态空间转换后的两个带前馈反馈功能的多变量控制问题基于 当前时刻拆分成过去和将来作用两部分，并转化为求带有约束条件的二次规划最优解问题；

4)通过锅炉排速度区间，划分不同的层燃锅炉模型组，通过模型组自动切换解决层燃锅炉的强非线性问题。

在一较佳实施例中：所述3变量输入和4变量输出的多变量控制问题具体表达为：

$\left[\begin{array}{c}{y}_1\\ y_2\\ y_3\\ y_4\end{array}\right]=\left[\begin{array}{ccc}{g}_{11}& g_{12}& g_{13}\\ g_{21}& g_{22}& g_{23}\\ g_{31}& g_{32}& g_{33}\\ g_{41}& g_{42}& g_{43}\end{array}\right]\left[\begin{array}{c}{u}_1\\ u_2\\ u_3\end{array}\right]---\left(1\right)$

g_ij表示第j个操控变量与第i个被控变量间的传递函数模型；由于炉排速度u3对炉膛负压y2、含氧量y3的影响小，因此g₂₃、g₃₃这里取值为0，即：

$\left[\begin{array}{c}{y}_1\\ y_2\\ y_3\\ y_4\end{array}\right]=\left[\begin{array}{ccc}{g}_{11}& g_{12}& g_{13}\\ g_{21}& g_{22}& 0\\ g_{31}& g_{32}& 0\\ g_{41}& g_{42}& g_{43}\end{array}\right]\left[\begin{array}{c}{u}_1\\ u_2\\ u_3\end{array}\right]---\left(2\right)$

在一较佳实施例中：所述两个带前馈反馈功能的多变量控制器控制问题具体表达为：

$\left[\begin{array}{c}{y}_1\\ y_2\\ y_3\end{array}\right]=\left[\begin{array}{ccc}{g}_{11}& g_{12}& g_{13}\\ g_{21}& g_{22}& 0\\ g_{31}& g_{32}& g_0\end{array}\right]\left[\begin{array}{c}{u}_1\\ u_2\\ u_3\end{array}\right]=\left[\begin{array}{cc}{g}_{11}& g_{12}\\ g_{21}& g_{22}\\ g_{31}& g_{32}\end{array}\right]\left[\begin{array}{c}{u}_1\\ u_2\end{array}\right]+\left[\begin{array}{c}{g}_{13}\\ 0\\ 0\end{array}\right]u_3---\left(3\right)$

即：Y₁＝G₁U₁+G_d1D₁     (4) 

$y_4=\left[\begin{array}{ccc}{g}_{41}& g_{42}& g_{43}\end{array}\right]\left[\begin{array}{c}{u}_1\\ u_2\\ u_3\end{array}\right]=g_{43}{u}_3+\begin{array}{}\end{array}\left[\begin{array}{cc}{g}_{41}& g_{42}\end{array}\right]\left[\begin{array}{c}{u}_1\\ u_2\end{array}\right]---\left(5\right)$

即：Y₂＝G₂U₂+G_d2D₂    (6) 

其中，Y₁、Y₂为被控矢量，U₁、U₂为输入矢量，D₁、D₂为扰动矢量； 

比较上式可见，引风转速u₁、鼓风转速u₂与炉膛温度y₁、炉膛负压y₂、含 氧量y₃形成闭环反馈控制回路，而炉排速度u₃与炉膛温度y₁形成的前馈预测控制回路；同理，炉排速度u₃与火焰熄灭位置y₄形成闭环反馈控制回路，而引风转速u₁、鼓风转速u₂与火焰熄灭位置y₄形成的前馈预测控制回路。

将式4和式6进行归纳，可以写出一类控制的通式，

Y＝GU+G_dD     (7)

其中，Y是m维被控矢量，G是m*n的控制矩阵，U是n维输入矢量，G_d是m*q的扰动矩阵，D是q维扰动矢量，这里的m>0，n>0，q≥0，且都为整数；

首先将G和G_d转换为空间状态表达式:

Y(k)＝GU(k-T_d1)+G_dD(k-T_d2)   (8) 

进一步写为：Y(k)＝Y′+Y″    (9) 

其中T_d1为输入最大时延，T_d2为扰动最大时延；Y′为操控响应，Y″为扰动响应；

对Y′、Y″分别进行状态空间转换：

首先，对Y′进行状态空间转换，可得：

$\left\{\begin{array}{c}x(k+1)=\mathrm{Ax}\left(k\right)+\mathrm{BU}(k-T_{d1})\\ Y^{\′}\left(k\right)=\mathrm{Cx}\left(k\right)\end{array}\right.---\left(10\right)$

其中x(k)为状态变量，U(k-T_d1)为输入变量，A为系统矩阵、B为输入矩阵、C为输出矩阵。本控制方案的特点在于，对式10进行了进一步状态升阶操作，状态空间表达式10可以转化为：

$\left\{\begin{array}{c}\stackrel{\‾}{X}(k+1)=\stackrel{\‾}{A}\stackrel{\‾}{X}\left(k\right)+\stackrel{\‾}{B}\mathrm{\ΔU}(k-T_{d1})\\ {\stackrel{\‾}{Y}}^{\′}\left(k\right)=\stackrel{\‾}{C}\stackrel{\‾}{X}\left(k\right)\end{array}\right.---\left(11\right)$

其中状态变量x(k)被升阶为： 

$\stackrel{\‾}{X}\left(k\right)=\left[\begin{array}{c}\mathrm{\Δx}\left(k\right)\\ x\left(k\right)\end{array}\right],\mathrm{\Δx}\left(k\right)=x\left(k\right)-x(k-1).$

升阶状态矩阵由原状态空间表达式10中的系数矩阵(A,B,C)构成；

同理，对Y″进行状态空间转换，可得：

$\left\{\begin{array}{c}{x}_d(k+1)=A_d{x}_d\left(k\right)+B_{d}D(k-T_{d2})\\ Y^{\′\′}=C_d{x}_d\left(k\right)\end{array}\right.---\left(12\right)$

状态空间表达式12可被升阶为： 

$\left\{\begin{array}{c}{\stackrel{\‾}{X}}_d(k+1)={\stackrel{\‾}{A}}_d{\stackrel{\‾}{X}}_d\left(k\right)+{\stackrel{\‾}{B}}_d\mathrm{\ΔD}(k-T_{d2})\\ {\stackrel{\‾}{Y}}^{\′\′}\left(k\right)={\stackrel{\‾}{C}}_d{\stackrel{\‾}{X}}_d\left(k\right)\end{array}\right.---\left(13\right)$

其中状态x(k)被升阶为： 

${\stackrel{\‾}{X}}_d=\left[\begin{array}{c}\mathrm{\Δ}{x}_d\left(k\right)\\ x_d\left(k\right)\end{array}\right],\mathrm{\Δ}{x}_d\left(k\right)=x_d\left(k\right)-x_d(k-1).$

升阶状态矩阵由原状态空间表达式12中的系数矩阵(A,B,C)构成；是一个满足因果关系的扰动,在预测周期内该扰动保持不变；

这里将增量ΔU(k)＝U(k)-U(k-1)系统作为输入信号(包括控制变量与扰动变量)，并结合状态变量进行控制。

在一较佳实施例中：所述步骤3的具体过程为:

设定预测过程中的当前时刻为k，给定输入的控制步长为H_u，即输入序列的时间跨度为：ΔU_k→ΔU_k+Hu-1；输出的预测步长H_p，即从预测当前时刻的下一步开始，输出轨迹的时间跨度为Y_k+1→Y_k+Hp；

将式11的操控输入信号和式13的扰动输入信号基于当前时刻k拆分成过去和将来作用两部分，转换后的状态空间表达式如下所示：

$\left\{\begin{array}{c}\tilde{X}\left(k\right)=\tilde{A}\stackrel{\‾}{X}\left(k\right)+{\tilde{B}}_p\mathrm{\Δ}{\tilde{U}}_p\left(k\right)+{\tilde{B}}_f\mathrm{\Δ}{\tilde{U}}_f\left(k\right)\\ {\tilde{X}}_d\left(k\right)={\tilde{A}}_d{\stackrel{\‾}{X}}_d\left(k\right)+{\tilde{B}}_{\mathrm{dp}}\mathrm{\Δ}{\tilde{D}}_p\left(k\right)+{\tilde{B}}_{\mathrm{df}}\mathrm{\Δ}{\tilde{D}}_f\left(k\right)\end{array}\right.---\left(14\right)$

和分别指过去的操控信号输入和扰动信号输入的作用部分，和分别指未来预测的操控信号输入和扰动信号输入的作用部分；

其中，操控信号输入的各个系数详细表达为：

$\tilde{X}\left(k\right)=\left[\begin{array}{c}\stackrel{\‾}{X}(k+1)\\ \stackrel{\‾}{X}(k+2)\\ .\\ .\\ .\\ \stackrel{\‾}{X}(k+\mathrm{Hp})\end{array}\right]$

$\stackrel{\‾}{A}=\left[\begin{array}{c}\stackrel{\‾}{A}\\ {\stackrel{\‾}{A}}^2\\ .\\ .\\ .\\ {\stackrel{\‾}{A}}^{\mathrm{Hp}}\end{array}\right]$

${\tilde{B}}_p=\left[\begin{array}{cccccc}\stackrel{\‾}{B}& 0& ...& ...& ...& 0\\ \stackrel{\‾}{A}\stackrel{\‾}{B}& \stackrel{\‾}{B}& 0& ...& ...& 0\\ {\stackrel{\‾}{A}}^2& \stackrel{\‾}{A}\stackrel{\‾}{B}& \stackrel{\‾}{B}& 0& ...& 0\\ .& .& .& .& .& .\\ .& .& .& .& .& .\\ .& .& .& .& .& .\\ {\stackrel{\‾}{A}}^{\mathrm{Td}-1}\stackrel{\‾}{B}& {\stackrel{\‾}{A}}^{\mathrm{Td}-2}\stackrel{\‾}{B}& ...& ...& ...& \stackrel{\‾}{B}\\ {\stackrel{\‾}{A}}^{\mathrm{Td}}\stackrel{\‾}{B}& {\stackrel{\‾}{A}}^{\mathrm{Td}-1}\stackrel{\‾}{B}& ...& ...& ...& \stackrel{\‾}{A}\stackrel{\‾}{B}\\ .& .& .& .& .& .\\ .& .& .& .& .& .\\ .& .& .& .& .& .\\ {\stackrel{\‾}{A}}^{\mathrm{Hp}-1}\stackrel{\‾}{B}& {\stackrel{\‾}{A}}^{\mathrm{Hp}-2}\stackrel{\‾}{B}& ...& ...& ...& {\stackrel{\‾}{A}}^{\mathrm{Hp}-\mathrm{Td}-2}\stackrel{\‾}{B}\end{array}\right]$

$\mathrm{\Δ}{\tilde{U}}_p\left(k\right)=\left[\begin{array}{c}\mathrm{\ΔU}(k-T_d)\\ \mathrm{\ΔU}(k-T_d+1)\\ .\\ .\\ .\\ \mathrm{\ΔU}(k-1)\end{array}\right]$

$\mathrm{\Δ}{\tilde{U}}_f\left(k\right)=\left[\begin{array}{c}\mathrm{\ΔU}\left(k\right)\\ \mathrm{\ΔU}(k+1)\\ .\\ .\\ .\\ \mathrm{\ΔU}(k+\mathrm{Hu}-1)\end{array}\right]$

扰动信号输入的各个系数详细表达为：

${\tilde{X}}_d\left(k\right)=\left[\begin{array}{c}{\stackrel{\‾}{X}}_d(k+1)\\ {\stackrel{\‾}{X}}_d(k+2)\\ .\\ .\\ .\\ {\stackrel{\‾}{X}}_d(k+\mathrm{Hp})\end{array}\right]$

$\tilde{A}=\left[\begin{array}{c}{\stackrel{\‾}{A}}_d\\ {{\stackrel{\‾}{A}}_d}^2\\ .\\ .\\ .\\ {{\stackrel{\‾}{A}}_d}^{\mathrm{Hp}}\end{array}\right]$

$\mathrm{\Δ}{\tilde{D}}_p\left(k\right)=\left[\begin{array}{c}\mathrm{\ΔD}(k-T_d)\\ \mathrm{\ΔD}(k-T_d+1)\\ .\\ .\\ .\\ \mathrm{\ΔD}(k-1)\end{array}\right]$

${\tilde{B}}_{\mathrm{dp}}=\left[\begin{array}{cccccc}{\stackrel{\‾}{B}}_d& 0& ...& ...& ...& 0\\ {\stackrel{\‾}{A}}_d{\stackrel{\‾}{B}}_d& {\stackrel{\‾}{B}}_d& 0& ...& ...& 0\\ {{\stackrel{\‾}{A}}_d}^2{\stackrel{\‾}{B}}_d& {\stackrel{\‾}{A}}_d{\stackrel{\‾}{B}}_d& {\stackrel{\‾}{B}}_d& 0& ...& 0\\ .& .& .& .& .& .\\ .& .& .& .& .& .\\ .& .& .& .& .& .\\ {{\stackrel{\‾}{A}}_d}^{\mathrm{Td}-1}{\stackrel{\‾}{B}}_d& {{\stackrel{\‾}{A}}_d}^{\mathrm{Td}-2}{\stackrel{\‾}{B}}_d& ...& ...& ...& {\stackrel{\‾}{B}}_d\\ {{\stackrel{\‾}{A}}_d}^{\mathrm{Td}}{\stackrel{\‾}{B}}_d& {{\stackrel{\‾}{A}}_d}^{\mathrm{Td}-1}{\stackrel{\‾}{B}}_d& ...& ...& ...& {\stackrel{\‾}{A}}_d{\stackrel{\‾}{B}}_d\\ .& .& .& .& .& .\\ .& .& .& .& .& .\\ .& .& .& .& .& .\\ {{\stackrel{\‾}{A}}_d}^{\mathrm{Hp}-1}{\stackrel{\‾}{B}}_d& {{\stackrel{\‾}{A}}_d}^{\mathrm{Hp}-2}{\stackrel{\‾}{B}}_d& ...& ...& ...& {{\stackrel{\‾}{A}}_d}^{\mathrm{Hp}-\mathrm{Td}-2}{\stackrel{\‾}{B}}_d\end{array}\right]$

$\mathrm{\Δ}{\tilde{D}}_f\left(k\right)=\left[\begin{array}{c}\mathrm{\ΔD}\left(k\right)\\ \mathrm{\ΔD}(k+1)\\ .\\ .\\ .\\ \mathrm{\ΔD}(k+\mathrm{Hu}-1)\end{array}\right]$

在一较佳实施例中：将式14中的设为0，则式11、13可以简化为：

$\left\{\begin{array}{c}{\tilde{Y}}^{\′}\left(k\right)=\tilde{C}\tilde{A}\stackrel{\‾}{X}\left(k\right)+\tilde{C}{\tilde{B}}_p\mathrm{\Δ}{\tilde{U}}_p\left(k\right)+\tilde{C}{\tilde{B}}_f\mathrm{\Δ}{\tilde{U}}_f\left(k\right)\\ {\tilde{Y}}^{\′\′}\left(k\right)={\tilde{C}}_d{\tilde{A}}_d{\stackrel{\‾}{X}}_d\left(k\right)+{\tilde{C}}_d{\tilde{B}}_{\mathrm{dp}}\mathrm{\Δ}{\tilde{D}}_p\left(k\right)\end{array}\right.---\left(15\right)$

其中，分别是由构成的输出系数矩阵，可分别表示成：  $\tilde{C}=\mathrm{diag}\{\stackrel{\‾}{C},\stackrel{\‾}{C},...,\stackrel{\‾}{C}\}$ 和 ${\tilde{C}}_d=\mathrm{diag}\{{\stackrel{\‾}{C}}_d,{\stackrel{\‾}{C}}_d,...,{\stackrel{\‾}{C}}_d\}.$

通过式7和式15，可以总结得出：

$\begin{array}{c}\tilde{Y}\left(k\right)={\tilde{Y}}^{\′}\left(k\right)+{\tilde{Y}}^{\′\′}\left(k\right)\\ =\tilde{C}{\tilde{B}}_f\mathrm{\Δ}{\tilde{U}}_f\left(k\right)+[\tilde{C}\tilde{A}\stackrel{\‾}{X}\left(k\right)+\tilde{C}{\tilde{B}}_p\mathrm{\Δ}{\tilde{U}}_p\left(k\right)+{\tilde{C}}_d{\tilde{A}}_d{\stackrel{\‾}{X}}_d\left(k\right)+{\tilde{C}}_d{\tilde{B}}_{\mathrm{dp}}\mathrm{\Δ}{\tilde{D}}_p\left(k\right)]\end{array}---\left(16\right)$

式16由未知的未来输入响应项和已知的过去输入响应项组成，对这两部分进行整理，可将优化问题转化成式17所示的求带有约束条件的二次规划最优解问题：

$J={(\tilde{Y}\left(k\right)-Y_{\mathrm{tgt}})}^{T}Q(\tilde{Y}\left(k\right)-Y_{\mathrm{tgt}})+\mathrm{\Δ}{\tilde{U}}_f{\left(k\right)}^T\mathrm{R\Δ}{\tilde{U}}_f\left(k\right)$

$\begin{array}{cc}s.t.& \mathrm{\Σ\Δ}{\tilde{U}}_f\end{array}\left(k\right)\≤T---\left(17\right)$

即： $\min J={\min }_{{\tilde{U}}_f\left(k\right)}\frac{1}{2}\mathrm{\Δ}{\tilde{U}}_f{\left(k\right)}^T{Q}^{\′}\mathrm{\Δ}{\tilde{U}}_f\left(k\right)+R^{\′}\mathrm{\Δ}{\tilde{U}}_f\left(k\right)$

$\begin{array}{cc}s.t.& \mathrm{\Σ\Δ}{\tilde{U}}_f\end{array}\left(k\right)\≤T---\left(18\right)$

这里的Q是输出的权重矩阵，R是输入的权重矩阵；Y_tgt是预测信号的目标值，Q′、R′是由式17推到所得的等价权重矩阵。线性约束条件是指输入/输出的物理约束，以及其线性变化速率。

相较于现有技术，本发明的技术方案具备以下有益效果：

1.本发明提供的一种层燃锅炉的多变量控制与实时优化的方法，将一个3变量输入4变量输出的控制问题拆分为两个前馈反馈相互耦合的控制问题，使得控制过程更加迅速、稳定、抗干扰能力强。

2.本发明提供的一种层燃锅炉的多变量控制与实时优化的方法，将控制变量和扰动变量的增量作为输入信号，并结合状态变量进行控制，大大减少了MPC预测与优化的计算量，使控制更加精确快捷。

3.本发明提供的一种层燃锅炉的多变量控制与实时优化的方法，将过去的控制输入和扰动输入整合到对未来的预测过程中，使得控制器在作用过程中不局限于当前时刻开始预测，而是将延时部分一并考虑进预测过程中来。

4.本发明提供的一种层燃锅炉的多变量控制与实时优化的方法，将未来扰动信号输入部分设为0，简化了预测过程。

5.本发明提供了一种线性控制器解决非线性问题的有效策略，由于层燃锅炉的强非线性特点，传统的线性控制器在实际应用中控制效果一般都不理想，这是由于锅炉在不同工作点，具有不同的动静态系统特征。为了很好的反应系统非线性，该发明在基本控制算法不变的情况下，增加了控制组切换机制。根据锅炉的工作速度不同，利用不同的多变量模型进行系统动静态表达，而且实现了不同模型组切换的自适应，全智能，有效解决了层燃锅炉非线性控制难题。

具体实施方式

下文结合具体实施方式对本发明做进一步说明：

一种层燃锅炉的多变量控制与实时优化的方法，包括如下步骤：

1)将层燃锅炉的控制转化为3变量输入和4变量输出的多变量控制问题；所述3变量输入包括引风转速u₁、鼓风转速u₂、炉排速度u₃；所述4变量输出包括炉膛温度y₁、炉膛负压y₂、含氧量y₃、火焰熄灭位置y₄；

具体表达为： 

$\left[\begin{array}{c}{y}_1\\ y_2\\ y_3\\ y_4\end{array}\right]=\left[\begin{array}{ccc}{g}_{11}& g_{12}& g_{13}\\ g_{21}& g_{22}& g_{23}\\ g_{31}& g_{32}& g_{33}\\ g_{41}& g_{42}& g_{43}\end{array}\right]\left[\begin{array}{c}{u}_1\\ u_2\\ u_3\end{array}\right]---\left(1\right)$

G_ij表示第j个操控变量与第i个被控变量间的传递函数模型；由于炉排速度u3 对炉膛负压y2、含氧量y3的影响小，因此G₂₃、G₃₃这里取值为0，即：

$\left[\begin{array}{c}{y}_1\\ y_2\\ y_3\\ y_4\end{array}\right]=\left[\begin{array}{ccc}{g}_{11}& g_{12}& g_{13}\\ g_{21}& g_{22}& 0\\ g_{31}& g_{32}& 0\\ g_{41}& g_{42}& g_{43}\end{array}\right]\left[\begin{array}{c}{u}_1\\ u_2\\ u_3\end{array}\right]---\left(2\right)$

2)将所述多变量控制问题拆分为两个带前馈反馈功能的多变量控制器相互耦合的控制问题；其中，第一控制器的反馈控制变量在第二控制器中以前馈预测变量的形式出现；第二控制器的反馈控制变量在第一控制器中以前馈预测变量的形式出现；

具体表达为： 

$\left[\begin{array}{c}{y}_1\\ y_2\\ y_3\end{array}\right]=\left[\begin{array}{ccc}{g}_{11}& g_{12}& g_{13}\\ g_{21}& g_{22}& 0\\ g_{31}& g_{32}& g_0\end{array}\right]\left[\begin{array}{c}{u}_1\\ u_2\\ u_3\end{array}\right]=\left[\begin{array}{cc}{g}_{11}& g_{12}\\ g_{21}& g_{22}\\ g_{31}& g_{32}\end{array}\right]\left[\begin{array}{c}{u}_1\\ u_2\end{array}\right]+\left[\begin{array}{c}{g}_{13}\\ 0\\ 0\end{array}\right]u_3---\left(3\right)$

即：Y₁＝G₁U₁+G_d1D₁    (4) 

$y_4=\left[\begin{array}{ccc}{g}_{41}& g_{42}& g_{43}\end{array}\right]\left[\begin{array}{c}{u}_1\\ u_2\\ u_3\end{array}\right]=g_{43}{u}_3+\begin{array}{}\end{array}\left[\begin{array}{cc}{g}_{41}& g\end{array},42\right]\left[\begin{array}{c}{u}_1\\ u_2\end{array}\right]---\left(5\right)$

即：Y₂＝G₂U₂+G_d2D₂   (6) 

其中，Y₁、Y₂为被控矢量，U₁、U₂为输入矢量，D₁、D₂为扰动矢量； 

比较上式可见，引风转速u₁、鼓风转速u₂与炉膛温度y₁、炉膛负压y₂、含氧量y₃形成闭环反馈控制回路，而炉排速度u₃与炉膛温度y₁形成的前馈预测控制回路；同理，炉排速度u₃与火焰熄灭位置y₄形成闭环反馈控制回路，而引风转速u₁、鼓风转速u₂与火焰熄灭位置y₄形成的前馈预测控制回路。

接下来，对所述两个带前馈反馈功能的多变量控制问题分别进行状态空间转换；

所述状态空间转换的具体过程为：

将式5和式6进行归纳，可以写出一类控制的通式，

Y＝GU+G_dD    (7)

其中，Y是m维被控矢量，G是m*n的控制矩阵，U是n维输入矢量，G_d是m*q的扰动矩阵，D是q维扰动矢量，这里的m>0，n>0，q≥0，且都为整数；

首先将G和G_d转换为空间状态表达式:

Y(k)＝GU(k-T_d1)+G_dD(k-T_d2)     (8) 

进一步写为：Y(k)＝Y′+Y″     (9) 

其中T_d1为输入最大时延，T_d2为扰动最大时延；Y′为操控响应，Y″为扰动响应；

对Y′、Y″分别进行状态空间转换：

首先，对Y′进行状态空间转换，可得：

$\left\{\begin{array}{c}x(k+1)=\mathrm{Ax}\left(k\right)+\mathrm{BU}(k-T_{d1})\\ Y^{\′}\left(k\right)=\mathrm{Cx}\left(k\right)\end{array}\right.---\left(10\right)$

其中x(k)为状态变量，u(k-T_d1)为输入变量，A为系统矩阵、B为输入矩阵、C为输出矩阵。本控制方案的特点在于，对式10进行了进一步状态升阶操作，状态空间表达式10可以转化为：

$\left\{\begin{array}{c}\stackrel{\‾}{X}(k+1)=\stackrel{\‾}{A}\stackrel{\‾}{X}\left(k\right)+\stackrel{\‾}{B}\mathrm{\ΔU}(k-T_{d1})\\ {\stackrel{\‾}{Y}}^{\′}\left(k\right)=\stackrel{\‾}{C}\stackrel{\‾}{X}\left(k\right)\end{array}\right.---\left(11\right)$

其中状态变量x(k)被升阶为： 

$\stackrel{\‾}{X}\left(k\right)=\left[\begin{array}{c}\mathrm{\Δx}\left(k\right)\\ x\left(k\right)\end{array}\right],\mathrm{\Δx}\left(k\right)=x\left(k\right)-x(k-1).$

升阶状态矩阵由原状态空间表达式10中的系数矩阵(A,B,C)构成；

同理，对Y″进行状态空间转换，可得：

$\left\{\begin{array}{c}{x}_d(k+1)=A_d{x}_d\left(k\right)+B_{d}D(k-T_{d2})\\ Y^{\′\′}=C_d{x}_d\left(k\right)\end{array}\right.---\left(12\right)$

状态空间表达式12可被升阶为： 

$\left\{\begin{array}{c}{\stackrel{\‾}{X}}_d(k+1)={\stackrel{\‾}{A}}_d{\stackrel{\‾}{X}}_d\left(k\right)+{\stackrel{\‾}{B}}_d\mathrm{\ΔD}(k-T_{d2})\\ {\stackrel{\‾}{Y}}^{\′\′}\left(k\right)={\stackrel{\‾}{C}}_d{\stackrel{\‾}{X}}_d\left(k\right)\end{array}\right.---\left(13\right)$

其中状态x(k)被升阶为： 

${\stackrel{\‾}{X}}_d=\left[\begin{array}{c}\mathrm{\Δ}{x}_d\left(k\right)\\ x_d\left(k\right)\end{array}\right],\mathrm{\Δ}{x}_d\left(k\right)=x_d\left(k\right)-x_d(k-1).$

升阶状态矩阵由原状态空间表达式12中的系数矩阵(A,B,C)构成；是一个满足因果关系的扰动,在预测周期内该扰动保持不变。

3)对经过状态空间转换后的两个带前馈反馈功能的多变量控制问题基于当前时刻拆分成过去和将来作用两部分，并转化为求带有约束条件的二次规划最优解问题。

具体过程为： 

设定预测过程中的当前时刻为k，给定输入的控制步长为H_u，即输入序列的时间跨度为：ΔU_k→ΔU_k+Hu-1；输出的预测步长H_p，即从预测当前时刻的下一步开始，输出轨迹的时间跨度为Y_k+1→Y_k+Hp；

将式11的操控输入信号和式13的扰动输入信号基于当前时刻k拆分成过去和将来作用两部分，转换后的状态空间表达式如下所示：

$\left\{\begin{array}{c}\tilde{X}\left(k\right)=\tilde{A}\stackrel{\‾}{X}\left(k\right)+{\tilde{B}}_p\mathrm{\Δ}{\tilde{U}}_p\left(k\right)+{\tilde{B}}_f\mathrm{\Δ}{\tilde{U}}_f\left(k\right)\\ {\tilde{X}}_d\left(k\right)={\tilde{A}}_d{\stackrel{\‾}{X}}_d\left(k\right)+{\tilde{B}}_{\mathrm{dp}}\mathrm{\Δ}{\tilde{D}}_p\left(k\right)+{\tilde{B}}_{\mathrm{df}}\mathrm{\Δ}{\tilde{D}}_f\left(k\right)\end{array}\right.---\left(14\right)$

和分别指过去的操控信号输入和扰动信号输入的作用部分，和分别指未来预测的操控信号输入和扰动信号输入的作用部分；

其中，操控信号输入的各个系数详细表达为：

$\tilde{X}\left(k\right)=\left[\begin{array}{c}\stackrel{\‾}{X}(k+1)\\ \stackrel{\‾}{X}(k+2)\\ .\\ .\\ .\\ \stackrel{\‾}{X}(k+\mathrm{Hp})\end{array}\right]$

$\stackrel{\‾}{A}=\left[\begin{array}{c}\stackrel{\‾}{A}\\ {\stackrel{\‾}{A}}^2\\ .\\ .\\ .\\ {\stackrel{\‾}{A}}^{\mathrm{Hp}}\end{array}\right]$

${\tilde{B}}_p=\left[\begin{array}{cccccc}\stackrel{\‾}{B}& 0& ...& ...& ...& 0\\ \stackrel{\‾}{A}\stackrel{\‾}{B}& \stackrel{\‾}{B}& 0& ...& ...& 0\\ {\stackrel{\‾}{A}}^2& \stackrel{\‾}{A}\stackrel{\‾}{B}& \stackrel{\‾}{B}& 0& ...& 0\\ .& .& .& .& .& .\\ .& .& .& .& .& .\\ .& .& .& .& .& .\\ {\stackrel{\‾}{A}}^{\mathrm{Td}-1}\stackrel{\‾}{B}& {\stackrel{\‾}{A}}^{\mathrm{Td}-2}\stackrel{\‾}{B}& ...& ...& ...& \stackrel{\‾}{B}\\ {\stackrel{\‾}{A}}^{\mathrm{Td}}\stackrel{\‾}{B}& {\stackrel{\‾}{A}}^{\mathrm{Td}-1}\stackrel{\‾}{B}& ...& ...& ...& \stackrel{\‾}{A}\stackrel{\‾}{B}\\ .& .& .& .& .& .\\ .& .& .& .& .& .\\ .& .& .& .& .& .\\ {\stackrel{\‾}{A}}^{\mathrm{Hp}-1}\stackrel{\‾}{B}& {\stackrel{\‾}{A}}^{\mathrm{Hp}-2}\stackrel{\‾}{B}& ...& ...& ...& {\stackrel{\‾}{A}}^{\mathrm{Hp}-\mathrm{Td}-2}\stackrel{\‾}{B}\end{array}\right]$

$\mathrm{\Δ}{\tilde{U}}_p\left(k\right)=\left[\begin{array}{c}\mathrm{\ΔU}(k-T_d)\\ \mathrm{\ΔU}(k-T_d+1)\\ .\\ .\\ .\\ \mathrm{\ΔU}(k-1)\end{array}\right]$

$\mathrm{\Δ}{\tilde{U}}_f\left(k\right)=\left[\begin{array}{c}\mathrm{\ΔU}\left(k\right)\\ \mathrm{\ΔU}(k+1)\\ .\\ .\\ .\\ \mathrm{\ΔU}(k+\mathrm{Hu}-1)\end{array}\right]$

扰动信号输入的各个系数详细表达为：

${\tilde{X}}_d\left(k\right)=\left[\begin{array}{c}{\stackrel{\‾}{X}}_d(k+1)\\ {\stackrel{\‾}{X}}_d(k+2)\\ .\\ .\\ .\\ {\stackrel{\‾}{X}}_d(k+\mathrm{Hp})\end{array}\right]$

$\tilde{A}=\left[\begin{array}{c}{\stackrel{\‾}{A}}_d\\ {{\stackrel{\‾}{A}}_d}^2\\ .\\ .\\ .\\ {{\stackrel{\‾}{A}}_d}^{\mathrm{Hp}}\end{array}\right]$

$\mathrm{\Δ}{\tilde{D}}_p\left(k\right)=\left[\begin{array}{c}\mathrm{\ΔD}(k-T_d)\\ \mathrm{\ΔD}(k-T_d+1)\\ .\\ .\\ .\\ \mathrm{\ΔD}(k-1)\end{array}\right]$

${\tilde{B}}_{\mathrm{dp}}=\left[\begin{array}{cccccc}{\stackrel{\‾}{B}}_d& 0& ...& ...& ...& 0\\ {\stackrel{\‾}{A}}_d{\stackrel{\‾}{B}}_d& {\stackrel{\‾}{B}}_d& 0& ...& ...& 0\\ {{\stackrel{\‾}{A}}_d}^2{\stackrel{\‾}{B}}_d& {\stackrel{\‾}{A}}_d{\stackrel{\‾}{B}}_d& {\stackrel{\‾}{B}}_d& 0& ...& 0\\ .& .& .& .& .& .\\ .& .& .& .& .& .\\ .& .& .& .& .& .\\ {{\stackrel{\‾}{A}}_d}^{\mathrm{Td}-1}{\stackrel{\‾}{B}}_d& {{\stackrel{\‾}{A}}_d}^{\mathrm{Td}-2}{\stackrel{\‾}{B}}_d& ...& ...& ...& {\stackrel{\‾}{B}}_d\\ {{\stackrel{\‾}{A}}_d}^{\mathrm{Td}}{\stackrel{\‾}{B}}_d& {{\stackrel{\‾}{A}}_d}^{\mathrm{Td}-1}{\stackrel{\‾}{B}}_d& ...& ...& ...& {\stackrel{\‾}{A}}_d{\stackrel{\‾}{B}}_d\\ .& .& .& .& .& .\\ .& .& .& .& .& .\\ .& .& .& .& .& .\\ {{\stackrel{\‾}{A}}_d}^{\mathrm{Hp}-1}{\stackrel{\‾}{B}}_d& {{\stackrel{\‾}{A}}_d}^{\mathrm{Hp}-2}{\stackrel{\‾}{B}}_d& ...& ...& ...& {{\stackrel{\‾}{A}}_d}^{\mathrm{Hp}-\mathrm{Td}-2}{\stackrel{\‾}{B}}_d\end{array}\right]$

$\mathrm{\Δ}{\tilde{D}}_f\left(k\right)=\left[\begin{array}{c}\mathrm{\ΔD}\left(k\right)\\ \mathrm{\ΔD}(k+1)\\ .\\ .\\ .\\ \mathrm{\ΔD}(k+\mathrm{Hu}-1)\end{array}\right]$

将式14中的设为0，则式11、13可以简化为：

$\left\{\begin{array}{c}{\tilde{Y}}^{\′}\left(k\right)=\tilde{C}\tilde{A}\stackrel{\‾}{X}\left(k\right)+\tilde{C}{\tilde{B}}_p\mathrm{\Δ}{\tilde{U}}_p\left(k\right)+\tilde{C}{\tilde{B}}_f\mathrm{\Δ}{\tilde{U}}_f\left(k\right)\\ {\tilde{Y}}^{\′\′}\left(k\right)={\tilde{C}}_d{\tilde{A}}_d{\stackrel{\‾}{X}}_d\left(k\right)+{\tilde{C}}_d{\tilde{B}}_{\mathrm{dp}}\mathrm{\Δ}{\tilde{D}}_p\left(k\right)\end{array}\right.---\left(15\right)$

其中，分别是由构成的输出系数矩阵，可分别表示成：  $\tilde{C}=\mathrm{diag}\{\stackrel{\‾}{C},\stackrel{\‾}{C},...,\stackrel{\‾}{C}\}$ 和 ${\tilde{C}}_d=\mathrm{diag}\{{\stackrel{\‾}{C}}_d,{\stackrel{\‾}{C}}_d,...,{\stackrel{\‾}{C}}_d\}.$

通过式7和式15，可以总结得出：

$\begin{array}{c}\tilde{Y}\left(k\right)={\tilde{Y}}^{\′}\left(k\right)+{\tilde{Y}}^{\′\′}\left(k\right)\\ =\tilde{C}{\tilde{B}}_f\mathrm{\Δ}{\tilde{U}}_f\left(k\right)+[\tilde{C}\tilde{A}\stackrel{\‾}{X}\left(k\right)+\tilde{C}{\tilde{B}}_p\mathrm{\Δ}{\tilde{U}}_p\left(k\right)+{\tilde{C}}_d{\tilde{A}}_d{\stackrel{\‾}{X}}_d\left(k\right)+{\tilde{C}}_d{\tilde{B}}_{\mathrm{dp}}\mathrm{\Δ}{\tilde{D}}_p\left(k\right)]\end{array}---\left(16\right)$

式16由未知的未来输入响应项和已知的过去输入响应项组成，对这两部分进行整理，可将优化问题转化成式17所示的求带有约束条件的二次规划最优解问题：

$J={(\tilde{Y}\left(k\right)-Y_{\mathrm{tgt}})}^{T}Q(\tilde{Y}\left(k\right)-Y_{\mathrm{tgt}})+\mathrm{\Δ}{\tilde{U}}_f{\left(k\right)}^T\mathrm{R\Δ}{\tilde{U}}_f\left(k\right)$

$\begin{array}{cc}s.t.& \mathrm{\Σ\Δ}{\tilde{U}}_f\end{array}\left(k\right)\≤T---\left(17\right)$

即： $\min J={\min }_{{\tilde{U}}_{f\left(k\right)}}\frac{1}{2}\mathrm{\Δ}{\tilde{U}}_f{\left(k\right)}^T{Q}^{\′}\mathrm{\Δ}{\tilde{U}}_f\left(k\right)+R^{\′}\mathrm{\Δ}{\tilde{U}}_f\left(k\right)$

$\begin{array}{cc}s.t.& \mathrm{\Σ\Δ}{\tilde{U}}_f\end{array}\left(k\right)\≤T---\left(18\right)$

这里的Q是输出的权重矩阵，R是输入的权重矩阵；Y_tgt是预测信号的目标值，Q′、R′是由式17推到所得的等价权重矩阵。线性约束条件是指输入/输出的物理约束，以及其线性变化速率。

4)通过锅炉排速度区间，划分不同的层燃锅炉模型组，通过模型组自动切换解决层燃锅炉的强非线性问题。

以上所述，仅为本发明专利较佳的具体实施方式，但本发明专利的保护范围并不局限于此，任何熟悉本技术领域的技术人员在本发明揭露的技术范围内，可轻易想到的变化或替换，都应涵盖在本发明的保护范围之内。因此，本发明专利的保护范围应该以权利要求的保护范围为准。

Claims (5)

1.一种层燃锅炉的多变量控制与实时优化的方法，其特征在于包括如下步骤：

1)将层燃锅炉的控制转化为3变量输入和4变量输出的多变量控制问题；所述3变量输入包括引风转速u₁、鼓风转速u₂、炉排速度u₃；所述4变量输出包括炉膛温度y₁、炉膛负压y₂、含氧量y₃、火焰熄灭位置y₄；

2)将所述多变量控制问题拆分为两个带前馈反馈功能的多变量控制器相互耦合的控制问题；其中，第一控制器的反馈控制变量在第二控制器中以前馈预测变量的形式出现；第二控制器的反馈控制变量在第一控制器中以前馈预测变量的形式出现；对所述两个带前反馈功能的多变量控制问题分别进行状态空间升阶转换，将当前时刻输入信号的增量作为新的输入信号进行控制；

3)对经过状态空间转换后的两个带前馈反馈功能的多变量控制问题基于当前时刻拆分成过去和将来作用两部分，并转化为求带有约束条件的二次规划最优解问题；

4)通过锅炉排速度区间，划分不同的层燃锅炉模型组，通过模型组自动切换解决层燃锅炉的强非线性问题。

2.根据权利要求1所述的一种层燃锅炉的多变量控制与实时优化的方法，其特征在于：所述3变量输入和4变量输出的多变量控制问题具体表达为：

g_ij表示第j个操控变量与第i个被控变量间的传递函数模型；由于 炉排速度u3对炉膛负压y2、含氧量y3的影响小，因此g₂₃、g₃₃这里取值为0，即：

3.根据权利要求2所述的一种层燃锅炉的多变量控制与实时优化的方法，其特征在于：所述两个带前馈反馈功能的多变量控制器控制问题具体表达为：

即：Y₁＝G₁U₁+G_d1D₁              (4) 

即：Y₂＝G₂U₂+G_d2D₂         (6) 

其中，Y₁、Y₂为被控矢量，U₁、U₂为输入矢量，D₁、D₂为扰动矢量；

比较上式可见，引风转速u₁、鼓风转速u₂与炉膛温度y₁、炉膛负压y₂、含氧量y₃形成闭环反馈控制回路，而炉排速度u₃与炉膛温度y₁形成的前馈预测控制回路；同理，炉排速度u₃与火焰熄灭位置y₄形成闭环反馈控制回路，而引风转速u₁、鼓风转速u₂与火焰熄灭位置y₄形成的前馈预测控制回路。

将式4和式6进行归纳，可以写出一类控制的通式，

Y＝GU+G_dD             (7) 

其中，Y是m维被控矢量，G是m*n的控制矩阵，U是n维输入矢量，G_d是m*q的扰动矩阵，D是q维扰动矢量，这里的m>0，n>0，q≥0， 且都为整数；

首先将G和G_d转换为空间状态表达式:

Y(k)＝GU(k-T_d1)+G_dD(k-T_d2)          (8) 

进一步写为：Y(k)＝Y′+Y″                (9) 

其中T_d1为输入最大时延，T_d2为扰动最大时延；Y′为操控响应，Y″为扰动响应；

对Y′、Y″分别进行状态空间转换：

首先，对Y′进行状态空间转换，可得：

其中x(k)为状态变量，U(k-T_d1)为输入变量，A为系统矩阵、B为输入矩阵、C为输出矩阵。本控制方案的特点在于，对式10进行了进一步状态升阶操作，状态空间表达式10可以转化为：

其中状态变量x(k)被升阶为：

升阶状态矩阵由原状态空间表达式10中的系数矩阵(A,B,C)构成；

同理，对Y″进行状态空间转换，可得：

状态空间表达式12可被升阶为：

其中状态x(k)被升阶为：

升阶状态矩阵由原状态空间表达式12中的系数矩阵(A,B,C)构成；是一个满足因果关系的扰动,在预测周期内该扰动保持不变；

这里将增量ΔU(k)＝U(k)-U(k-1)作为系统的输入信号，包括控制变量与扰动变量，并结合状态变量进行控制。

4.根据权利要求3所述的一种层燃锅炉的多变量控制与实时优化的方法，其特征在于：所述步骤3的具体过程为:

设定预测过程中的当前时刻为k，给定输入的控制步长为H_u，即输入序列的时间跨度为：ΔU_k→ΔU_k+Hu-1；输出的预测步长H_p，即从预测当前时刻的下一步开始，输出轨迹的时间跨度为Y_k+1→Y_k+Hp；

将式11的操控输入信号和式13的扰动输入信号基于当前时刻k拆分成过去和将来作用两部分，转换后的状态空间表达式如下所示：

和分别指过去的操控信号输入和扰动信号输入的作用部分，和分别指未来预测的操控信号输入和扰动信号输入的作用部分；

其中，操控信号输入的各个系数详细表达为：

扰动信号输入的各个系数详细表达为：

5.根据权利要求4所述的一种层燃锅炉的多变量控制与实时优化的方法，其特征在于：将式14中的设为0，则式11、13可以简化为：

其中，分别是由构成的输出系数矩阵，可分别表示成： 和

通过式7和式15，可以总结得出：

式16由未知的未来输入响应项和已知的过去输入响应项组成，对这两部分进行整理，可将优化问题转化成式17所示的求带有约束条件的二次规划最优解问题：

即：

这里的Q是输出的权重矩阵，R是输入的权重矩阵；Y_tht是预测信号的目标值，Q′、R′是由式17推到所得的等价权重矩阵。线性约束条件是指输入/输出的物理约束，以及其线性变化速率。