CN104765350A - 基于组合模型预测控制技术的水泥分解炉控制方法及系统 - Google Patents

Abstract

本发明涉及基于组合模型预测控制技术的水泥分解炉控制方法及系统，其方法为：(1)中控机利用数据通讯接口采集水泥分解炉生料预分解过程数据；(2)对数据进行分类，分别进行模型辨识，得到LSSVM稳态模型和ARMAX动态模型，并将其结合为组合模型；(3)用含有时滞环节的组合模型预测出生料预分解过程的未来输出状态，并设定输出值的参考轨迹；(4)非线性控制器采用序列二次规划法对目标函数滚动求解，得到水泥分解炉的最优控制值。其系统包括分解炉智能测量仪表、执行器、通讯接口及中控机，组合模型预测控制算法嵌入在中控机中。本发明模型辨识准确度高，算法鲁棒性强，可适应分解炉多变量间耦合性、非线性、时滞性，具有良好的控制效果。

Description

基于组合模型预测控制技术的水泥分解炉控制方法及系统

技术领域

本发明涉及水泥生料预分解过程的先进控制领域，特别是涉及一种基于组合模型预测控制技术的水泥分解炉控制方法和系统。

背景技术

水泥分解炉是新型干法水泥生产过程中的关键设备，生料中大部分碳酸盐在其中进行分解，生料稳定有效的分解可保证整条水泥生产线的长期稳定运行。随着国家对经济结构调整和对节能减排等方面的要求，使得水泥分解炉智能控制的研究显得越来越重要。

目前，我国水泥生产针对分解炉控制主要依靠操作员凭其经验来操作控制。由于操作员的操作经验参差不齐，极易导致分解炉温度忽高忽低、分解炉系统气体含量波动较大以及生料预分解率不稳定，造成整个水泥烧成系统工况波动较大以及能源浪费。还由于分解炉控制具有明显非线性、大滞后和大扰动的特性并且难以建立其数理模型，导致传统控制方法难以适用。

众多国内外工艺、自控专家对此做了大量的研究工作，提出了许多先进控制理论，也取得了一些令人鼓舞的进展，但是却难以有效地建立分解炉系统的模型，并且相关算法计算量大、应用性差。因此，有必要从水泥分解炉的自身特性出发寻求一种精确度高的建模方法及控制算法，实现水泥分解炉系统的智能控制。

发明内容

针对现有分解炉控制方案的控制变量单一和不能适应多变量间强耦合性、非线性特性的不足，本发明要解决的技术问题是提供一种能够克服非线性、强耦合和大干扰的影响，并能解决分解炉多变量控制的基于组合模型预测控制技术的水泥分解炉控制方法和系统。

为解决上述技术问题，本发明的目的之一是通过以下方案实现的：

一种基于组合模型预测控制技术的水泥分解炉控制方法，其内容包括以下步骤：

步骤一：控制变量选取，分解过程的主要操作变量是分解炉喂煤量、高温风机转速和三次风管挡板开度，被控变量是分解炉出口温度、一级筒氧气含量和一级筒一氧化碳含量；

步骤二：数据采集及分类，通过中控机的数据通讯接口采集步骤一所述各操作变量和被控变量的现场数据，并按照阶跃响应特性将数据分类为稳态数据和动态数据；

步骤三：LSSVM-ARMAX组合模型辨识，利用基于粒子群优化参数的最小二乘支持向量机对分类出的稳态数据进行辨识得到水泥分解炉系统的非线性LSSVM稳态模型，利用递推增广最小二乘法对动态数据辨识得到水泥分解炉系统的线性ARMAX动态模型，并将稳态模型与动态模型有机并联组合，得到水泥分解炉系统的LSSVM-ARMAX组合模型；

步骤四：输出状态预测，在当前时刻k，利用水泥分解炉的输入输出信息和预测的未来输入信息，通过步骤三辨识得到的LSSVM-ARMAX组合模型，将其进行迭代预测出分解炉出口温度y_1n(k+j)、一级筒氧气含量y_2n(k+j)及一级筒一氧化碳含量y_3n(k+j)未来预测时域P的输出状态，式中j＝1,2,···,P；

步骤五：时滞及误差校正，根据步骤四中得到的各变量未来预测时域输出状态记为y_in(k+j)，进行时滞处理，得到各变量输出时滞预测状态y_in(k|k-1)，将其与k+1时刻的实际输出值相减求得输出误差e_i(k)，然后将y_in(k+j)与误差相加得到基于k时刻水泥分解炉未来预测输出y_ic(k+j)，上述各式中i＝1,2,3；

步骤六：滚动优化计算，引入参考输入轨迹分解炉出口温度y_1r(k+j)、一级筒氧气含量y_2r(k+j)以及一级筒一氧化碳含量y_3r(k+j)，将其与步骤五求得分解炉各未来预测输出y_ic(k+j)进行比较，构建带约束的二次型目标函数，应用序列二次规划法对其进行滚动优化求解，计算出当前时刻水泥分解炉被控变量分解炉喂煤量u₁(k)、高温风机转速u₂(k)和三次风管挡板开度u₃(k)；

步骤七：实时输出控制，根据步骤六计算出的当前时刻水泥分解炉被控变量控制现场执行器，通过滚动实时控制，使现场水泥分解炉被控变量实际输出与设定值y_sp相等，实现水泥分解炉的自动控制。

在步骤一中，所述的操作变量：分解炉喂煤量、高温风机转速和三次风管挡板开度是影响水泥分解炉系统的主要输入控制变量，它们是分解炉内热量和气体的主要来源；所述的被控变量：分解炉出口温度、一级筒氧气含量和一级筒一氧化碳含量是影响水泥分解炉系统的主要输出控制变量，选取分解炉出口温度作为输出变量即能直接反应生料预分解所需热量，也能间接反应水泥生料预分解情况；分解炉出口温度的稳定可保证生料预分解过程的稳定有效运行；另选取一级筒氧气含量和一级筒一氧化碳含量既能反应分解炉内燃煤的燃烧状况，也兼顾了分解炉单产能耗和安全运转。

