经典预测控制程序

%DMC控制算法
clc;
clear all;
G=input（‘输入传递函数G=‘）;%输入传递函数
%判断是否为稳定系统，若是可以控制，若不是，则无法用DMC算法进行控制
den=G.den{1};%取传函的分母
p=real（roots（den））;%求传函的极点的实部
for i=1:length（p）
       r=p（i）;
   if r>0   %若有某一个极点的实部的实部大于零，则为不稳定系统，DMC无法控制
       pG   %在命令窗口显示极点和传函
Error=（‘您要控制的对象为不稳定系统，DMC算法只适用于稳定系统!‘）
    return
  end
end
%设置DMC参数
Ts=input（‘采样周期 Ts= ‘）;%采样时间
P=input（‘预测时域 P= ‘）;%预测步长
M=input（‘控制时域 M= ‘）;%控制步长
N=80;%截断步长
%设定参考值
yr=10;
%建立系统阶跃响应模型
[y0t0]=step（G0:5:500）;
%初始化DMC
A=zeros（PM）;%动态矩阵
a=zeros（N1）;
for i=1:N
    a（i）=y0（i）;
end
for i=1:P
    for j=1:M
        if i-j+1>0
            A（ij）=a（i-j+1）; %构造矩阵A
        end
    end
end
%初始化向量ys，yue和矩阵A0
ys=ones（N1）;
y=zeros（N1）;
u=zeros（N1）;
e=zeros（N1）;
A0=zeros（PN-1）;
for i=1:P
    for j=N-2:-1:1
        if N-j+1+i-1<=N
            A0（ij）=a（N-j+1+i-1）-a（N-j+i-1）;%构造矩阵A0
        else
            A0（ij）=0;
        end
    end
    A0（iN-1）=a（i+1）;
end
%DMC程序
for k=2:N
    Uk_1=zeros（N-11）;
    for i