一种氨法脱硫系统的扰动抑制预测控制方法与流程

本发明涉及热工自动控制领域，特别是涉及一种氨法脱硫系统的扰动抑制预测控制方法。

背景技术：

烟气氨法脱硫作为一项较新的环保技术，虽然具有很多其它脱硫方法无可比拟的优点，但同样也面临着许多问题：氨作为危险品，若加氨量或加氨的方式不当，会造成二次污染并影响脱硫效率；增加吸收剂氨浓度有利于提高脱硫效率，但同时也会增加氨逃逸量；浓缩液过饱和将不利于so2的吸收且容易结晶堵塞装置；氨法脱硫系统是一个典型的具有非线性及大惯性特性，干扰量较多的系统。因此，操作和控制好电厂烟气脱硫装置，可以保证电厂的达标排放和经济运行。

氨法脱硫系统中，需要操作和控制的变量较多，其中最为突出的就是循环浆液ph值。ph值反应了吸收塔中脱硫剂的酸碱度，它的大小影响到反应浆液的品质，间接影响到脱硫效率及副产物品质。整个脱硫过程非常复杂，加氨后对脱硫效率和循环液ph值的影响具有很大的惯性和滞后；同时氨法脱硫过程是一个化学反应过程，具有强烈的非线性特征，易受不可测扰动的影响。总之，氨法脱硫系统中，惯性大，且烟气流量、烟气so2浓度、氨水浓度等干扰众多。在这种情况下采用传统的pid控制方案，根本无法解决控制系统稳定性和控制系统品质之间的矛盾，必然会引起控制系统的不稳定和振荡。尤其是在受到未知扰动影响时，常规控制方案更是难以取得令人满意的控制效果。

技术实现要素：

发明目的：本发明的目的是提供一种氨法脱硫系统的扰动抑制预测控制方法，能够解决系统稳定性和系统品质之间的矛盾。

技术方案：为达到此目的，本发明采用以下技术方案：

本发明所述的氨法脱硫系统的扰动抑制预测控制方法，包括以下步骤：

s1：以脱硫塔的氨水流量作为控制量，以循环浆液ph作为被控量，在稳态工况下进行开环阶跃响应试验，获取控制量与被控量之间的传递函数，得到被控对象的辨识模型；

s2：将循环浆液ph对氨水流量的传递函数转换为带扰动项的系统状态空间模型；

s3：选取阶跃输出扰动模型、阶跃输入扰动模型、斜坡输入扰动模型和周期输入扰动模型这四个扰动模型按顺序组成扰动模型集；针对扰动模型集中的每一个扰动模型，将步骤s2得到的系统状态空间模型转换为增量状态空间模型，再将扰动增量项和系统输出项作为扩增量对状态变量进行扩增，得到扩增状态变量和扩增状态空间模型；

s4：根据步骤s3中每一个扰动模型的特点，得出所述四个扰动模型的扩增状态空间模型中的系数矩阵；

s5：利用标准kalman滤波器对每个扩增状态空间模型进行状态估计；

s6：计算脱硫系统的实际输出与kalman滤波器得到的预估输出之间的偏差；

s7：采用递推贝叶斯公式计算出各个扰动模型与当前脱硫系统的匹配概率；

s8：根据各个扰动模型与当前脱硫系统的匹配概率计算各个扩增状态空间模型的权重；

s9：根据各个扩增状态空间模型的权重，对各个扩增状态空间模型进行加权平均，计算出脱硫系统的p步平均预测输出值；

s10：根据扩增状态空间模型构建平均预测模型，计算最优控制增量向量；

s11：计算得到脱硫系统的输出。

进一步，所述步骤s2中的系统状态空间模型如式(1)所示：

式(1)中，xd(k)表示k时刻的系统状态变量，xd(k+1)表示k+1时刻的系统状态变量，u(k)表示k时刻的控制量，y(k)表示k时刻的被控量，d(k)表示k时刻的扰动量，ad表示xd(k)的系数矩阵，bd表示u(k)的系数矩阵，brd表示d(k)的系数矩阵，cd表示xd(k)的系数矩阵。

进一步，所述步骤s3中的扩增状态空间模型如式(2)所示：

