📄 mpc.txt
字号:
程序一：跟踪阶跃函数
clear all;
close all;
Ts=0.5;%采样周期0.5s
Tr=1;H=10;%参考轨迹时间常数为1,预测时域优化长度为10
Km=3;Tm=3;%预测模型的参数
beta=exp(-Ts/Tr);
alph=exp(-Ts/Tm);
%对象离散化
sys=tf(3,[3,1]);
dsys=c2d(sys,Ts,'z');
[num,den]=tfdata(dsys,'v');
%预测模型离散化
sysm=tf(Km,[Tm,1]);
dsysm=c2d(sysm,Ts,'z');
[numm,denm]=tfdata(dsysm,'v');
u_1=0;
y_1=0;
ym_1=0;
c=1;%设定值
for k=1:1:100
time(k)=k*Ts;
%经z变换后的离散化对象
y(k)=-den(2)*y_1+num(2)*u_1;
ym(k)=-denm(2)*ym_1+numm(2)*u_1;
%修正后的过程输出值ypav(k)=y(k)
ypav(k)=y(k)
%控制器输出
u(k)=((c-ypav(k))*(1-beta^H)+ym(k)*(1-alph^H))/(Km*(1-alph^H));
%参数更新
u_1=u(k);
y_1=y(k);
ym_1=ym(k);
end
figure(1);
plot(time,c,'b',time,y,'r');
axis([0 50 0 1.2])
xlabel('time(s)');ylabel('c,y');
程序二：跟踪阶跃加扰动后响应
clear all;
close all;
Ts=0.5;%采样周期0.5s
Tr=1;H=10;%参考轨迹时间常数为1,预测时域优化长度为10
Km=3;Tm=3;%预测模型的参数
beta=exp(-Ts/Tr);
alph=exp(-Ts/Tm);
%对象离散化
sys=tf(3,[3,1]);
dsys=c2d(sys,Ts,'z');
[num,den]=tfdata(dsys,'v');
%预测模型离散化
sysm=tf(Km,[Tm,1]);
dsysm=c2d(sysm,Ts,'z');
[numm,denm]=tfdata(dsysm,'v');
u_1=0;
y_1=0;
ym_1=0;
c=1;%设定值
for k=1:1:40
time(k)=k*Ts;
%经z变换后的离散化对象
y(k)=-den(2)*y_1+num(2)*u_1;
ym(k)=-denm(2)*ym_1+numm(2)*u_1;
%控制器输出
u(k)=((c-y(k))*(1-beta^H)+ym(k)*(1-alph^H))/(Km*(1-alph^H));
%参数更新
u_1=u(k);
y_1=y(k);
ym_1=ym(k);
end
for k=41:1:44
time(k)=k*Ts;
c=1.2;
%经z变换后的离散化对象
y(k)=-den(2)*y_1+num(2)*u_1;
ym(k)=-denm(2)*ym_1+numm(2)*u_1;
%控制器输出
u(k)=((c-y(k))*(1-beta^H)+ym(k)*(1-alph^H))/(Km*(1-alph^H));
%参数更新
u_1=u(k);
y_1=y(k);
ym_1=ym(k);
end
for k=45:1:140
time(k)=k*Ts;
c=1;
%经z变换后的离散化对象
y(k)=-den(2)*y_1+num(2)*u_1;
ym(k)=-denm(2)*ym_1+numm(2)*u_1;
%控制器输出
u(k)=((c-y(k))*(1-beta^H)+ym(k)*(1-alph^H))/(Km*(1-alph^H));
%参数更新
u_1=u(k);
y_1=y(k);
ym_1=ym(k);
end
figure(1);
plot(time,y,'r');
axis([0 70 0 1.2]);
xlabel('time(s)');ylabel('c,y');
程序三：跟踪斜坡函数
clear all;
close all;
Ts=0.5;%采样周期0.5s
Tr=1;H1=4;H2=5;%参考轨迹时间常数为1,预测步长H1为4，H2为5
Km=15;Tm=15;%预测模型的参数
beta=exp(-Ts/Tr);
alph=exp(-Ts/Tm);
%对象离散化
sys=tf(5,[5,1]);
dsys=c2d(sys,Ts,'z');
[num,den]=tfdata(dsys,'v');
%预测模型离散化
sysm=tf(Km,[Tm,1]);
dsysm=c2d(sysm,Ts,'z');
[numm,denm]=tfdata(dsysm,'v');
u_1=0;
y_1=0;
ym_1=0;
s1=0;
s2=0;
G1(H1)=Km*(1-alph^H1);
G1(H2)=Km*(1-alph^H2);
for i=1:1:H1-1
s1=s1+i*alph^(H1-1-i);
end;
for j=1:1:H2-1
s2=s2+i*alph^(H2-1-i);
end;
G2(H1)=Km*(1-alph)*s1;
G2(H2)=Km*(1-alph)*s2;
for k=1:1:200
time(k)=k*Ts;
c(k)=time(k);%设定值
%c(k+H1)=(k+H1)*Ts;
%c(k+H2)=(k+H2)*Ts;
%经z变换后的离散化对象
y(k)=-den(2)*y_1+num(2)*u_1;
ym(k)=-denm(2)*ym_1+numm(2)*u_1;
%修正后的过程输出值ypav(k)=y(k)
%ypav(k)=y(k);
%yr(k+H1)=e(k+H1)-beta^H1*(c(k)-ypav(k))其中的e(k+H1)应该为c(k+H1)
%即yr(k+H1)=c(k+H1)-beta^H1*(c(k)-ypav(k))
%X1=yr(k+H1)-alph^H1*ym(k)-e(k+H1)
yr(k+H1)=(k+H1)*Ts-beta^H1*((k+0)*Ts-y(k));
yr(k+H2)=(k+H2)*Ts-beta^H2*((k+0)*Ts-y(k));
e(k+H1)=y(k)-ym(k);
e(k+H2)=y(k)-ym(k);
X1=yr(k+H1)-alph^H1*ym(k)-e(k+H1);
X2=yr(k+H2)-alph^H2*ym(k)-e(k+H2);
%控制器输出
u(k)=(G2(H2)*X1-G2(H1)*X2)/(G1(H1)*G2(H2)-G1(H2)*G2(H1));
%参数更新
u_1=u(k);
y_1=y(k);
ym_1=ym(k);
end
figure(1);
plot(time,c,'b-',time,y,'r--');
%axis([0 50 0 1.2])
xlabel('time(s)');ylabel('c,y');
程序四：跟踪设定曲线
clear all;
close all;
Ts=0.5;%采样周期0.5s
Tr=1;H1=4;H2=5;%参考轨迹时间常数为1,预测步长H1为4，H2为5
Km=15;Tm=15;%预测模型的参数
beta=exp(-Ts/Tr);
alph=exp(-Ts/Tm);
%对象离散化
sys=tf(5,[5,1]);
dsys=c2d(sys,Ts,'z');
[num,den]=tfdata(dsys,'v');
%预测模型离散化
sysm=tf(Km,[Tm,1]);
dsysm=c2d(sysm,Ts,'z');
[numm,denm]=tfdata(dsysm,'v');
u_1=0;
y_1=0;
ym_1=0;
s1=0;
s2=0;
G1(H1)=Km*(1-alph^H1);
G1(H2)=Km*(1-alph^H2);
for i=1:1:H1-1
s1=s1+i*alph^(H1-1-i);
end;
for j=1:1:H2-1
s2=s2+i*alph^(H2-1-i);
end;
G2(H1)=Km*(1-alph)*s1;
G2(H2)=Km*(1-alph)*s2;
for k=1:1:80
time(k)=k*Ts;
c(k)=time(k);%设定值
%c(k+H1)=(k+H1)*Ts;
%c(k+H2)=(k+H2)*Ts;
%经z变换后的离散化对象
y(k)=-den(2)*y_1+num(2)*u_1;
ym(k)=-denm(2)*ym_1+numm(2)*u_1;
%修正后的过程输出值ypav(k)=y(k)
%ypav(k)=y(k);
%yr(k+H1)=e(k+H1)-beta^H1*(c(k)-ypav(k))其中的e(k+H1)应该为c(k+H1)
%即yr(k+H1)=c(k+H1)-beta^H1*(c(k)-ypav(k))
%X1=yr(k+H1)-alph^H1*ym(k)-e(k+H1)
yr(k+H1)=(k+H1)*Ts-beta^H1*((k+0)*Ts-y(k));
yr(k+H2)=(k+H2)*Ts-beta^H2*((k+0)*Ts-y(k));
e(k+H1)=y(k)-ym(k);
e(k+H2)=y(k)-ym(k);
X1=yr(k+H1)-alph^H1*ym(k)-e(k+H1);
X2=yr(k+H2)-alph^H2*ym(k)-e(k+H2);
%控制器输出
u(k)=(G2(H2)*X1-G2(H1)*X2)/(G1(H1)*G2(H2)-G1(H2)*G2(H1));
%参数更新
u_1=u(k);
y_1=y(k);
ym_1=ym(k);
end
for k=81:1:160
time(k)=k*Ts;
c(k)=40;%设定值
y(k)=-den(2)*y_1+num(2)*u_1;
ym(k)=-denm(2)*ym_1+numm(2)*u_1;
yr(k+H1)=40-beta^H1*(40-y(k));
yr(k+H2)=40-beta^H2*(40-y(k));
e(k+H1)=y(k)-ym(k);
e(k+H2)=y(k)-ym(k);
X1=yr(k+H1)-alph^H1*ym(k)-e(k+H1);
X2=yr(k+H2)-alph^H2*ym(k)-e(k+H2);
u(k)=(G2(H2)*X1-G2(H1)*X2)/(G1(H1)*G2(H2)-G1(H2)*G2(H1));
u_1=u(k);
y_1=y(k);
ym_1=ym(k);
end
for k=161:1:200
time(k)=k*Ts;
c(k)=-1.5*k*Ts+160;
y(k)=-den(2)*y_1+num(2)*u_1;
ym(k)=-denm(2)*ym_1+numm(2)*u_1;
yr(k+H1)=-1.5*(k+H1)*Ts+160-beta^H1*(c(k)-y(k));
yr(k+H2)=-1.5*(k+H2)*Ts+160-beta^H2*(c(k)-y(k));
e(k+H1)=y(k)-ym(k);
e(k+H2)=y(k)-ym(k);
X1=yr(k+H1)-alph^H1*ym(k)-e(k+H1);
X2=yr(k+H2)-alph^H2*ym(k)-e(k+H2);
u(k)=(G2(H2)*X1-G2(H1)*X2)/(G1(H1)*G2(H2)-G1(H2)*G2(H1));
u_1=u(k);
y_1=y(k);
ym_1=ym(k);
end
for k=200:1:300
time(k)=k*Ts;
c(k)=10;
y(k)=-den(2)*y_1+num(2)*u_1;
ym(k)=-denm(2)*ym_1+numm(2)*u_1;
yr(k+H1)=10-beta^H1*(c(k)-y(k));
yr(k+H2)=10-beta^H2*(c(k)-y(k));
e(k+H1)=y(k)-ym(k);
e(k+H2)=y(k)-ym(k);
X1=yr(k+H1)-alph^H1*ym(k)-e(k+H1);
X2=yr(k+H2)-alph^H2*ym(k)-e(k+H2);
u(k)=(G2(H2)*X1-G2(H1)*X2)/(G1(H1)*G2(H2)-G1(H2)*G2(H1));
u_1=u(k);
y_1=y(k);
ym_1=ym(k);
end
figure(1);
plot(time,c,'b-',time,y,'r--');
axis([0 150 0 45])
xlabel('time(s)');ylabel('c,y');
⌨️ 快捷键说明
复制代码
Ctrl + C
搜索代码
Ctrl + F
全屏模式
F11
切换主题
Ctrl + Shift + D
显示快捷键
?
增大字号
Ctrl + =
减小字号
Ctrl + -