在现代控制系统设计中，模型预测控制（MPC）因其优异的控制性能而受到广泛关注。本文介绍了一种结合了在线模型拟合和MPC的Simulink子系统，该子系统能够适应ARIMA模型，并通过MPC算法实现对原始系统的控制。
在控制系统设计中，模型预测控制（MPC）是一种先进的控制策略，它通过预测未来系统的行为来优化当前的控制输入。然而，传统的MPC方法通常需要对系统有一个精确的数学模型。在许多实际应用中，系统的精确模型可能难以获得或者过于复杂。因此，数据驱动的控制方法应运而生，它利用在线数据来估计系统的动态行为，从而实现控制目标。
本文介绍的子系统在Simulink环境中实现，它能够同时进行在线模型拟合和MPC控制。首先，子系统考虑了ARIMA模型来描述原始系统的输出和输入之间的关系：
y(t) = A1*y(t-1) + A2*y(t-2) + ... + Ap*y(t-p) + B1*u(t-1) + B2*u(t-2) + ... + Bq*u(t-q)
其中，y(t)和u(t)分别代表原始系统的输出和输入。参数p和q分别代表ARIMA模型的阶数。这些参数可以通过修改子系统的参数来手动调整。子系统将上述模型拟合到原始系统的输入输出数据上，即使原始系统的内部特性完全未知。
在第二步中，子系统将适当的MPC应用于获得的ARIMA模型。控制器尝试通过合理设计控制序列，使系统输出y(t)达到期望的参考信号yd(t)。控制器参数，包括预测和控制范围，分别用Ny和Nu表示。控制器预测系统输出y(t)在未来的步骤中，如下所示：
y(t+1) = A1*y(t) + A2*y(t-1) + ... + Ap*y(t-p+1) + B1*U(t) + B2*u(t-1) + ... + Bq*u(t-q+1)
控制器精确地确定控制序列U=[U(t) U(t+1) ... U(t+Nu-1)]，使得预测输出Y=[y(t+1) y(t+2) ... y(t+Ny)]接近参考信号yd(t)。为此，定义并解决了一个优化问题，以最优地确定控制序列U。在每一步t，获得最优控制序列，其第一个条目作为控制输入应用于原始系统。
在本节中，通过一些示例来说明如何使用开发的子系统。
在第一个示例中，假设希望控制LTI系统1/(s+1)以跟随正弦波参考信号。为此，原始系统（即1/(s+1)）被实现，其输出连接到数据驱动MPC的"y"端口，数据驱动MPC的"uc"端口连接到其输入。期望的参考信号被视为正弦波，其参数可以修改。然后，控制器尝试将原始系统的输出信号关闭到期望的信号。
通过单击数据驱动MPC，显示以下属性：
其中一些参数是预先定义的，其他参数在这里定义。参数"T"代表时间步长。参数umin、umax和udmax分别代表控制输入信号的最小值、最大值和时间导数最大值。所有这些参数都可以根据原始系统的限制进行更改。
开发的子系统利用一个名为myfun的脚本，该脚本模拟在线估计和MPC设计。以下是代码：
// Matlab function
function Z=myfun(y,yref,u,t,umin,umax,udmax,p,q,my,mu,Ny,Nu,T,th0,P0)
global U Y P th Uo yh;
U=[U(2:end,:);u'];
if t==0
P=P0;th=th0;
U=zeros(q,mu);Y=zeros(p,my);
Uo=zeros(Nu,mu);
end
H=zeros(my,my^2*p+my*mu*q);
ind=1;
for i=1:p
H(:,ind:ind+my^2-1)=kron(Y(p-i+1,:),eye(my));ind=ind+my^2;
end
for i=1:q
H(:,ind:ind+my*mu-1)=kron(U(q-i+1,:),eye(my));ind=ind+my*mu;
end
H=H'
;
m=size(H,
2
);
L=(P*H)/(eye(m)+H
'
*P*H);
P=P-L*H'
*P;
th=th+L*(y
'
-H'
*th);
yh=zeros(my,
1
);
for
i=
1
:p
ind=1+(i-1)*my^2;
Ai=reshape(th(ind:ind+my*my-1),my,my);
yh=yh+Ai*Y(p-i+1,:)
'
;
end
for i=1:q
ind=1+p*my^2+(i-1)*my*mu;
Bi=reshape(th(ind:ind+my*mu-1),my,mu);
yh=yh+Bi*U(q-i+1,:)'
;
end
Y=[Y(
2
:end,:);y
'
];
umin=reshape(umin,length(umin),1);
umax=reshape(umax,length(umax),1);
udmax=reshape(udmax,length(udmax),1);
G=zeros(Ny*my,Nu*mu);g=zeros(Ny*my,1);
for i=1:Ny
Ge=zeros(my,size(G,2));ge=zeros(my,1);
for j=1:p
ind=1+(j-1)*my^2;
Aj=reshape(th(ind:ind+my*my-1),my,my);
if j>=i
ge=ge+Aj*Y(end-j+i,:)'
;
else
Gw=G((i-j-1)*my+1:(i-j)*my,:);
gw=g((i-j-1)*my+1:(i-j)*my,:);
Ge=Ge+Aj*Gw;
ge=ge+Aj*gw;
end
end
for
j=
1
:q
ind=1+p*my^2+(j-1)*my*mu;
Bj=reshape(th(ind:ind+my*mu-1),my,mu);
if
j>i
ge=ge+Bj*U(end-j+i+1,:)
'
;
elseif i-j+1<=Nu
Ge(:,(i-j)*mu+1:(i-j+1)*mu)=Ge(:,(i-j)*mu+1:(i-j+1)*mu)+Bj;
end
end
G((i-1)*my+1:i*my,:)=Ge;
g((i-1)*my+1:i*my,:)=ge;
end
z1=ones(Nu,1);
if ~isempty(udmax)
A1=eye(Nu*mu);b1=Uo+T*kron(z1,udmax);
A2=-eye(Nu*mu);b2=-Uo+T*kron(z1,udmax);
A=[A1;A2];b=[b1;b2];
else
A=[];b=[];
end
if ~isempty(umin)
Umin=kron(z1,umin);
else
Umin=[];
end
if ~isempty(umax)
Umax=kron(z1,umax);
else
Umax=[];
end
Yref=kron(ones(Ny,1),yref);
gb=g-Yref;
opt=optimset('
disp
'
,'
none
'
);
Uo=quadprog(2*(G'
*G),2*G
'
*gb,A,b,[],[],Umin,Umax,Uo,opt);
uc=Uo(1:mu);
Z=[uc;yh];
end
//
上述代码在开发的子系统中内部使用，解释了如何实现数据驱动的MPC。