式(2)中，xi(k+1)(i＝1,2,3,4)表示k+1时刻第i个扰动模型扩增的状态向量，包括三个子向量δxd(k+1)、δd(k+1)和y(k+1)，δxd(k+1)表示k+1时刻系统原来的状态增量，δd(k+1)表示k+1时刻当前扰动模型的扰动增量，y(k+1)表示k+1时刻脱硫系统的被控量；xi(k)表示k时刻第i个扰动模型扩增的状态向量，包括三个子向量δxd(k)、δd(k)和y(k)，δxd(k)表示k时刻系统原来的状态增量，δd(k)表示k时刻当前扰动模型的扰动增量，y(k)表示k时刻脱硫系统的被控量；δu(k)表示k时刻脱硫系统的控制量，aai表示第i个扰动模型的扩增状态空间模型的第一系数矩阵，a1,...,a9表示矩阵aai中的元素，bai表示第i个扰动模型的扩增状态空间模型的第二系数矩阵，b1,...,b3表示矩阵bai中的元素，cai表示第i个扰动模型的扩增状态空间模型的第三系数矩阵，c1,...,c3表示矩阵cai中的元素。

进一步，所述步骤s4中，四个扰动模型的扩增状态空间模型中的系数矩阵aai、bai和cai通过式(3)-(6)得到：

其中，aa1表示阶跃输出扰动模型的扩增状态空间模型的第一系数矩阵，aa2表示阶跃输入扰动模型的扩增状态空间模型的第一系数矩阵，aa3表示斜坡输入扰动模型的扩增状态空间模型的第一系数矩阵，aa4表示周期输入扰动模型的扩增状态空间模型的第一系数矩阵，ba1表示阶跃输出扰动模型的扩增状态空间模型的第二系数矩阵，ba2表示阶跃输入扰动模型的扩增状态空间模型的第二系数矩阵，ba3表示斜坡输入扰动模型的扩增状态空间模型的第二系数矩阵，ba4表示周期输入扰动模型的扩增状态空间模型的第二系数矩阵，ca1表示阶跃输出扰动模型的扩增状态空间模型的第三系数矩阵，ca2表示阶跃输入扰动模型的扩增状态空间模型的第三系数矩阵，ca3表示斜坡输入扰动模型的扩增状态空间模型的第三系数矩阵，ca4表示周期输入扰动模型的扩增状态空间模型的第三系数矩阵。

进一步，所述步骤s5中，状态估计通过式(7)-(12)实现：

p(k+1|k)＝aaip(k|k)aai^t+q(8)

l(k+1)＝p(k+1|k)cai^t[caip(k+1|k)cai^t+r]^-1(9)

p(k+1|k+1)＝[i-l(k+1)cai]p(k+1|k)(11)

其中，表示k时刻对k+1时刻的第i个扰动模型的扩增状态空间模型状态量的预估值，表示k时刻对k时刻的第i个扰动模型的扩增状态空间模型状态量的预估值，表示k+1时刻对k+1时刻的第i个扰动模型的扩增状态空间模型状态量的预估值，表示k时刻对k+1时刻的第i个扰动模型的扩增状态空间模型输出量的预估值，p(k+1|k)表示k时刻对k+1时刻的协方差矩阵，p(k|k)表示k时刻对k时刻的协方差矩阵，p(k+1|k+1)表示k+1时刻对k+1时刻的协方差矩阵，l(k+1)表示k+1时刻的状态估计增益矩阵，y(k+1)表示k+1时刻的脱硫系统的实际输出，q表示输入噪声项的协方差，r表示输出噪声项的协方差。

进一步，所述步骤s6中，脱硫系统的实际输出与kalman滤波器的预估输出之间的偏差通过式(13)得到：

式(13)中，ei(k+1)表示k+1时刻脱硫系统的实际输出与第i个扰动模型的kalman滤波器的预估输出之间的偏差。

进一步，所述步骤s7中，各个扰动模型与当前脱硫系统的匹配概率通过式(14)得到：

式(14)中，ρi(k+1)表示k+1时刻第i个扰动模型与当前脱硫系统的匹配概率，ρi(k)表示k时刻第i个扰动模型与当前脱硫系统的匹配概率，ei(k+1)表示表示k+1时刻脱硫系统的实际输出与第i个扰动模型的kalman滤波器的预估输出之间的偏差，λ表示贝叶斯概率协方差矩阵；设定匹配概率的下限值δ，当k时刻的匹配概率ρi(k)低于下限值δ时，先以该下限值δ置换ρi(k)后，再代入式(14)计算出k+1时刻的匹配概率ρi(k+1)。

