%Clarke廣義預測控制（C=1）（對象參數已知）?
%N1=d、N、Nu取不同的值?
clear?all;?close?all;?
?
a=[1?-2?1.1];
b=[1?2];?
c=1;?
d=4;?%對象參數?
na=length（a）-1;?
b=[zeros（1d-1）?b];?
nb=length（b）-1;?%na、nb為多項式A、B階次（因d!=1，對b添0）?
aa=conv（a[1?-1]）;?
naa=na+1;?%aa的階次?
?
N1=d;?N=8;?Nu=5;?%最小輸出長度、預測長度、控制長度?
gamma=1*eye（Nu）;?alpha=0.7;?%控制加權矩陣、輸出柔化系數?
?
L=400;?%控制步數?
uk=zeros（d+nb1）;?%輸入初值：uk（i）表示u（k-i）?
duk=zeros（d+nb1）;?%控制增量初值?
yk=zeros（naa1）;?%輸出初值?
%w=10*[ones（L/41）;-ones（L/41）;ones（L/41）;-ones（L/4+d1）];?%設定值?
w=10*[zeros（L/21）;ones（L/2+d1）];?
xi=sqrt（0.01）*randn（L1）;?%白噪聲序列?
?
%求解多步Diophantine方程并構建F1、F2、G?
[EFG]=multidiophantine（aabcN）;?
G=G（N1:N:）;?
F1=zeros（N-N1+1Nu）;?F2=zeros（N-N1+1nb）;?
for?i=1:N-N1+1?
????for?j=1:min（iNu）;????
????????F1（ij）=F（i+N1-1i+N1-1-j+1）;??
????end?
????for?j=1:nb;???????????
????????F2（ij）=F（i+N1-1i+N1-1+j）;??
????end?
end?
?
for?k=1:L?
????time（k）=k;?
????y（k）=-aa（2:naa+1）*yk+b*duk（1:nb+1）+xi（k）;?%采集輸出數據?
????Yk=[y（k）;?yk（1:na）];?%構建向量Y（k）?
????dUk=duk（1:nb）;?%構建向量ΔU（k-j）?
?????
????%參考軌跡?
????yr（k）=y（k）;?
????for?i=1:N?
????????yr（k+i）=alpha*yr（k+i-1）+（1-alpha）*w（k+d）;?
????end?
????Yr=[yr（k+N1:k+N）]‘;?%構建向量Yr（k）?
????????????
????%求控制量?
????dU=inv（F1‘*F1+gamma）*F1‘*（Yr-F2*dUk-G*Yk）;?%ΔU?
????du（k）=dU（1）;?
????u（k）=uk（1）+du（k）;?
?????
????%更新數據?
????for?i=1+nb:-1:2?
????????uk（i）=uk（i-1）;?
????????duk（i）=duk（i-1）;?
????end?
????uk（1）=u（k）;?
????duk（1）=du（k）;?
?????
????for?i=naa:-1:2?
????????yk（i）=yk（i-1）;?
????end?
????yk（1）=y（k）;?
end?
subplot（211）;?
plot（timew（1:L）‘r:‘timey）;?
xlabel（‘k‘）;?ylabel（‘w（k）、y（k）‘）;?
legend（‘w（k）‘‘y（k）‘）;?
subplot（212）;?
plot（timeu）;?
xlabel（‘k‘）;?ylabel（‘u（k）‘）;

?屬性????????????大小?????日期????時間???名稱
-----------?---------??----------?-----??----

?????文件???????1876??2010-05-05?10:41??GPC\GPC1.m

?????文件????????922??2010-04-27?12:36??GPC\multidiophantine.m

?????目錄??????????0??2010-05-05?11:40??GPC

-----------?---------??----------?-----??----

?????????????????2798????????????????????3