在步骤二中，所述的数据采集是通过中控机的数据通讯接口采集步骤一所述各控制变量和被控变量的现场数据并进行存储；所述的数据分类是依据现场数据的阶跃响应特性将其分类为稳态数据和动态数据。

在步骤三中，所述的水泥分解炉系统的LSSVM-ARMAX组合模型是通过非线性稳态模型的增益在线调节线性动态模型参数实现模型的有机并联组合，水泥分解炉的非线性LSSVM稳态模型是通过粒子群优化参数的最小二乘支持向量机根据分解炉的稳态数据辨识获得，水泥分解炉系统的线性ARMAX动态模型是通过递推增广最小二乘法对分解炉的动态数据辨识获得。

在步骤六中所述的滚动优化计算是在滚动时域下采用序列二次规划法对带有约束的二次型目标函数进行求解，保证求出的分解炉预测输入解是在约束条件下的全局最优解。

本发明的另一目的是提供一种基于组合模型预测控制技术的水泥分解炉控制系统，该控制系统包括测量仪表、数据通讯接口、中控机和执行器；

所述的测量仪表用于测量水泥分解炉系统的输出量，即测量分解炉出口温度、一级筒气体含量和一级筒一氧化碳含量；

所述的数据通讯接口用于将现场测量仪表测量的数据传输到中控机，并将中控机处理后下达的输入量变化增量传输到现场执行器；

所述的中控机用于运行水泥分解炉组合模型预测控制算法，根据水泥分解炉系统的实际输出量计算出当前时刻应加入分解炉内燃煤用量、高温风机转速以及三次风管挡板开度的控制量；

所述的执行器通过数据通讯接口接收中控机下达的变量调节量指令，通过其自动调节，实现现场分解炉的优化自动控制。

本发明具有以下有益效果：

1、选择分解炉出口温度、一级筒氧气含量和一级筒一氧化碳含量作为水泥分解炉预测控制器的三个输出量，选择分解炉喂煤量、高温风机转速和三次风管挡板开度作为其三个输入量，这既完全准确的反应了分解炉实际运行状况，又综合考虑了分解炉的单位能耗及安全运行的问题；

2、无须使用经验知识，只需使用输入输出分类的数据分别进行非线性LSSVM稳态模型和线性ARMAX动态模型的离线辨识，再通过动态增益即可将其有机结合为反映水泥分解炉系统的组合模型，辨识过程简单，模型精度高；

3、充分利用预测控制技术的优点，引入在线实时预测输出值、参考轨迹和滚动优化计算技术，设计的水泥分解炉组合模型预测控制算法克服了时滞因素影响，平稳达到各变量设定值，跟踪性能好，抗干扰能力强；

4、利用序列二次规划法对带有约束的二次型目标函数进行滚动优化求解，即易得到系统全局最优解，又充分考虑了各变量的约束，防止了系统操作变量发生跳变。

附图说明

图1为基于组合模型预测控制技术的水泥分解炉控制方法流程图；

图2为本发明的LSSVM-ARMAX组合模型框图；

图3为本发明提出的基于组合模型预测控制技术的水泥分解炉控制系统的方框图；

图4为本发明提出的基于组合模型预测控制技术的水泥分解炉控制系统的现场接线图。

具体实施方式

下面结合附图对本发明做进一步的详细描述。

从图1所示的基于组合模型预测控制技术的水泥分解炉控制方法的流程图和图3所示的基于组合模型预测控制技术的水泥分解炉控制系统方框图可知，本发明的一种基于组合模型预测控制技术的水泥分解炉控制方法，其内容包括如下步骤：

步骤一控制变量选取

分解炉喂煤量、高温风机转速和三次风管挡板开度是影响水泥分解炉系统的主要输入控制变量，它们是分解炉内热量和气体的主要来源；分解炉出口温度、一级筒氧气含量和一级筒一氧化碳含量是影响水泥分解炉系统的主要输出控制变量，选取分解炉出口温度作为输出变量即可直接反应生料预分解所需热量，也可间接反应水泥生料预分解情况；分解炉出口温度的稳定可保证生料预分解过程的稳定有效运行；另选取一级筒氧气含量和一级筒一氧化碳含量即可反应分解炉内燃煤的燃烧状况，也兼顾了分解炉单产能耗和安全运转。

因此，本发明具体实施方式中选取分解炉喂煤量、高温风机转速和三次风管挡板开度作为水泥分解炉分解过程的主要操作变量，选取分解炉出口温度、一级筒氧气含量和一级筒一氧化碳含量作为主要被控变量。

步骤二数据采集及分类

本发明具体实施方式中，通过图3所示的数据通讯接口10采集某水泥生产线的分解炉生料预分解过程中步骤一所述六个变量的现场测量仪表数据，采样周期为30s，共采集30000组数据。

将所采集数据分成动态数据和稳态数据两种数据，其划分依据为：当输入量分解炉喂煤量u₁、高温风机转速u₂和三次风管挡板开度u₃中任意一个发生阶跃变化时，输出量分解炉出口温度y₁、一级筒O₂含量y₂和一级筒CO含量y₃从阶跃时刻到最终都达到稳定状态时的数据选为动态数据；由于辨识算法的特殊性，只选取输出量从稳定状态到下一个阶跃发生前的少部分数据为稳态数据。经数据分类及不合理数据删除后得到动态数据9823组，稳态数据1632组。

步骤三LSSVM-ARMAX组合模型辨识

LSSVM-ARMAX组合模型是一种由非线性LSSVM稳态模型和线性ARMAX动态模型并联而成的组合模型，其详细的模型结构关系如图2所示。其模型结构分为三部分：LSSVM稳态模型6、动态增益K7和ARMAX动态模型8，用LSSVM稳态模型求出当前时刻系统的动态增益K，再利用动态增益K实时调整ARMAX动态模型参数，完成稳态模型和动态模型的有机并联组合，建立代表水泥分解炉系统综合特性的LSSVM-ARMAX组合模型。

非线性LSSVM稳态模型辨识