进一步，所述步骤s8中，通过式(15)计算各个扩增状态空间模型的权重：

式(15)中，wi(k)表示k时刻第i个扰动模型的扩增状态空间模型的权重。

进一步，所述步骤s9中，脱硫系统的p步平均预测输出值通过式(16)得到：

式(16)中，表示k时刻对未来j个时刻的第i个扩增状态空间模型的预测输出值，wi(k)表示k时刻第i个扰动模型的扩增状态空间模型的权重，表示k时刻对未来j个时刻加权后的平均扩增状态空间模型的预测输出值。

进一步，所述步骤s10中的平均预测模型如式(17)所示：

式(17)中，为k时刻循环浆液ph值的p步平均预测输出向量，δu(k)为k时刻控制量向量，f为加权后的平均预测模型的状态变量的系数矩阵，wi(k)表示k时刻第i个扰动模型的扩增状态空间模型的权重，fi表示第i个扰动模型的扩增状态空间模型的状态向量的系数矩阵，φ为加权后的平均预测模型的控制量向量的系数矩阵，φi表示第i个扰动模型的扩增状态空间模型的控制量向量的系数矩阵；p表示预测步长，m表示控制步长；

所述最优控制增量向量通过以下过程计算得到：基于性能指标函数计算k时刻控制量增量δu(k)，进而基于k时刻控制量增量δu(k)构造控制增量向量δu，制定最终性能指标函数来计算最优控制增量向量，最终性能指标函数和最优控制增量公式如式(18)：

式(18)中，yr为浆液ph值设定值，wy为误差权矩阵，j为性能指标，x(k)为k时刻的平均预测模型的状态量。

有益效果：本发明公开了一种氨法脱硫系统的扰动抑制预测控制方法，利用扩增状态空间模型，建立扰动模型集，结合贝叶斯概率加权的方法，得到平均预测模型，通过预测控制方法，采用新的目标函数，实现进氨水流量对循环浆液ph值的良好控制，能够很好地应对氨法脱硫过程的大迟延、大惯性特性。同时考虑不可测扰动对脱硫系统的影响，使系统获得了更好的抗干扰能力，能够在受不可测扰动影响下，快速识别扰动并抑制扰动对系统的影响，总体提高了脱硫系统的控制品质，进一步促进了脱硫系统经济环保运行。

附图说明

图1为本发明的氨法脱硫系统的扰动抑制预测控制方法的结构框图；

图2为本发明与其它扰动抑制控制算法在循环浆液ph值受阶跃输出不可测扰动时的控制效果对比图；

图3为本发明与其它扰动抑制控制算法在循环浆液ph值受阶跃输入不可测扰动时的控制效果对比图；

图4为本发明与其它扰动抑制控制算法在循环浆液ph值受斜坡输入不可测扰动时的控制效果对比图；

图5为本发明与其它扰动抑制控制算法在循环浆液ph值受周期输入不可测扰动时的控制效果对比图；

图6为本发明与其它扰动抑制控制算法在循环浆液ph值同时受阶跃和周期输入不可测扰动时的控制效果对比图。

具体实施方式

下面结合具体实施方式对本发明的技术方案作进一步的介绍。

本具体实施方式公开了一种氨法脱硫系统的扰动抑制预测控制方法，包括以下步骤：

s1：以脱硫塔的氨水流量作为控制量，以循环浆液ph作为被控量，在稳态工况下进行开环阶跃响应试验，获取控制量与被控量之间的传递函数，得到被控对象的辨识模型；

s2：将循环浆液ph对氨水流量的传递函数转换为带扰动项的系统状态空间模型；

s3：选取阶跃输出扰动模型、阶跃输入扰动模型、斜坡输入扰动模型和周期输入扰动模型这四个扰动模型按顺序组成扰动模型集；针对扰动模型集中的每一个扰动模型，将步骤s2得到的系统状态空间模型转换为增量状态空间模型，再将扰动增量项和系统输出项作为扩增量对状态变量进行扩增，得到扩增状态变量和扩增状态空间模型；

s4：根据步骤s3中每一个扰动模型的特点，得出所述四个扰动模型的扩增状态空间模型中的系数矩阵；

s5：利用标准kalman滤波器对每个扩增状态空间模型进行状态估计；

s6：计算脱硫系统的实际输出与kalman滤波器得到的预估输出之间的偏差；

s7：采用递推贝叶斯公式计算出各个扰动模型与当前脱硫系统的匹配概率；

