广义预测控制在机组控制系统中的仿真研究
摘要：由于单元机组是一个复杂的系统，一般都存在着很大的延迟和惯性，被控对象具有不确定性和强耦合，常规的PID 控制器在处理大延迟对象上很难获得较好的控制效果。我们将广义预测算法应用到机组控制中，利用MA TLAB 仿真，对响应速度和鲁棒性进行分析，通过比较仿真效果表明广义预测控制方法对系统具有良好的适应性和鲁棒性。 关键词：预测控制 控制系统 仿真
0 引言
广义预测控制（GPC ）是在最小方差控制的基础上，汲取动态矩阵控制（DMC ）和模型算法控制（MAC ）中的多步预测优化策略，随着自适应控制的研究而发展起来的一种预测控制方法。此种控制算法的数学模型是通过推导控制律参数而获取的，在理论研究和实际应用中都得到了重视[1]。本文采用广义预测控制算法对500MW 机组80%负荷点处定压运行进行仿真研究。
1 多变量广义预测控制算法
预测控制算法是在工业生产过程和相关的理论知识结合的基础上提出来的，与传统的PID 控制方法相比较，预测控制算法具有自身的优越性。
广义预测控制具有预测控制算法的基本特征，因而具有优良的控制性能和鲁棒性，使系统的设计更加灵活，被认为是在复杂工业控制领域具有代表性的预测控制算法之一。广义预测控制算法采用的预测模型是具有一定结构和参数的最小化模型，因此对系统输出进行预测时需要已知模型结构，这给预测模型的建立增加了难度，但所需确定的参数较少，减少了计算量。
在GPC 中，采用了最小方差控制中所用的受控自回归积分滑动平均模型来描述收到随机干扰的对象[2]
C (q -1) ξ(t )
（1） A （q ）y (t ) =B (q ) u (t -1) +
∆
-1
-1
其中
A (q -1) =1+a 1q -1+... +a n a q -n a B (q -1) =b 0+b 1q -1+... +b n b q -n b A (q -1) =c 0+c 1q -1+... +c n c q -n c
式中，q 是后移算子，表示后退一个采样周期的相应的量；∆=1-q -1为差分算子；
-1
ξ（t ）是一个不相关的随机序列，表示一类随机噪声的影响。A 、B 、C 都是q 的多项式，
-1
其中多项式B (q ) 的若干首项元素b 0, b 1,... 可以是零，以表示对象相应的时滞数。为了突出方法原理，这里假设C (q -1) =1。
式（1）实际上是用Z 传递函数给出了对象的描述，由输入u 到输出y 间的Z 传递函数为
-1
z -1B (z -1)
（2） G (z ) =
A (z -1)
-1
为了利用模型（1）导出j 步后输出y(t+j)的预测值，首先考虑下述丢翻图（Diophantine ）方程
1=E j (q -1) A ∆+q -1F j (q -1) （3）
其中，E j ，F j 是由A (q ) 和预测长度j 唯一确定的
-1
E j (q -1) =e j . 0+e j . 1q -1+... +e j . j -1q -(j -1) F j (q ) =f j . 0+f j . 1q +... +f j . na q
这样，在式（1）两端乘以E j ∆q 后可得
j -1
-1
-na
E j A ∆y (t +j ) =E j B ∆u (t +j -1) +E j ξ(t +j )
利用式（3），可以写出t+j时刻的输出量
y (t +j ) =E j B ∆u (t +j -1) +F j y (t ) +E j ξ(t +j )
注意到E j ，F j 的形式，可以知道E j B ∆u (t +j -1) 与u (t +j -1) ，u (t +j -2),... 有关，
F j y (t ) 只与y(t)，y(t-1),...有关，而E j ξ(t +j ) 与ξ(t +j ),..., ξ(t +1) 有关[3]。由于在t 时刻
未来的噪声ξ(t +i ), i =1,... j 都是未知的，所以对y(t+j)最合适的预测值可以由下式得到
Λ
y (t +j |t ) =E j B ∆u (t +j -1) +F j y (t ) （4）
这样，根据已知的输入输出信息及未来的的输入值，就可以预测对象未来的输出，式（4）就是GPC 的预测模型。
在GPC 中，t 时刻的优化性能指标具有以下形式：
NU ⎧N 2
2⎫2
min J (t ) =E ⎨∑(y (t +j ) -ω(t +j )) +∑λ(j ) [∆u (t +j -1) ]⎬ （5）
j =1⎭⎩j =N 1
其中，E 为数学期望；ω为对象输出的期望值；N 1和N 2分别为优化时域的识值和终值；NU 为控制时域，即在NU 步后控制量不再变化：u (t +j -1) =u (t +NU -1), j >NU
λ(j ) 为控制加权系数，为简化设计一般常假设其为常数λ。
2 机组受控对象模型
在单元机组协调控制系统中，包括主控系统和子控系统，它们最终的控制对象为锅炉、汽轮机和发电机[4]，它具有大时滞、强耦合、非线性，同时难以平匹配的特点。在整个机组的发电过程中，锅炉和汽轮机的动态特性差别很大，它们对负荷的响应速度一慢一快，所以单元机组内外的两个能量平衡关系互相制约，而它内部的负荷响应性能又与内部参数稳定性存在固有固有的矛盾。因此，在单元机组实施负荷控制的时候，如何协调控制相当重要。为了便于对系统的控制进行分析，在保证满足控制要求的前提下，我们将受控对象的模型简化为两输入两输出的系统，如下图1所示：
图1 单元机组受控对象模型
上图由输入控制和传递函数以及过程输出组成，并且此模型的动态特性由线性常系数多变量传递函数式（6）和（7）来描述[5]：
Y (s ) =G (s ) ⋅U (s ) （6）
⎡y 1(s ) ⎤⎡u 1(s ) ⎤⎡W 11(s ) W 12(s ) ⎤
Y (s ) =U (s ) = 其中⎢y (s ) ⎥，⎢u (s ) ⎥，G (s ) =⎢⎥ （7）
W (s ) W (s ) ⎣2⎦⎣2⎦22⎣21⎦
式中Y (s ) 是被控制矢量，U (s ) 是控制矢量，G (s ) 是传递函数矩阵。
我们选取的动态模型是某厂300MW 直流锅炉燃煤机组在100%负荷点上的数学模型[6]，模型如下[7]：
4. 665s (99s +1) ⎡
⎡N e ⎤⎢582s 2+50s +1(4. 1s +1) ⎢P ⎥=⎢
⎣t ⎦⎢-1. 42(0. 04+0. 96)
⎢70s +1⎣
2. 069(311s +1) ⎤
⎡u t ⎤149s +12(22. 4s +1) ⎥
⎥⋅⎥
1. 265(205s +1) ⎥⎢⎣u b ⎦ （8） (128s +1) 2(11. 7s +1) ⎥⎦
由图1可以知道单元机组复杂的多输入多输出的控制系统简化成了有汽轮机调节门指令u t 和锅炉的燃烧率指令u b 为输入，实发电功率N e 和机前压力P t 为输出的模型。在公式（8）中，u t 和u b 对应于图1 中的u 1和u 2，N e 和P t 对应于图1中的y1(t)和y2(t)。
3 仿真
控制时域长度M 越大，则控制的机动性变强，控制的灵敏度变高，但系统地稳定性和鲁棒性就会下降，为了提高系统的稳定性和鲁棒性，要求M 小，M 越小，远程跟踪控制能力虽有所削弱，但可得到一个稳定的控制[8]，因此，M 的选择应兼顾快速性和稳定性两种，平衡考虑，本文仿真时选取M=6。预测时域P=10，仿真时间为500s ，利用如下程序[9]： clear all
g11=poly2tfd([461.835 4.665 0],conv([58 50 1],[4.1 1]),0,0);
g12=poly2tfd([643.459 2.069],conv(conv([149 1],[149 1]),[22.4 1]),0,0); g22=poly2tfd([259.325 1.265],conv(conv([128 1],[128 1]),[11.7 1]),0,0); g21=poly2tfd([-3.92 -1.42],[70 1],0,0); gd1=poly2tfd(3.8,[14.9 1],0,0); gd2=poly2tfd(4.9,[13.2 1],0,0);
tfinal=500; %Model horizon N delt=10; %sampling period ny=2; %Number of outputs model=tfd2step(tfinal,delt,ny,g11,g21,g12,g22);
plant=model; %No plant/model mismatch dmodel=[]; %Default disturbance model dplant=tfd2step(tfinal,delt,ny,gd1,gd2); P=10;M=6; %Horizons
ywt=[1 1]; uwt=[1 1]; %Q and R
tend=450; %final time for simulation r=[10 0]; %Set-point change inXD XB a=zeros([1,tend]); for i=100:1:tend
a(i)=1; % 3 %step in F at t=50 min. end dstep=[a']; ulim=[];%u limits ylim=[]; tfilter=[]; dplant=[]; dmodel=[]; dstep=[];
[y,u,ym]=cmpc(plant,model,ywt,uwt,M,P,tend,r,ulim,ylim,tfilter,dplant,dmodel,dstep); figure(1) subplot(211) plot(y) grid on legend('Ne','Pt')
ylabel('Ne /% Pt /%') xlabel('k(采样时间10s)') subplot(212) plot(u) legend('Ut','Ub') ylabel('Ub /% Ut /%') xlabel('k(采样时间10s)') grid on
仿真结果如下图2所示：
15N e /% P t /%
1050-5
51015
202530k(采样时间10s)
35404550
15U b /% U t /%
10
5
00
51015
202530k(采样时间10s)
35404550
图2 100%负荷下功率响应曲线和气压响应曲线（上图）、汽机调门指令和锅炉燃烧率指
令响应曲线（下图）
4 结论
随着电力工业的发展，如何将一些先进的控制技术应用于单元机组负荷控制，以提高单元机组负荷控制系统的性能，是一个具有很大应用研究的价值的课题。本文采用了多变量广义预测控制的算法，将其应用于单元机组负荷控制中，进行仿真研究，结果表明该方法在线计算量小，实时性好，参数调整灵活，改善了系统的稳定性和动态特性。
参考文献
[1] 舒迪前. 预测控制系统及其应用[M].北京:机械工业出版社，1996.
[2] 曹楠. 多变量广义预测控制在单元机组协调控制系统中的应用研究[D] 北京:华北电力 大学,2000.
[3] 王伟. 一种广义预测自适应控制的直接方法[J]. 自动化学报 , 1996, 22(3) : 270- 276. [4] 凌呼君, 白雄怀. 火电单元机组机炉协调的广义预测控制策略[J].控制与决策, 2003, 18(增刊 1 ) : 31- 33.
[5] 吴宝琴. 广义预测控制在协调控制系统中的应用[D].保定：华北电力大学，2006. [6] 杜之正. 预测控制在电厂协调控制系统中的研究[D].保定：华北电力大学，2012. [7] 雎刚, 韦红旗, 陈绍炳, 徐治皋. 单元机组负荷多变量模型预测控制[J].中国电机工程学报,2002,04:145-149.
[8] 罗辑，曹建国等. 基于广义预测控制的综合控制算法设计与仿真[J].自动化技术与应用 ， 2004,23（2）：10-14.
[9] 李国勇. 智能控制及其MATLAB 实现[M]北京：电子工业出 版社，2005. 第8组 分工
[1**********] 黄锋：查找相关资料
[1**********] 孟霄：写论文 仿真部分请教的学长
广义预测控制在机组控制系统中的仿真研究
摘要：由于单元机组是一个复杂的系统，一般都存在着很大的延迟和惯性，被控对象具有不确定性和强耦合，常规的PID 控制器在处理大延迟对象上很难获得较好的控制效果。我们将广义预测算法应用到机组控制中，利用MA TLAB 仿真，对响应速度和鲁棒性进行分析，通过比较仿真效果表明广义预测控制方法对系统具有良好的适应性和鲁棒性。 关键词：预测控制 控制系统 仿真
0 引言
广义预测控制（GPC ）是在最小方差控制的基础上，汲取动态矩阵控制（DMC ）和模型算法控制（MAC ）中的多步预测优化策略，随着自适应控制的研究而发展起来的一种预测控制方法。此种控制算法的数学模型是通过推导控制律参数而获取的，在理论研究和实际应用中都得到了重视[1]。本文采用广义预测控制算法对500MW 机组80%负荷点处定压运行进行仿真研究。
1 多变量广义预测控制算法
预测控制算法是在工业生产过程和相关的理论知识结合的基础上提出来的，与传统的PID 控制方法相比较，预测控制算法具有自身的优越性。
广义预测控制具有预测控制算法的基本特征，因而具有优良的控制性能和鲁棒性，使系统的设计更加灵活，被认为是在复杂工业控制领域具有代表性的预测控制算法之一。广义预测控制算法采用的预测模型是具有一定结构和参数的最小化模型，因此对系统输出进行预测时需要已知模型结构，这给预测模型的建立增加了难度，但所需确定的参数较少，减少了计算量。
在GPC 中，采用了最小方差控制中所用的受控自回归积分滑动平均模型来描述收到随机干扰的对象[2]
C (q -1) ξ(t )
（1） A （q ）y (t ) =B (q ) u (t -1) +
∆
-1
-1
其中
A (q -1) =1+a 1q -1+... +a n a q -n a B (q -1) =b 0+b 1q -1+... +b n b q -n b A (q -1) =c 0+c 1q -1+... +c n c q -n c
式中，q 是后移算子，表示后退一个采样周期的相应的量；∆=1-q -1为差分算子；
-1
ξ（t ）是一个不相关的随机序列，表示一类随机噪声的影响。A 、B 、C 都是q 的多项式，
-1
其中多项式B (q ) 的若干首项元素b 0, b 1,... 可以是零，以表示对象相应的时滞数。为了突出方法原理，这里假设C (q -1) =1。
式（1）实际上是用Z 传递函数给出了对象的描述，由输入u 到输出y 间的Z 传递函数为
-1
z -1B (z -1)
（2） G (z ) =
A (z -1)
-1
为了利用模型（1）导出j 步后输出y(t+j)的预测值，首先考虑下述丢翻图（Diophantine ）方程
1=E j (q -1) A ∆+q -1F j (q -1) （3）
其中，E j ，F j 是由A (q ) 和预测长度j 唯一确定的
-1
E j (q -1) =e j . 0+e j . 1q -1+... +e j . j -1q -(j -1) F j (q ) =f j . 0+f j . 1q +... +f j . na q
这样，在式（1）两端乘以E j ∆q 后可得
j -1
-1
-na
E j A ∆y (t +j ) =E j B ∆u (t +j -1) +E j ξ(t +j )
利用式（3），可以写出t+j时刻的输出量
y (t +j ) =E j B ∆u (t +j -1) +F j y (t ) +E j ξ(t +j )
注意到E j ，F j 的形式，可以知道E j B ∆u (t +j -1) 与u (t +j -1) ，u (t +j -2),... 有关，
F j y (t ) 只与y(t)，y(t-1),...有关，而E j ξ(t +j ) 与ξ(t +j ),..., ξ(t +1) 有关[3]。由于在t 时刻
未来的噪声ξ(t +i ), i =1,... j 都是未知的，所以对y(t+j)最合适的预测值可以由下式得到
Λ
y (t +j |t ) =E j B ∆u (t +j -1) +F j y (t ) （4）
这样，根据已知的输入输出信息及未来的的输入值，就可以预测对象未来的输出，式（4）就是GPC 的预测模型。
在GPC 中，t 时刻的优化性能指标具有以下形式：
NU ⎧N 2
2⎫2
min J (t ) =E ⎨∑(y (t +j ) -ω(t +j )) +∑λ(j ) [∆u (t +j -1) ]⎬ （5）
j =1⎭⎩j =N 1
其中，E 为数学期望；ω为对象输出的期望值；N 1和N 2分别为优化时域的识值和终值；NU 为控制时域，即在NU 步后控制量不再变化：u (t +j -1) =u (t +NU -1), j >NU
λ(j ) 为控制加权系数，为简化设计一般常假设其为常数λ。
2 机组受控对象模型
在单元机组协调控制系统中，包括主控系统和子控系统，它们最终的控制对象为锅炉、汽轮机和发电机[4]，它具有大时滞、强耦合、非线性，同时难以平匹配的特点。在整个机组的发电过程中，锅炉和汽轮机的动态特性差别很大，它们对负荷的响应速度一慢一快，所以单元机组内外的两个能量平衡关系互相制约，而它内部的负荷响应性能又与内部参数稳定性存在固有固有的矛盾。因此，在单元机组实施负荷控制的时候，如何协调控制相当重要。为了便于对系统的控制进行分析，在保证满足控制要求的前提下，我们将受控对象的模型简化为两输入两输出的系统，如下图1所示：
图1 单元机组受控对象模型
上图由输入控制和传递函数以及过程输出组成，并且此模型的动态特性由线性常系数多变量传递函数式（6）和（7）来描述[5]：
Y (s ) =G (s ) ⋅U (s ) （6）
⎡y 1(s ) ⎤⎡u 1(s ) ⎤⎡W 11(s ) W 12(s ) ⎤
Y (s ) =U (s ) = 其中⎢y (s ) ⎥，⎢u (s ) ⎥，G (s ) =⎢⎥ （7）
W (s ) W (s ) ⎣2⎦⎣2⎦22⎣21⎦
式中Y (s ) 是被控制矢量，U (s ) 是控制矢量，G (s ) 是传递函数矩阵。
我们选取的动态模型是某厂300MW 直流锅炉燃煤机组在100%负荷点上的数学模型[6]，模型如下[7]：
4. 665s (99s +1) ⎡
⎡N e ⎤⎢582s 2+50s +1(4. 1s +1) ⎢P ⎥=⎢
⎣t ⎦⎢-1. 42(0. 04+0. 96)
⎢70s +1⎣
2. 069(311s +1) ⎤
⎡u t ⎤149s +12(22. 4s +1) ⎥
⎥⋅⎥
1. 265(205s +1) ⎥⎢⎣u b ⎦ （8） (128s +1) 2(11. 7s +1) ⎥⎦
由图1可以知道单元机组复杂的多输入多输出的控制系统简化成了有汽轮机调节门指令u t 和锅炉的燃烧率指令u b 为输入，实发电功率N e 和机前压力P t 为输出的模型。在公式（8）中，u t 和u b 对应于图1 中的u 1和u 2，N e 和P t 对应于图1中的y1(t)和y2(t)。
3 仿真
控制时域长度M 越大，则控制的机动性变强，控制的灵敏度变高，但系统地稳定性和鲁棒性就会下降，为了提高系统的稳定性和鲁棒性，要求M 小，M 越小，远程跟踪控制能力虽有所削弱，但可得到一个稳定的控制[8]，因此，M 的选择应兼顾快速性和稳定性两种，平衡考虑，本文仿真时选取M=6。预测时域P=10，仿真时间为500s ，利用如下程序[9]： clear all
g11=poly2tfd([461.835 4.665 0],conv([58 50 1],[4.1 1]),0,0);
g12=poly2tfd([643.459 2.069],conv(conv([149 1],[149 1]),[22.4 1]),0,0); g22=poly2tfd([259.325 1.265],conv(conv([128 1],[128 1]),[11.7 1]),0,0); g21=poly2tfd([-3.92 -1.42],[70 1],0,0); gd1=poly2tfd(3.8,[14.9 1],0,0); gd2=poly2tfd(4.9,[13.2 1],0,0);
tfinal=500; %Model horizon N delt=10; %sampling period ny=2; %Number of outputs model=tfd2step(tfinal,delt,ny,g11,g21,g12,g22);
plant=model; %No plant/model mismatch dmodel=[]; %Default disturbance model dplant=tfd2step(tfinal,delt,ny,gd1,gd2); P=10;M=6; %Horizons
ywt=[1 1]; uwt=[1 1]; %Q and R
tend=450; %final time for simulation r=[10 0]; %Set-point change inXD XB a=zeros([1,tend]); for i=100:1:tend
a(i)=1; % 3 %step in F at t=50 min. end dstep=[a']; ulim=[];%u limits ylim=[]; tfilter=[]; dplant=[]; dmodel=[]; dstep=[];
[y,u,ym]=cmpc(plant,model,ywt,uwt,M,P,tend,r,ulim,ylim,tfilter,dplant,dmodel,dstep); figure(1) subplot(211) plot(y) grid on legend('Ne','Pt')
ylabel('Ne /% Pt /%') xlabel('k(采样时间10s)') subplot(212) plot(u) legend('Ut','Ub') ylabel('Ub /% Ut /%') xlabel('k(采样时间10s)') grid on
仿真结果如下图2所示：
15N e /% P t /%
1050-5
51015
202530k(采样时间10s)
35404550
15U b /% U t /%
10
5
00
51015
202530k(采样时间10s)
35404550
图2 100%负荷下功率响应曲线和气压响应曲线（上图）、汽机调门指令和锅炉燃烧率指
令响应曲线（下图）
4 结论
随着电力工业的发展，如何将一些先进的控制技术应用于单元机组负荷控制，以提高单元机组负荷控制系统的性能，是一个具有很大应用研究的价值的课题。本文采用了多变量广义预测控制的算法，将其应用于单元机组负荷控制中，进行仿真研究，结果表明该方法在线计算量小，实时性好，参数调整灵活，改善了系统的稳定性和动态特性。
参考文献
[1] 舒迪前. 预测控制系统及其应用[M].北京:机械工业出版社，1996.
[2] 曹楠. 多变量广义预测控制在单元机组协调控制系统中的应用研究[D] 北京:华北电力 大学,2000.
[3] 王伟. 一种广义预测自适应控制的直接方法[J]. 自动化学报 , 1996, 22(3) : 270- 276. [4] 凌呼君, 白雄怀. 火电单元机组机炉协调的广义预测控制策略[J].控制与决策, 2003, 18(增刊 1 ) : 31- 33.
[5] 吴宝琴. 广义预测控制在协调控制系统中的应用[D].保定：华北电力大学，2006. [6] 杜之正. 预测控制在电厂协调控制系统中的研究[D].保定：华北电力大学，2012. [7] 雎刚, 韦红旗, 陈绍炳, 徐治皋. 单元机组负荷多变量模型预测控制[J].中国电机工程学报,2002,04:145-149.
[8] 罗辑，曹建国等. 基于广义预测控制的综合控制算法设计与仿真[J].自动化技术与应用 ， 2004,23（2）：10-14.
[9] 李国勇. 智能控制及其MATLAB 实现[M]北京：电子工业出 版社，2005. 第8组 分工
[1**********] 黄锋：查找相关资料
[1**********] 孟霄：写论文 仿真部分请教的学长