本发明具体实施方式中，基于粒子群优化参数的最小二乘支持向量机辨识算法采用高斯径向基核函数作为最小二乘支持向量机的核函数，采用粒子群优化算法对支持向量机正则化参数C与核宽度σ进行自动寻优，初始化粒子群中最大遗传代数maxgen＝20，种群规模sizepop＝30，学习因子c₁＝c₂＝2，算法中w为惯性权重，eps为适应度预设值，粒子群优化LSSVM参数算法(PSO)的具体步骤如下：

1)初始化粒子群；

2)给定正则化系数C与核宽度σ的取值区间：

C∈[C_max,C_min]，σ∈[σ_max,σ_min]  (1)

3)确定正则化系数C与核宽度σ迭代速度的最大值与最小值，k为常数；

V_Cmax＝k·C_max

V_Cmin＝-k·C_min  (2)

V_dmax＝k·σ_max

V_dmin＝-k·σ_min

4)随机产生种群的初始位置pop和初始速度V：

pop(1,i)＝(C_max-C_min)*rand()+C_min

pop(2,i)＝(deta_max-deta_min)*rand()+deta_min  (3)

V(1,i)＝Vcmax*rand()

V(2,i)＝Vdmax*rand()

5)LSSVM辨识，以测试集的预测输出与实际输出间的均方误差作为每个粒子的适应度fitness，全局极值点p_g和个体极值点p为：

p_g＝min(fitness)

p＝fitness  (4)

6)检查是否满足结束条件：p_g≤eps或已经达到最大迭代次数；如果满足，则跳转到步骤10)；

7)分别根据式(5)和式(6)更新粒子的速度与位置，如果超出范围则可利用程序将其限定在规定的范围内；

$v_{\mathrm{id}}^{(t+1)}=w*v_{\mathrm{id}}^{\left(t\right)}+c_1*r_1*(p_{\mathrm{id}}^{\left(t\right)}-x_{\mathrm{id}}^{\left(t\right)})+c_2*r_2*(p_{\mathrm{gd}}^{\left(t\right)}-x_{\mathrm{id}}^{\left(t\right)})---\left(5\right)$

$x_{\mathrm{id}}^{(t+1)}=x_{\mathrm{id}}^{\left(t\right)}+v_{\mathrm{id}}^{(t+1)}---\left(6\right)$

8)LSSVM辨识，以测试集的均方误差作为当前每个粒子的适应度值fitness，如果某个粒子的适应度值小于其个体极值，则对其个体极值进行更新；如果某个粒子的适应度值小于全局极值，则对全局极值进行更新，同时更新最优解正则化系数C与核宽度σ组成的向量global_x；

9)迭代次数加1，跳转到步骤6)；

10)得到最优解global_x；

11)采用已优化的正则化系数C与核宽度σ进行LSSVM辨识，优化算法结束。

本发明具体实施方式中采用的最小二乘支持向量机在一定程度上提高了模型辨识的学习速度并在模式识别中得到了广泛应用。其具体实施方式中，训练数据集{(x₁,y₁),(x₂,y₂),···,(x_N,y_N)}为现场采集的稳态数据集，其中，N为训练样本数N＝1632，x_i∈R^N为输入变量数据集，y_i∈R分别为输出变量数据集，i＝1,2,···,N。

具体实施方式过程如下：

首先确定非线性回归函数如式(7)所示：

其中，w和b为待确定的参数，为核函数。

LSSVM的优化问题被定义为：

$\underset{w,e}{\min }J(w,e)=\frac{1}{2}{\left|\right|w\left|\right|}^2+\frac{C}{2}\underset{i=1}{\overset{N}{\mathrm{\Σ}}}{e}_i^2,C>0---\left(8\right)$

其中，目标函数的第一项反映了模型的泛化能力，目标函数的第二项体现了模型的准确度，正则化参数C调节二者的平衡，e_i表示第i个数据实际输出与预测输出间的误差。

假定矩阵：

$\mathrm{\Φ}=\left[\begin{array}{cc}0& E^T\\ E& \mathrm{\Ω}+{\mathrm{\γ}}^{-1}I\end{array}\right]---\left(9\right)$

则参数α和b的解可由式(10)得到：

$\left[\begin{array}{c}b\\ \mathrm{\α}\end{array}\right]={\mathrm{\Φ}}^{-1}\left[\begin{array}{c}0\\ y\end{array}\right]---\left(10\right)$

最终得到LSSVM模型的表达式：

$y\left(x\right)=\underset{i=1}{\overset{N}{\mathrm{\Σ}}}{\mathrm{\α}}_{i}K(x,x_i)+b---\left(11\right)$

在本具体实施方式中，核函数选取为高斯径向基核函数(RBF)：

$K(x,x_i)=\exp (-\frac{{\left|\right|x-x_i\left|\right|}^2}{2{\mathrm{\σ}}^2})---\left(12\right)$

在本发明具体实施方式LSSVM模型中正则化参数C与核宽度σ由上述粒子群优化算法PSO进行自动寻优求解。

综上所述，可以利用上述的PSO-LSSVM算法对给定的输入、输出样本数据进行训练，最终训练出如下反映水泥分解炉系统输入到输出关系的非线性稳态模型：

$y\left(x\right)=\underset{i=1}{\overset{N}{\mathrm{\Σ}}}{\mathrm{\α}}_{i}K(x,x_i)+b---\left(13\right)$

式中：x是系统的稳态输入(分解炉喂煤量、高温风机转速三次风挡板开度)，y(x)是系统的稳态输出(分解炉出口温度、一级筒O2含量和一级筒CO含量)，K(·,·)为选取的核函数，α_i和b是LSSVM的模型参数。

线性ARMAX时滞动态模型辨识

本发明的线性动态模型辨识是利用递推增广最小二乘法的辨识方法辨识得到带外部输入的自回归滑动平均模型ARMAX(autoregressive moving averagewith exogenous inputs)。带有色噪声的多变量ARMAX模型最常用表达式如下式所示：

A(z)y(t)＝B(z)u(t)+D(z)v(t)  (14)

