gpc预测控制CARIMA模型

%产生周期为200时间为T，幅值为l的方波信号的给定值
T=600;
[yr0t]=gensig（‘square‘200T1）;
%plot（tyr0）

nn=input（‘时域长度 nn=‘）;
n=input（‘预测长度 n=‘）;
m=input（‘控制长度 m=‘）;
lmd=input（‘控制加权系数 lmd=‘）;
a=input（‘柔化系数 a=‘）;
P=（1e+5）*eye（n+1）;
%参数初始值
uuu=0;yyy=0;
uu=zeros（n1）;u=zeros（m1）;
yy=zeros（n1）;y1=zeros（n1）;
Q=zeros（n+11）;
Q（11）=1;
Q（n+11）=1;
t1=1;
nm=length（t）;
for ij=2:nm
    yr=yr0（ij）;
%根据系统模型，计算k时刻的输出值y（k）
%     A=yy（n-2）*yy（n-1）;
%     B=1.2*yy（n-2）^2+0.8*yy（n-1）^2+1;
%     y=A./B+0.5*uu（n1）;
y=0.1*（2.571*yy（n-11）-2.806*yy（n-21）+2.04*yy（n-31）-1.095*yy（n-41）+...
   0.2818*yy（n-51）+0.001283+0.0243*uu（n1）+0.0473*uu（n-11）+...
    0.01576*uu（n-21）+0.005498*uu（n-31））;
    %产生均匀分布的白噪声
    %保存k时刻及以前的n个输出值y（k）y（k-1）...y（k-n）
    for i=1:n-1
        yy（i1）=yy（i+11）;
    end
    yy（n1）=y;
    yyy=[yyy;y]; %保存各k时刻的n个输出量以便绘图
    %根据最小二乘公式．由y（k）计算G阵的各元素值g0g0...gn
    for i=1:n
        X（1i）=uu（i1）;
    end
    %     a9=0;
%     for i=1:1
%         a9=a9+rand;
%     end
%     a8=0.01*（a9-6）;

    X（1n+1）=1;
    K=P*X‘*inv（t1+X*P*X‘）;
    P=（eye（n+1