s8：根据各个扰动模型与当前脱硫系统的匹配概率计算各个扩增状态空间模型的权重；

s9：根据各个扩增状态空间模型的权重，对各个扩增状态空间模型进行加权平均，计算出脱硫系统的p步平均预测输出值；

s10：根据扩增状态空间模型构建平均预测模型，计算最优控制增量向量；

s11：计算得到脱硫系统的输出。

步骤s2中的系统状态空间模型如式(1)所示：

式(1)中，xd(k)表示k时刻的系统状态变量，xd(k+1)表示k+1时刻的系统状态变量，u(k)表示k时刻的控制量，y(k)表示k时刻的被控量，d(k)表示k时刻的扰动量，ad表示xd(k)的系数矩阵，bd表示u(k)的系数矩阵，brd表示d(k)的系数矩阵，cd表示xd(k)的系数矩阵。

步骤s3中的扩增状态空间模型如式(2)所示：

式(2)中，xi(k+1)(i＝1,2,3,4)表示k+1时刻第i个扰动模型扩增的状态向量，包括三个子向量δxd(k+1)、δd(k+1)和y(k+1)，δxd(k+1)表示k+1时刻系统原来的状态增量，δd(k+1)表示k+1时刻当前扰动模型的扰动增量，y(k+1)表示k+1时刻脱硫系统的被控量；xi(k)表示k时刻第i个扰动模型扩增的状态向量，包括三个子向量δxd(k)、δd(k)和y(k)，δxd(k)表示k时刻系统原来的状态增量，δd(k)表示k时刻当前扰动模型的扰动增量，y(k)表示k时刻脱硫系统的被控量；δu(k)表示k时刻脱硫系统的控制量，aai表示第i个扰动模型的扩增状态空间模型的第一系数矩阵，a1,...,a9表示矩阵aai中的元素，bai表示第i个扰动模型的扩增状态空间模型的第二系数矩阵，b1,...,b3表示矩阵bai中的元素，cai表示第i个扰动模型的扩增状态空间模型的第三系数矩阵，c1,...,c3表示矩阵cai中的元素。

步骤s4中，四个扰动模型的扩增状态空间模型中的系数矩阵aai、bai和cai通过式(3)-(6)得到：

其中，aa1表示阶跃输出扰动模型的扩增状态空间模型的第一系数矩阵，aa2表示阶跃输入扰动模型的扩增状态空间模型的第一系数矩阵，aa3表示斜坡输入扰动模型的扩增状态空间模型的第一系数矩阵，aa4表示周期输入扰动模型的扩增状态空间模型的第一系数矩阵，ba1表示阶跃输出扰动模型的扩增状态空间模型的第二系数矩阵，ba2表示阶跃输入扰动模型的扩增状态空间模型的第二系数矩阵，ba3表示斜坡输入扰动模型的扩增状态空间模型的第二系数矩阵，ba4表示周期输入扰动模型的扩增状态空间模型的第二系数矩阵，ca1表示阶跃输出扰动模型的扩增状态空间模型的第三系数矩阵，ca2表示阶跃输入扰动模型的扩增状态空间模型的第三系数矩阵，ca3表示斜坡输入扰动模型的扩增状态空间模型的第三系数矩阵，ca4表示周期输入扰动模型的扩增状态空间模型的第三系数矩阵。

步骤s5中，状态估计通过式(7)-(12)实现：

p(k+1|k)＝aaip(k|k)aai^t+q(8)

l(k+1)＝p(k+1|k)cai^t[caip(k+1|k)cai^t+r]^-1(9)

p(k+1|k+1)＝[i-l(k+1)cai]p(k+1|k)(11)

其中，表示k时刻对k+1时刻的第i个扰动模型的扩增状态空间模型状态量的预估值，表示k时刻对k时刻的第i个扰动模型的扩增状态空间模型状态量的预估值，表示k+1时刻对k+1时刻的第i个扰动模型的扩增状态空间模型状态量的预估值，表示k时刻对k+1时刻的第i个扰动模型的扩增状态空间模型输出量的预估值，p(k+1|k)表示k时刻对k+1时刻的协方差矩阵，p(k|k)表示k时刻对k时刻的协方差矩阵，p(k+1|k+1)表示k+1时刻对k+1时刻的协方差矩阵，l(k+1)表示k+1时刻的状态估计增益矩阵，y(k+1)表示k+1时刻的脱硫系统的实际输出，q表示输入噪声项的协方差，r表示输出噪声项的协方差。