式中 $y\left(t\right)={[y_1\left(t\right),y_2\left(t\right),...y_{n_y}\left(t\right)]}^T;u\left(t\right)={[u_1\left(t\right),u_2\left(t\right),...u_{n_u}\left(t\right)]}^T;v\left(t\right)={[v_1\left(t\right),v_2\left(t\right),...v_{n_y}\left(t\right)]}^T;$ n_u，n_y分别是输入变量的个数和输出变量的个数，本发明中n_u＝3，n_y＝3。

为了方便推导多变量ARMAX模型辨识方法，定义以下参数向量：

$\begin{array}{c}\mathrm{\φ}\left(t\right)=[-y^T(t-1),{-y}^T(t-2),...-y^T(t-n_a),u^T(t-1),u^T(t-2),\\ {...u^T(t-n_b),v^T(t-1),v^T(t-2),...v^T(t-n_d)]}^T\∈R^{(n_y{n}_a+n_u{n}_b+n_y{n}_d)\×1}\end{array}---\left(15\right)$

则输出变量可以写成如下向量形式：

递推增广最小二乘法算法步骤如下：

1)采集输入输出数据，设数据共L组，本发明具体实施方式中现场采集的水泥分解炉动态输入输出数据共计9823组

2)设定初值，对于协方差阵初始值P(0)可以先选择一个大正定单位矩阵本发明实例中设定P(0)＝10⁶*I₀，参数估计初值可选取一个很小的实相量，如其中是单位均为1的n_y×(n_yn_a+n_un_b+n_yn_d)维向量。

3)按照式(18)计算协方差矩阵P(t)由式(19)求得信息向量φ(t)。

$P\left(t\right)=[I-L\left(t\right)\widehat{{\mathrm{\φ}}^T}\left(t\right)]P(t-1)\∈R^{(n_y{n}_a+n_u{n}_b+n_y{n}_d)\×(n_y{n}_a+n_u{n}_b+n_y{n}_d)}---\left(18\right)$

$\begin{array}{c}\widehat{\mathrm{\φ}}\left(t\right)=[-y(t-1),-y(t-2),...-y(t-n_a),u(t-1),u(t-2),...u(t-n_b),\\ \hat{v}(t-1),\hat{v}(t-2),...\hat{v}(t-n_d){]}^T\end{array}---\left(19\right)$

4)利用式(20)获得递推估计参数。

5)如果t＝L+1则辨识完成，得到最终辨识结果为否则t+1，重复3-5步。

通过上述算法可得到用白噪声信号代替噪声增量信号得到最终ARMAX模型：

Δy(k)＝A₁Δy(k-1)+A₂Δy(k-2)+B₁Δu(k-1)+B₂Δu(k-2)+Dv(k)  (21)

本发明辨识的ARMAX线性模型比ARX模型增加了噪声增量，带噪声ARMAX模型更能精确地表达被控非线性水泥分解炉系统的动态特性。

LSSVM-ARMAX组合模型建立

LSSVM-ARMAX组合模型是先将系统的稳态和动态特性进行单独辨识，最后再通过增益将LSSVM稳态模型与动态ARMAX模型进行有机结合，如图2所示，形成一个并联的组合模型。最终用LSSVM-ARMAX组合模型来完全体现出水泥分解炉系统的特性。稳态模型和动态模型具体结合过程如下，

ARMAX模型如式(21)所示，设：

A₀＝I-A₁-A₂  (22)

式中：I是n_y×n_y的单位矩阵。在增益问题上，可以先不考虑噪声的影响，式(21)中ARMAX增量模型的稳态增益矩阵可以表示为：

$K=\frac{B_1+B_2}{I-A_1-A_2}={(I-A_1-A_2)}^{-1}\·(B_1+B_2)=A_0^{-1}\·(B_1+B_2)---\left(23\right)$

ARMAX增量模型的稳态增益是描述实际水泥分解炉系统各变量之间变化的LSSVM稳态模型的动态增益。首先将动态ARMAX模型的增益设为1，具体的变换过程如下：

令

$B_i^{\′}=A_0\·\frac{B_i}{B_1+B_2}=(I-A_1-A_2)\·[B_i./(B_1+B_2)]---\left(24\right)$

式中：./代表两个矩阵对应位置的元素相除，i＝1,2,3。保持ARMAX模型(21)中的用A_i不变，用式(24)中的B'_i代替模型(21)中的B_i，忽略噪声项的变化，则可得到变换后的稳态增益为1的ARMAX模型如式(25)所示：

Δy(k)＝A₁Δy(k-1)+A₂Δy(k-2)+B'₁Δu(k-1)+B'₂Δu(k-2)+Dv(k)  (25)

下面利用数值微分法求LSSVM稳态模型的增益K，设当前时刻的输入为u_i(k)，u_i(k)发生微小变化Δu_i时的输入为u_i(k)+Δu_i，将输入量分别代入到LSSVM稳态模型(13)式中，可求得其对应的第j个输出分别为：

$y_{\mathrm{sj}}=f(u_1\left(k\right),u_2\left(k\right),...,u_i\left(k\right),...,u_{n_u}\left(k\right))---\left(26\right)$

$y_{\mathrm{sj}}^{\′}=f(u_1\left(k\right),u_2\left(k\right),...,u_i\left(k\right)+\mathrm{\Δ}{u}_i,...,u_{n_u}\left(k\right))---\left(27\right)$

式中：f(·)代表LSSVM稳态模型；Δu_i是一个很小的设定值，要根据输入量的值来选取；j＝1,2,···,n_y；i＝1,2,···,n_u；n_u、n_y分别表示系统输入变量个数和输出变量的个数，本实施方式中n_u＝n_y＝3。

根据数值微分法可以根据系统稳态输入输出可求得系统输入变量u_i(k)到输出变量y_j(k)的动态增益为：

$K_{i,j}^s=\frac{{\∂y}_j\left(k\right)}{{\∂u}_i\left(k\right)}=\frac{y_{\mathrm{sj}}^{\′}-y_{\mathrm{sj}}}{\mathrm{\Δ}{u}_i}---\left(28\right)$

