load V7_1.txt;%阀门开度数据
u1=V7_1;
u=u1(1:150);
load T_1.txt;%温度数据
z1=T_1;
z=z1(1:150);
%d=zeros(10,1);
%RLS递推最小二乘法
c0=[0.001 0.001 0.001 0.001 0.001 0.001 0.001 0.001 0.001 0.001]';
p0=10^6*eye(10,10);
E=0.00000000005;%相对误差E=0.00000000005
c=[c0,zeros(10,14)];
e=zeros(10,15);
for k=6:150
h1=[-z(k-1),-z(k-2),-z(k-3),-z(k-4),-z(k-5),u(k-1),u(k-2),u(k-3),u(k-4),u(k-5)]';
x=h1'*p0*h1+1;
x1=inv(x); % 矩阵求逆
k1=p0*h1*x1;%求出K值 K(k)=p(k-1)h(k)inv(h'(k)p(k-1)h(k)+1)
m1=z(k)-h1'*c0; % c=c(k-1)+K(k)[z(k)-h'(k)c(k-1)]
c1=c0+k1*m1;%辨识参数c
e1=c1-c0; %求参数当前值与上一次的值的差值
e2=e1./c0; %求参数的相对变化
e(:,k)=e2; %把当前相对变化的列向量加入误差矩阵的最后一列
c0=c1; %新获得的参数作为下一次递推的旧参数
c(:,k)=c1; %把递推出的C列向量加入辨识参数矩阵的最后一列
p1=p0-k1*h1'*p0; %求出p(k) p(k)=p(k-1)-K(k)h'(k)p(k-1)
p0=p1;%给下次使用
if e2<=E break;%若参数收敛满足要求，终止计算
end%小循环结束
end%大循环结束
d=[c1];
d%显示被辨识参数
a1=d(1);a2=d(2);a3=d(3);a4=d(4);a5=d(5);b1=d(6);b2=d(7);b3=d(8);b4=d(9);b5=d(10);
%y(k)=-a1*y(k-1)-a2*y(k-2)-a3*y(k-3)-a4*y(k-4)+b1*u(k-1)+b2*u(k-2)+b3*u(k-3)+b4*u(k-4);
%进行所建的模型与实际模型进行误差分析
y(1)=z(1);y(2)=z(2);y(3)=z(3);y(4)=z(4);y(5)=z(5);
Y(1)=z(1);Y(2)=z(2);Y(3)=z(3);Y(4)=z(4);Y(5)=z(5);
for k=6:150
y(k-1)=z(k-1);y(k-2)=z(k-2);y(k-3)=z(k-3);y(k-4)=z(k-4);y(k-5)=z(k-5);
y(k)=-a1*y(k-1)-a2*y(k-2)-a3*y(k-3)-a4*y(k-4)-a5*y(k-5)+b1*u(k-1)+b2*u(k-2)+b3*u(k-3)+b4*u(k-4)+b5*u(k-5);
Y(k)=y(k);
end
Y;
e3=Y'-z;
subplot(1,3,1);
plot(1:150,z,'r');
subplot(1,3,2);
plot(1:150,Y,'g');
subplot(1,3,3);
plot(1:150,e3,'b');
评论0