步骤s6中，脱硫系统的实际输出与kalman滤波器的预估输出之间的偏差通过式(13)得到：

式(13)中，ei(k+1)表示k+1时刻脱硫系统的实际输出与第i个扰动模型的kalman滤波器的预估输出之间的偏差。

步骤s7中，各个扰动模型与当前脱硫系统的匹配概率通过式(14)得到：

式(14)中，ρi(k+1)表示k+1时刻第i个扰动模型与当前脱硫系统的匹配概率，ρi(k)表示k时刻第i个扰动模型与当前脱硫系统的匹配概率，ei(k+1)表示表示k+1时刻脱硫系统的实际输出与第i个扰动模型的kalman滤波器的预估输出之间的偏差，λ表示贝叶斯概率协方差矩阵；设定匹配概率的下限值δ，当k时刻的匹配概率ρi(k)低于下限值δ时，先以该下限值δ置换ρi(k)后，再代入式(14)计算出k+1时刻的匹配概率ρi(k+1)。

步骤s8中，通过式(15)计算各个扩增状态空间模型的权重：

式(15)中，wi(k)表示k时刻第i个扰动模型的扩增状态空间模型的权重。

步骤s9中，脱硫系统的p步平均预测输出值通过式(16)得到：

式(16)中，表示k时刻对未来j个时刻的第i个扩增状态空间模型的预测输出值，wi(k)表示k时刻第i个扰动模型的扩增状态空间模型的权重，表示k时刻对未来j个时刻加权后的平均扩增状态空间模型的预测输出值。

步骤s10中的平均预测模型如式(17)所示：

式(17)中，为k时刻循环浆液ph值的p步平均预测输出向量，δu(k)为k时刻控制量向量，f为加权后的平均预测模型的状态变量的系数矩阵，wi(k)表示k时刻第i个扰动模型的扩增状态空间模型的权重，fi表示第i个扰动模型的扩增状态空间模型的状态向量的系数矩阵，φ为加权后的平均预测模型的控制量向量的系数矩阵，φi表示第i个扰动模型的扩增状态空间模型的控制量向量的系数矩阵；p表示预测步长，m表示控制步长；

所述最优控制增量向量通过以下过程计算得到：基于性能指标函数计算k时刻控制量增量δu(k)，进而基于k时刻控制量增量δu(k)构造控制增量向量δu，制定最终性能指标函数来计算最优控制增量向量，最终性能指标函数和最优控制增量公式如式(18)：

式(18)中，yr为浆液ph值设定值，wy为误差权矩阵，j为性能指标，x(k)为k时刻的平均预测模型的状态量。

利用本发明中氨法脱硫系统的扰动抑制预测控制方法对氨法脱硫系统进行控制，在系统受不可测扰动影响时，循环浆液ph值控制效果与其它扰动抑制方法对比如图2-6所示。所有图示中，循环浆液ph值的设定值为5.5，整个过程持续1000s。图2中，在1000s时循环浆液ph侧有扰动进入系统情况下，本发明能够快速识别扰动类型为输出扰动，提前改变氨水流量，快速抑制扰动，达到很好的控制效果。与其它扰动抑制方法相比，本发明有更快的调节能力，且控制量波动较小。图3中，在1000s时氨水流量侧有扰动进入系统情况下，本发明能够快速识别扰动类型为阶跃输入扰动，并快速抑制扰动。与其它扰动抑制方法相比，本发明有更快的调节能力和更小的动态偏差。图4中，在1000s时氨水流量侧有扰动进入系统情况下，本发明能够快速识别扰动类型为斜坡输入扰动，并在扰动持续过程中将其抑制消除，取得很好的控制效果。但是在斜坡扰动影响系统情况下，其它扰动抑制方法并不能有效抑制扰动，有较大的动态偏差。图5中，氨水流量侧受持续扰动影响，本发明能够快速识别出扰动类型为周期输入扰动，改变氨水流量，有效抑制扰动。与其它扰动抑制方法相比，本发明有较小的动态偏差。图6中，氨水流量侧受多种扰动影响，本发明能识别出扰动类型主要为阶跃输入和周期输入扰动，并首先完全抑制阶跃扰动的影响，再抑制周期扰动，有较好的抑制效果。与其它扰动抑制方法相比，有较小的动态偏差。