通过式(27)可以求多变量系统任意输入点u(k)所对应的任意输出变量的LSSVM稳态模型的动态增益矩阵K^s如下式表示：

通过式(28)求得的实际系统的动态增益K^s在线修正增量ARMAX模型的动态增益，这样能保证ARMAX模型的动态增益和系统的实际增益保持一致，即可得到具体的调整方法是用K^s与B'_i相乘作为调整输入变量的参数：

Δy(k)＝A₁Δy(k-1)+A₂Δy(k-2)+B₁₁Δu(k-1)+B₁₂Δu(k-2)+

B₂₁Δu²(k-1)+B₂₂Δu²(k-2)+Dv(k)  (30)

式中：K'＝(Kⁿ-K^c).u_s/(k+(-u_s1k)，B₁₁＝A₀·[B'₁./(B'₁+B'₂).*K^c]，B₁₂＝A₀·[B'₂./(B'₁+B'₂).*K^c]，B₂₁＝A₀·[B'₁./(B'₁+B'₂).*K']，B₂₂＝A₀·[B'₂./(B'₁+B'₂).*K']，Kⁿ和K^c分别是u_s(k+1)和u_s(k)的稳态增益，B₁'和B₂'是ARMAX模型归一化的系数，./表示矩阵对应位置元素相除，.*表示两个矩阵中对应位置元素相乘。

步骤四预测未来输出状态

上述式(30)是基于k-1时刻以前的煅烧过程中输入输出量状态预测输出k时刻水泥分解炉的输出增量状态，得到预测时域P和控制时域M内的输出状态推导过程如下：

先预设如下三个递推表达式：

$\left(\begin{array}{c}S\left(n\right)=A_1\·S(n-1)+A_2\·S(n-2)\\ T\left(n\right)=S\left(n\right)\·B_{11}+S(n-1)\·B_{12}\\ T_g\left(n\right)=S\left(n\right)\·B_{21}+S(n-1)\·B_{22}\end{array}\right.---\left(31\right)$

式中：n＝3,4,···,P+1，S为n_y×n_y维矩阵，T、T_g为n_y×n_u维矩阵，亦即T、S(1)＝I，S(2)＝A₁，T(1)＝B₁₁，T_g(1)＝B₂₁，I为n_y×n_y维的单位矩阵，其他参数与式(30)相同。

将预设的式(31)参数应用到推导中，考虑到噪声的不可预测性，可得到如下表达式：

$\begin{array}{c}\mathrm{\Δ}\hat{y}(k+1|k)=A_1\mathrm{\Δy}\left(k\right)+A_2\mathrm{\Δy}(k-1)+B_{11}\mathrm{\Δu}\left(k\right)+B_{12}\mathrm{\Δu}(k-1)+B_{21}\mathrm{\Δ}{u}^2\left(k\right)+\\ B_{22}\mathrm{\Δ}{u}^2(k-1)+\mathrm{Dv}\left(k\right)\\ =S\left(2\right)\mathrm{\Δy}\left(k\right)+S\left(1\right)\·A_2\mathrm{\Δy}(k-1)+S\left(1\right)\·B_{12}\mathrm{\Δu}(k-1)+T\left(1\right)\mathrm{\Δu}\left(k\right)+\\ S\left(1\right)\·B_{22}\mathrm{\Δ}{u}^2(k-1)+T_g\left(1\right)\mathrm{\Δ}{u}^2\left(k\right)+\mathrm{Dv}\left(k\right)\end{array}---\left(32\right)$

$\begin{array}{c}\mathrm{\Δ}\hat{y}(k+2|k)=A_1\mathrm{\Δy}(k+1)+A_2\mathrm{\Δy}\left(k\right)+B_{11}\mathrm{\Δu}(k+1)+B_{12}\mathrm{\Δu}\left(k\right)+B_{21}\mathrm{\Δ}{u}^2(k+1)+\\ B_{22}\mathrm{\Δ}{u}^2\left(k\right)+\mathrm{Dv}\left(k\right)\\ =S\left(3\right)\mathrm{\Δy}\left(k\right)+S\left(2\right)\·A_2\mathrm{\Δy}(k-1)+S\left(2\right)\·B_{12}\mathrm{\Δu}(k-1)+T\left(2\right)\mathrm{\Δu}\left(k\right)+\\ T\left(1\right)\mathrm{\Δu}(k+1)+S\left(2\right)\·B_{22}\mathrm{\Δ}{u}^2(k-1)+T_g\left(2\right)\mathrm{\Δ}{u}^2\left(k\right)+\\ T_g\left(1\right)\mathrm{\Δ}{u}^2(k+1)+\mathrm{Dv}\left(k\right)\end{array}---\left(33\right)$

···

当代表预测步数的参数j满足j＜M(j＝2,3,···,P)时：

$\begin{array}{c}\mathrm{\Δ}\hat{y}(k+j|k)=S(j+1)\mathrm{\Δy}\left(k\right)+S\left(j\right)\·A_2\mathrm{\Δy}(k-1)+S\left(j\right)\·B_{12}\mathrm{\Δu}(k-1)+\\ \underset{i=1}{\overset{j}{\mathrm{\Σ}}}T\left(i\right)\mathrm{\Δu}(k+j-i)+S\left(j\right)\·B_{22}\mathrm{\Δ}{u}^2(k-1)+\\ \underset{i=1}{\overset{j}{\mathrm{\Σ}}}{T}_g\left(i\right)\mathrm{\Δ}{u}^2(k+j-i)+\mathrm{Dv}\left(k\right)\end{array}---\left(34\right)$

当j≥M(j＝2,3,···,P)时：

$\begin{array}{c}\mathrm{\Δ}\hat{y}(k+j|k)=S(j+1)\mathrm{\Δy}\left(k\right)+S\left(j\right)\·A_2\mathrm{\Δy}(k-1)+\\ S\left(j\right)\·B_{12}\mathrm{\Δu}(k-1)+\underset{i=j-M+1}{\overset{j}{\mathrm{\Σ}}}T\left(i\right)\mathrm{\Δu}(k+j-i)+\\ S\left(j\right)\·B_{22}\mathrm{\Δ}{u}^2(k-1)+\underset{i=j-M+1}{\overset{j}{\mathrm{\Σ}}}{T}_g\left(i\right)\mathrm{\Δ}{u}^2(k+j-i)+\mathrm{Dv}\left(k\right)\end{array}---\left(35\right)$

得到多变量系统多步预测输出增量矢量表达式：

ΔY(k)＝G₁₁ΔU(k)+G₁₂ΔU²(k)+G₂₁ΔU₂(k)+G₂₂ΔU₂ ²(k)+

F₁Δy(k)+F₂Δy(k-1)+Zv(k)  (36)

再由预测增量值经位移矩阵转换，得到输出量的预测值为：

Y(k)＝V·ΔY(k)+W·y(k)  (37)

式中：位移转换矩阵 $W={\left[\begin{array}{cccc}I& I& ...& I\end{array}\right]}_{(n_y\×P)\×n_y}^T,$ I为n_y×n_y的单位矩阵。

步骤五时滞误差校正

在图1中通过LSSVM-ARMAX组合模型1得到系统预测输出的值通过时滞环节2得到的实际预测输出值与水泥分解炉5的实际值存在一定误差，本发明通过预测输出与设定值的偏差来确定优化算法的表达式，为了使优化算法更准确通过预测误差修正预测输出，用柔化系数修正预测输出值。具体做法如下：

取第k时刻的预测误差：

$e\left(k\right)=y\left(k\right)-\hat{y}\left(k\right|k-1)---\left(38\right)$

用预测误差对系统未来预测时域内的预测输出值进行修正，为了能够调节修正程度，引入可以实时调节参数大小的校正参数矢量为H＝[h₁ h₂ ··· h_P]^T， $h_j=\mathrm{diag}\left\{\begin{array}{cccc}{h}_{j1}& h_{j2}& ...& h_{j{n}_y}\end{array}\right\},$ (j＝1,2,···,P)，对Y(k)进行修正。经过修正的预测输出值：

Y_c(k)＝Y(k)+He(k)  (39)

步骤六滚动优化计算

设水泥分解炉系统的参考轨迹如式(40)所示：

$Y_r={\left[\begin{array}{cccc}{y}_r^T(k+1)& y_r^T(k+2)& ...& y_r^T(k+P)\end{array}\right]}^T---\left(40\right)$

式中：y_r(k+j)＝Cy_r(k+j-1)+(I-C)y_sp，(j＝1,2,···,P)，y_r(k)＝y(k)代表多变量系统在k时刻的实际输出值， $C=\mathrm{diag}\left\{\begin{array}{cccc}{c}_1& c_2& ...& c_{n_y}\end{array}\right\}$ 为了减小设定值过分干扰而设定的柔化系数矩阵，柔化系数一般取值范围为：0≤c_i≤1，(i＝1,2,···,n_y)，y_sp＝[y_1sp y_2sp ··· y_nysp]^T代表多变量系统输出量的设定值。参考轨迹的引入如图1所示的参考轨迹3处。

在图1中所示组合模型预测控制4运算目的是求解目标函数，依据输入输出量的预测值和历史值构建带有约束的二次型目标函数，如式(41)所示：

$\left\{\begin{array}{c}\min J\left(k\right)={\left|\right|{\hat{Y}}_c\left(k\right)-Y_r\left(k\right)\left|\right|}_Q^2+{\left|\right|\mathrm{\ΔU}\left(k\right)\left|\right|}_R^2\\ s.\mathrm{t\Δ}{u}_{\min }\≤\mathrm{\Δu}(k+i)\≤\mathrm{\Δ}{u}_{\max }\\ u_{\min }\≤u(k+i)\≤u_{\max },(i=\mathrm{0,1},...,M-1)\\ y_{\min }\≤y(k+j)\≤y_{\max },(j=\mathrm{1,2},...,P)\end{array}\right.---\left(41\right)$

式中：Q是输出误差的权值矩阵取Q＝qI，R是输入增量的权值矩阵取R＝rI。

针对带有约束目标函数，本发明采用序列二次规划算法(SQP)求解目标函数。通过序列二次规划算法求出目标函数的最优解ΔU^*(k)，取该组解的第一组值即当前时刻的最优输入增量Δu(k)。将增量值Δu(k)作用到上一步的输入值上：

u(k)＝u(k-1)+Δu(k)  (42)

式中， $\mathrm{\Δu}\left(k\right)={\left[\begin{array}{cc}\mathrm{\Δ}{u}_1\left(k\right)& \mathrm{\Δ}{u}_2\left(k\right)...\mathrm{\Δ}{u}_{n_u}\left(k\right)\end{array}\right]}^T.$

得到下一步的最佳控制输入u(k)，将u(k)作用到系统中再进行下一步的优化求解，循环执行即可实现水泥分解炉系统组合模型预测控制的在线滚动优化求解。将每步得到的最优输入量输入到系统中，可以实现对输出值的预测控制，较好的控制效果是预测输出量向设定值靠拢，实际输出值围绕设定值上下波动，直到稳定在设定值。

从图3所示的水泥分解炉组合模型预测控制系统的方框图和图4所示的水泥分解炉组合模型预测控制系统的现场接线图可知，本发明的一种基于组合模型预测控制技术的水泥分解炉控制系统，其详细内容如下：

该水泥分解炉预测控制系统由测量仪表9、数据通讯接口10、中控机11和执行器12组成。测量仪表9和执行器12直接和水泥分解炉5相连，数据通讯接口10分别与测量仪表、执行器以及中控机相连，实现中控机与现场设备的数据通讯。

现场测量仪表9通过数据通讯接口10，将采集的水泥分解炉5的输出量分解炉出口温度13、一级筒O₂含量14和一级筒CO含量15的状态传输给中控机；中控机11中嵌入了上述的水泥分解炉组合模型的预测控制算法，通过此算法计算出水泥分解炉系统预测输入量；将预测出的水泥分解炉的输入量通过数据通讯接口传输给现场执行器，包括分解炉喂煤量16、高温风机转速17以及三次风管挡板开度18，作用到水泥分解炉5上；此时水泥分解炉的输出状态再由测量仪表9测量的分解炉状态反馈给中控机预测控制算法进行循环控制，构成水泥分解炉组合模型预测控制系统，使水泥分解炉在组合模型预测控制技术的调节下，实现实时智能控制的目的。

综上所述，整个控制方法和系统总结如下：

1.选取操作变量和被控变量；

2.通过中控机的数据通讯接口采集现场水泥分解炉系统的变量数据，并将其筛选分类；

3.选取LSSVM稳态模型和ARMAX动态模型的初始参数，然后分别通过基于PSO-LSSVM算法和递推增广最小二乘法对其进行离线模型训练，再将训练得到的LSSVM和ARMAX模型通过增益进行并联结合；

4.选择设定预测控制参数P，M，C，Q和输入输出变量的约束值Δu_min，Δu_max，u_min，u_max，y_min，y_max；在当前时刻k，计算出水泥分解炉系统中输出量分解炉出口温度y_1r(k+j)、一级筒O₂含量y_2r(k+j)和CO含量y_3r(k+j)；

5.计算水泥分解炉系统LSSVM-ARMAX组合模型输出，分解炉出口温度y_1c(k+j)、一级筒O₂含量y_2c(k+j)和CO含量y_3c(k+j)；

6.利用序列二次规划法求解二次型目标函数J，获得最优解Δu(k+j-1)，选取当前k时刻的输入量信号u(k)作为水泥分解炉系统中分解炉喂煤量、高温风机转速和三次风挡板开度的输入调节量。

7.将6中计算出的输入调解量通过水泥分解炉组合模型预测控制系统中的通讯接口下达到现场执行器，跳转到步骤4，完成系统的循环控制。

Claims (9)

1.一种基于组合模型预测控制技术的水泥分解炉控制方法，其特征在于该方法内容包括以下步骤：

步骤一：控制变量选取，分解过程的主要操作变量是分解炉喂煤量、高温风机转速和三次风管挡板开度，被控变量是分解炉出口温度、一级筒氧气含量和一级筒一氧化碳含量；

步骤二：数据采集及分类，通过中控机的数据通讯接口采集步骤一所述各操作变量和被控变量的现场数据，并按照阶跃响应特性将数据分类为稳态数据和动态数据；

步骤三：LSSVM-ARMAX组合模型辨识，利用基于粒子群优化参数的最小二乘支持向量机对分类出的稳态数据进行辨识得到水泥分解炉系统的非线性LSSVM稳态模型，利用递推增广最小二乘法对动态数据辨识得到水泥分解炉系统的线性ARMAX动态模型，并将稳态模型与动态模型有机并联组合，得到水泥分解炉系统的LSSVM-ARMAX组合模型；

步骤四：输出状态预测，在当前时刻k，利用水泥分解炉的输入输出信息和预测的未来输入信息，通过步骤三辨识得到的LSSVM-ARMAX组合模型，将其进行迭代预测出分解炉出口温度y_1n(k+j)、一级筒氧气含量y_2n(k+j)及一级筒一氧化碳含量y_3n(k+j)未来预测时域P的输出状态，式中j＝1,2,…,P；

步骤五：时滞及误差校正，根据步骤四中得到的各变量未来预测时域输出状态记为y_in(k+j)，进行时滞处理，得到各变量输出时滞预测状态y_in(k|k-1)，将其与k+1时刻的实际输出值相减求得输出误差e_i(k)，然后将y_in(k+j)与误差相加得到基于k时刻水泥分解炉未来预测输出y_ic(k+j)，上述各式中i＝1,2,3；

步骤六：滚动优化计算，引入参考输入轨迹分解炉出口温度y_1r(k+j)、一级筒氧气含量y_2r(k+j)以及一级筒一氧化碳含量y_3r(k+j)，将其与步骤五求得分解炉各未来预测输出y_ic(k+j)进行比较，构建带约束的二次型目标函数，应用序列二次规划法对其进行滚动优化求解，计算出当前时刻水泥分解炉被控变量分解炉喂煤量u₁(k)、高温风机转速u₂(k)和三次风管挡板开度u₃(k)；

步骤七：实时输出控制，根据步骤六计算出的当前时刻水泥分解炉被控变量控制现场执行器，通过滚动实时控制，使现场水泥分解炉被控变量实际输出与设定值y_sp相等，实现水泥分解炉的自动控制。

2.根据权利要求1所述的一种基于组合模型预测控制技术的水泥分解炉控制方法，其特征在于：步骤一中所述的操作变量：所述的操作变量：分解炉喂煤量、高温风机转速和三次风管挡板开度是影响水泥分解炉系统的主要输入控制变量，它们是分解炉内热量和气体的主要来源；所述的被控变量：分解炉出口温度、一级筒氧气含量和一级筒一氧化碳含量是影响水泥分解炉系统的主要输出控制变量，选取分解炉出口温度作为输出变量即能直接反应生料预分解所需热量，也能间接反应水泥生料预分解情况；分解炉出口温度的稳定可保证生料预分解过程的稳定有效运行；另选取一级筒氧气含量和一级筒一氧化碳含量既能反应分解炉内燃煤的燃烧状况，也兼顾了分解炉单产能耗和安全运转。

3.根据权利要求1所述的一种基于组合模型预测控制技术的水泥分解炉控制方法，其特征在于：步骤二中所述的数据采集是通过中控机的数据通讯接口采集步骤一所述各控制变量和被控变量的现场数据并进行存储；所述的数据分类是依据现场数据的阶跃响应特性将其分类为稳态数据和动态数据。

4.根据权利要求1所述的一种基于组合模型预测控制技术的水泥分解炉控制方法，其特征在于：步骤三中所述的水泥分解炉系统的LSSVM-ARMAX组合模型是通过非线性稳态模型的增益在线调节线性动态模型参数实现模型的有机并联组合，水泥分解炉的非线性LSSVM稳态模型是通过粒子群优化参数的最小二乘支持向量机根据分解炉的稳态数据辨识获得，水泥分解炉系统的线性ARMAX动态模型是通过递推增广最小二乘法对分解炉的动态数据辨识获得。

5.根据权利要求1所述的一种基于组合模型预测控制技术的水泥分解炉控制方法，其特征在于：在步骤六中所述的滚动优化计算是在滚动时域下采用序列二次规划法对带有约束的二次型目标函数进行求解，保证求出的分解炉预测输入解是在约束条件下的全局最优解。

6.根据权利要求1所述的一种基于组合模型预测控制技术的水泥分解炉控制方法，其特征在于：所述的基于粒子群优化参数的最小二乘支持向量机，其辨识算法采用高斯径向基核函数作为最小二乘支持向量机的核函数，采用粒子群优化算法对支持向量机正则化参数C与核宽度σ进行自动寻优，初始化粒子群中最大遗传代数maxgen＝20，种群规模sizepop＝30，学习因子c₁＝c₂＝2，算法中w为惯性权重，eps为适应度预设值，粒子群优化LSSVM参数算法的具体步骤如下：

1)初始化粒子群；

2)给定正则化系数C与核宽度σ的取值区间：

C∈[C_max,C_min]，σ∈[σ_max,σ_min]           (1)

3)确定正则化系数C与核宽度σ迭代速度的最大值与最小值，k为常数；

