bianshi.rar_RLS_模型辨识_过程辨识_递推最小二乘_递推最小二乘法

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');