%DMC控制算法
clc;
clear?all;
G=input（‘輸入傳遞函數(shù)G=‘）;%輸入傳遞函數(shù)
%判斷是否為穩(wěn)定系統(tǒng)，若是可以控制，若不是，則無法用DMC算法進(jìn)行控制
den=G.den{1};%取傳函的分母
p=real（roots（den））;%求傳函的極點(diǎn)的實(shí)部
for?i=1:length（p）
???????r=p（i）;
???if?r>0???%若有某一個(gè)極點(diǎn)的實(shí)部的實(shí)部大于零，則為不穩(wěn)定系統(tǒng)，DMC無法控制
???????pG???%在命令窗口顯示極點(diǎn)和傳函
Error=（‘您要控制的對(duì)象為不穩(wěn)定系統(tǒng)，DMC算法只適用于穩(wěn)定系統(tǒng)!‘）
????return
??end
end
%設(shè)置DMC參數(shù)
Ts=input（‘采樣周期?Ts=?‘）;%采樣時(shí)間
P=input（‘預(yù)測時(shí)域?P=?‘）;%預(yù)測步長
M=input（‘控制時(shí)域?M=?‘）;%控制步長
N=80;%截?cái)嗖介L
%設(shè)定參考值
yr=10;
%建立系統(tǒng)階躍響應(yīng)模型
[y0t0]=step（G0:5:500）;
%初始化DMC
A=zeros（PM）;%動(dòng)態(tài)矩陣
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）;?%構(gòu)造矩陣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）;%構(gòu)造矩陣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