V_C max＝k·C_max

V_C min＝-k·C_min

                     (2)

V_d man＝k·σ_min

V_d min＝-k·σ_min

4)随机产生种群的初始位置pop和初始速度V：

pop(1,i)＝(C_max-C_min)*rand()+C_min

pop(2，i)＝(deta_max-deta_min)*rand()+deta_min

                                     (3)

V(1，i)＝Vcmax*rand()

V(2,i)＝Vdmax*rand()

5)LSSVM辨识，以测试集的预测输出与实际输出间的均方误差作为每个粒子的适应度fitness，全局极值点p_g和个体极值点p为：

p_g＝min(fitness)

                             (4)

p＝fitness

6)检查是否满足结束条件：p_g≤eps或已经达到最大迭代次数；如果满足，则跳转到步骤10)；

7)分别根据式(5)和式(6)更新粒子的速度与位置，如果超出范围则可利用程序将其限定在规定的范围内；

$v_{\mathrm{id}}^{(t+1)}=w*v_{\mathrm{id}}^{\left(t\right)}+c_1*r_1*(p_{\mathrm{id}}^{\left(t\right)}-x_{\mathrm{id}}^{\left(t\right)})+c_2*r_2*(p_{\mathrm{gd}}^{\left(t\right)}-x_{\mathrm{id}}^{\left(t\right)})---\left(5\right)$

$x_{\mathrm{id}}^{(t+1)}=x_{\mathrm{id}}^{\left(t\right)}+v_{\mathrm{id}}^{(t+1)}---\left(6\right)$

8)LSSVM辨识，以测试集的均方误差作为当前每个粒子的适应度值fitness，如果某个粒子的适应度值小于其个体极值，则对其个体极值进行更新；如果某个粒子的适应度值小于全局极值，则对全局极值进行更新，同时更新最优解正则化系数C与核宽度σ组成的向量global_x；

9)迭代次数加1，跳转到步骤6)；

10)得到最优解global_x；

11)采用已优化的正则化系数C与核宽度σ进行LSSVM辨识，优化算法结束。

7.根据权利要求1所述的一种基于组合模型预测控制技术的水泥分解炉控制方法，其特征在于：所述的递推增广最小二乘法，其算法步骤如下：

1)采集输入输出数据，设数据共L组，本发明具体实施方式中现场采集的水泥分解炉动态输入输出数据共计9823组；

2)设定初值，对于协方差阵初始值P(0)可以先选择一个大正定单位矩阵本发明实例中设定P(0)＝10⁶*I₀，参数估计初值可选取一个很小的实相量，如其中是单位均为1的n_y×(n_yn_a+n_un_b+n_yn_d)维向量；

3)按照式(7)计算协方差矩阵P(t)由式(8)求得信息向量φ(t)；

$P\left(t\right)=[I-L\left(t\right){\widehat{\mathrm{\φ}}}^T\left(t\right)]P(t-1)\∈R^{(n_y{n}_a+n_u{n}_b+n_y{n}_d)\×(n_y{n}_a+n_u{n}_b+n_y{n}_d)}---\left(7\right)$

$\begin{array}{c}\widehat{\mathrm{\φ}}\left(t\right)=[-y(t-1),-y(t-2),...-y(t-n_a),u(t-1),u(t-2),...u(t-n_b),\\ \hat{v}(t-1),\hat{v}(t-2),...\hat{v}(t-n_d){]}^T\end{array}---\left(8\right)$

4)利用式(9)获得递推估计参数；

5)如果t＝L+1则辨识完成，得到最终辨识结果为否则t+1，重复3—5步。

8.实现权利要求1所述一种基于组合模型预测控制技术的水泥分解炉控制方法的控制系统，其特征在于：该控制系统包括测量仪表、数据通讯接口、中控机和执行器；

所述的测量仪表用于测量水泥分解炉系统的输出量，即测量分解炉出口温度、一级筒气体含量和一级筒一氧化碳含量；

所述的数据通讯接口用于将现场测量仪表测量的数据传输到中控机，并将中控机处理后下达的输入量变化增量传输到现场执行器；

所述的中控机用于运行水泥分解炉组合模型预测控制算法，根据水泥分解炉系统的实际输出量计算出当前时刻应加入分解炉内燃煤用量、高温风机转速以及三次风管挡板开度的控制量；

所述的执行器通过数据通讯接口接收中控机下达的变量调节量指令，通过其自动调节，实现现场分解炉的优化自动控制。

9.根据权利要求6所述的基于LSSVM-ARMAX组合模型的水泥分解炉预测控制系统，其特征在于：测量仪表和执行器直接和水泥分解炉相连，数据通讯接口分别与测量仪表、执行器以及中控机相连，实现中控机与现场设备的数据通讯。