系統(tǒng)辨識與自適應(yīng)控制作業(yè)

1、精選優(yōu)質(zhì)文檔-傾情為你奉上 系統(tǒng)辨識與自適應(yīng)控制學(xué) 院： 專 業(yè)： 學(xué) 號： 姓 名： 系統(tǒng)辨識與自適應(yīng)控制作業(yè)一、 對時變系統(tǒng)進(jìn)行參數(shù)估計。系統(tǒng)方程為：y(k)+a(k)y(k-1)=b(k)u(k-1)+e(k) 其中：e(k)為零均值噪聲，a(k)= b(k)=要求：1對定常系統(tǒng)（a=0.8,b=0.5）進(jìn)行結(jié)構(gòu)（階數(shù)）確定和參數(shù)估計； 2對時變系統(tǒng)，取不同值（0.90.99）時對系統(tǒng)辨識結(jié)果和過程進(jìn)行比較、討論 3對辨識結(jié)果必須進(jìn)行殘差檢驗解：一（1）：分析：采用最小二乘法（LS）：最小二乘的思想就是尋找一個的估計值，使得各次測量的與由估計確定的量測估計之差的平方和最小，由于此方法兼顧

2、了所有方程的近似程度，使整體誤差達(dá)到最小，因而對抑制誤差是有利的。在此，我應(yīng)用批處理最小二乘法，收斂較快，易于理解，在系統(tǒng)參數(shù)估計應(yīng)用中十分廣泛。作業(yè)程序： clear all;a=1 0.8' b= 0.5' d=3; %對象參數(shù)na=length(a)-1; nb=length(b)-1; %na、nb為A、B階次L=500; %數(shù)據(jù)長度uk=zeros(d+nb,1); %輸入初值：uk(i)表示u(k-i)yk=zeros(na,1); %輸出初值x1=1; x2=1; x3=1; x4=0; S=1; %移位寄存器初值、方波初值xi=randn(L,1); %白噪聲序

3、列theta=a(2:na+1);b; %對象參數(shù)真值for k=1:L phi(k,:)=-yk;uk(d:d+nb)' %此處phi(k,:)為行向量，便于組成phi矩陣 y(k)=phi(k,:)*theta+xi(k); %采集輸出數(shù)據(jù) IM=xor(S,x4); %產(chǎn)生逆M序列 if IM=0 u(k)=-1; else u(k)=1; end S=not(S); M=xor(x3,x4); %產(chǎn)生M序列 %更新數(shù)據(jù) x4=x3; x3=x2; x2=x1; x1=M; for i=d+nb:-1:2 uk(i)=uk(i-1); end uk(1)=u(k); for i=

4、na:-1:2 yk(i)=yk(i-1); end yk(1)=y(k);endthetae=inv(phi'*phi)*phi'*y' %計算參數(shù)估計值thetae結(jié)果：thetae =0.7787 ，0.5103 真值=0.8，0.5解：一（2）：采用遺忘因子遞推最小二乘參數(shù)估計；其仿真算法如下：Step1：設(shè)置初值、，及遺忘因子，輸入初始數(shù)據(jù)；Step2：采樣當(dāng)前輸入和輸出數(shù)據(jù)；Step3：利用含有遺忘因子的遞推公式計算、和；Step4：k=k+1，返回Step2繼續(xù)循環(huán)。其中：對時變系統(tǒng)，取不同值（0.90.99）時對系統(tǒng)辨識結(jié)果和過程進(jìn)行比較、討論作業(yè)程序：

5、clear all; close all;a=1 0.8;b=0.5;d=2; %對象參數(shù)na=length(a)-1;nb=length(b)-1; %na、nb為A、B階數(shù)L=500; %仿真長度uk=zeros(d+nb,1); %輸入初值yk=zeros(na,1); %輸出初值u=randn(L,1); %輸入采用白噪聲序列xi=sqrt(0.1)*randn(L,1); %白噪聲序列 thetae_1=zeros(na+nb+1,1);p=106*eye(na+nb+1);lambda=0.95; %遺忘因子for k=1:L; if k>300; a=1 0.6'b

6、=0.3' endthetae(:,k)=a(2:na+1);b; phi=-yk;uk(d:d+nb); y(k)=phi'*thetae(:,k)+xi(k); %采集輸出數(shù)據(jù) %遞推最小二乘法 K=p*phi/(lambda+phi'*p*phi); thetae(:,k)=thetae_1+K*(y(k)-phi'*thetae_1); p=(eye(na+nb+1)-K*phi')*p/lambda;    %更新數(shù)據(jù) thetae_1=thetae(:,k); for i=d+nb:-1:2 uk(i)=uk(i-1);

7、end uk(1)=u(k); for i=na:-1:2 yk(i)=yk(i-1); end yk(1)=y(k);endsubplot(1,2,1)plot(1:L,thetae(1:na,:);hold on; plot(1:L, thetae(1:na,:),'k:');xlabel('k');ylabel('參數(shù)估計a');axis(0 L -0.5 2);subplot(1,2,2)plot(1:L,thetae(na+1:na+nb+1,:);hold on;plot(1:L,thetae(na+1:na+nb+1,:),'

8、;k:');xlabel('k');ylabel('參數(shù)估計b');axis(0 L -0.5 2);仿真結(jié)果圖1-1遺忘因子遞推最小二乘法的參數(shù)估計結(jié)果（=0.95）圖1-2遺忘因子遞推最小二乘法的參數(shù)估計結(jié)果（=0.99）綜上所述，可知從仿真圖可以看出，當(dāng)0<k<300時，a的值大約在0.8左右波動，b的值大約在0.5左右波動。當(dāng)k>300時，a、b的值發(fā)生變化，此時a的值大約在0.6左右波動，b的值大約在0.3左右波動。所以在此種情況下，可以達(dá)到辨識目的，但是由于波動較大，辨識效果不理想。同時，當(dāng)=0.99時比前面所述=0.96時

9、的辨識結(jié)果更為明顯。這種情況下當(dāng)估計值變化穩(wěn)定時，估計值的變化率比較小，辨識效果較為理想。一（3）：根據(jù)殘差平方和估計模型的階次SISO過程的差分方程模型的輸出殘差為，數(shù)據(jù)長度L，為階時的數(shù)據(jù)矩陣，為階時的參數(shù)的估計量，為模型階次估計值，為真實階次，則殘差平方和函數(shù)：殘差平方和有這樣的性質(zhì)：當(dāng)L足夠大時，隨著增加先是顯著地下降，當(dāng)>時，值顯著下降的現(xiàn)象就終止。這就是損失函數(shù)法來定階的原理。程序：N=15;%4位移位寄存器產(chǎn)生的M序列的周期y1=1;y2=1;y3=1;y4=0; for i=1:N; x1=xor(y3,y4); x2=y1; x3=y2; x4=y3; y(i)=y4;

10、 if y(i)>0.5,u(i)=-1; else u(i)=1; end y1=x1;y2=x2;y3=x3;y4=x4; end%白噪聲的產(chǎn)生A=19;x0=12;M=1024;for k=1:N x=A*x0; x1=mod(x,M); v(k)=x1/512; x0=x1;end%辨識主程序z=zeros(7,N);zs=zeros(7,N);zm=zeros(7,N);zmd=zeros(7,N);z(1)=0;z(2)=0;zs(1)=0; zs(2)=0; zm(1)=0; zm(2)=0; zmd(1)=0;zmd(2)=0;theta=zeros(7,1); p0=1

11、06*eye(7,7);the=theta,zeros(7,14); e=zeros(7,15);for k=3:N; h=-z(k-1),-z(k-2),u(k-1),u(k-2),v(k),v(k-1),v(k-2)' x=inv(h'*p0*h+1); q=p0*h*x; d1=z(k)-h'*theta; theta1=theta+q*d1; zs(k)=-1.5*z(k-1)+0.7*z(k-2)+u(k-1)+0.5*u(k-2); zm(k)=-z(k-1),-z(k-2),u(k-1),u(k-2)*theta1(1);theta1(2);theta1(

12、3);theta1(4); zmd(k)=h'*theta1; e(:,k)=theta1-theta; theta=theta1; the(:,k)=theta1; p1=p0-q*q'*h'*p0*h+1; p0=p1; endfigure(1); i=1:N;plot(i,the(1,:),'r',i,the(2,:),'r:',i,the(3,:),'b',i,the(4,:),'b:',i,the(5,:),'g',i,the(6,:),'g:',i,the(7,:

13、),'g+')title('辨識曲線')figure(2); i=1:N;plot(i,e(1,:),'r',i,e(2,:),'r:',i,e(3,:),'b',i,e(4,:),'b:',i,e(5,:),'g',i,e(6,:),'g:',i,e(7,:),'r+') title('辨識誤差曲線')圖1-3-1辨識曲線圖1-3-2辨識誤差曲線由以上對比曲線可以看出，經(jīng)過最小二乘法估計得到的數(shù)據(jù)與實際數(shù)據(jù)之間雖然存在區(qū)別，但是基本符

14、合要求, 說明辨識效果較好。二、 已知系統(tǒng)方程：y(k)+1.2 y(k-1)+0.35 y(k-2)=u(k-2)+0.4 u(k-3)+e(k)+1.1e(k-1)+0.28e(k-2)其中e(k)為白噪聲序列。要求：1、對系統(tǒng)進(jìn)行參數(shù)遞推估計； 2、進(jìn)行最小方差控制器設(shè)計； 3、進(jìn)行基于最小方差控制的自適應(yīng)控制系統(tǒng)仿真分析（設(shè)輸入信號為方波信號）。二（1）：程序：clear all; close all;a=1 1.2 0.35' b=1 0.4' c=1 1.1 0.28' d=3; %對象參數(shù)na=length(a)-1; nb=length(b)-1; nc

15、=length(c)-1; %na、nb、nc為A、B、C階次L=1000; %仿真長度uk=zeros(d+nb,1); %輸入初值：uk(i)表示u(k-i)yk=zeros(na,1); %輸出初值xik=zeros(nc,1); %噪聲初值xiek=zeros(nc,1); %噪聲估計初值u=randn(L,1); %輸入采用白噪聲序列xi=sqrt(0.1)*randn(L,1); %白噪聲序列theta=a(2:na+1);b;c(2:nc+1); %對象參數(shù)thetae_1=zeros(na+nb+1+nc,1); %na+nb+1+nc為辨識參數(shù)個數(shù)P=106*eye(na+n

16、b+1+nc);for k=1:L phi=-yk;uk(d:d+nb);xik; y(k)=phi'*theta+xi(k); %采集輸出數(shù)據(jù) phie=-yk;uk(d:d+nb);xiek; %組建phie %遞推增廣最小二乘法 K=P*phie/(1+phie'*P*phie); thetae(:,k)=thetae_1+K*(y(k)-phie'*thetae_1); P=(eye(na+nb+1+nc)-K*phie')*P; xie=y(k)-phie'*thetae(:,k); %白噪聲的估計值 %更新數(shù)據(jù) thetae_1=thetae

17、(:,k); for i=d+nb:-1:2 uk(i)=uk(i-1); end uk(1)=u(k); for i=na:-1:2 yk(i)=yk(i-1); end yk(1)=y(k); for i=nc:-1:2 xik(i)=xik(i-1); xiek(i)=xiek(i-1); end xik(1)=xi(k); xiek(1)=xie;endfigure(1)plot(1:L,thetae(1:na,:);xlabel('k'); ylabel('參數(shù)估計a');legend('a_1','a_2'); axi

18、s(0 L -2 2);figure(2)plot(1:L,thetae(na+1:na+nb+1,:);xlabel('k'); ylabel('參數(shù)估計b');legend('b_0','b_1'); axis(0 L -2 2);figure(3)plot(1:L,thetae(na+nb+2:na+nb+nc+1,:);xlabel('k'); ylabel('參數(shù)估計c');legend('c_1','c_2'); axis(0 L -2 2);對系統(tǒng)進(jìn)行最

19、廣最小二乘的遞推估計（RELS），仿真結(jié)果如下圖2.1所示。當(dāng)仿真步數(shù)K大于50時，參數(shù)估計值就趨于穩(wěn)定。參數(shù) a1、a2、b1、b2、d1、d2、d3 估計結(jié)果：ans = 1.2060 0.3548 1.0000 0.4060 1.0000 1.1060 0.2843圖2-1系統(tǒng)辨識曲線二（2）：（1） 對系統(tǒng)施加最小方差控制。由于系統(tǒng)已知即：系統(tǒng)階次和系統(tǒng)參數(shù)A、B、C皆為已知，則最小方差控制：由上式得：最小方差控制預(yù)測方程為：預(yù)測誤差為:控制律為：由此搭建Simulink仿真圖形可得：未加最小方差控制：圖2-2-1系統(tǒng)控制效果最小方差控制后：圖2-2-2最小方差系統(tǒng)控制效果從響應(yīng)圖中可

20、以看出采用最小方差調(diào)節(jié)器構(gòu)成閉環(huán)系統(tǒng)后，穩(wěn)定后響應(yīng)的偏差減小，并且使控制過程得到了很大的改善，讓系統(tǒng)很快就達(dá)到穩(wěn)定。二（3）分析;利用遞推算法直接估計最小方差控制器的參數(shù)，即基于最小方差控制自適應(yīng)控制，計算量小。取初值p(0)=106I、(0)=0，f的上界為0，期望輸出有yr(k)為幅值為10的方波信號，并為此建立一個新的估計模型。程序：a=1 1.2 0.35; b=1 0.4; c=1 1.1 0.28; d=4; %對象參數(shù)na=length(a)-1; nb=length(b)-1; nc=length(c)-1; %na、nb、nc為多項式A、B、C階次nf=nb+d-1; %nf

21、為多項式F的階次L=1000; %控制步數(shù)uk=zeros(d+nb,1); %輸入初值：uk(i)表示u(k-i);yk=zeros(na,1); %輸出初值yrk=zeros(nc,1); %期望輸出初值xik=zeros(nc,1); %白噪聲初值xiek=zeros(nc,1); %白噪聲估計初值yr=10*ones(L/4,1);-ones(L/4,1);ones(L/4,1);-ones(L/4+d,1); %期望輸出xi=sqrt(0.1)*randn(L,1); %白噪聲序列%RELS初值設(shè)置thetae_1=0.001*ones(na+nb+1+nc,1); %非常小的正數(shù)（

22、這里不能為0）P=106*eye(na+nb+1+nc);for k=1:L time(k)=k; y(k)=-a(2:na+1)*yk+b*uk(d:d+nb)+c*xi(k);xik; %采集輸出數(shù)據(jù) %遞推增廣最小二乘法 phie=-yk;uk(d:d+nb);xiek; K=P*phie/(1+phie'*P*phie); thetae(:,k)=thetae_1+K*(y(k)-phie'*thetae_1); P=(eye(na+nb+1+nc)-K*phie')*P; xie=y(k)-phie'*thetae(:,k); %白噪聲的估計值 %提取

23、辨識參數(shù) ae=1 thetae(1:na,k)' be=thetae(na+1:na+nb+1,k)' ce=1 thetae(na+nb+2:na+nb+1+nc,k)' if abs(be(2)>0.9 be(2)=sign(ce(2)*0.9; %MVC算法要求B穩(wěn)定 end if abs(ce(2)>0.9 ce(2)=sign(ce(2)*0.9; end e,f,g=sindiophantine(ae,be,ce,d); %求解單步Diophantine方程 u(k)=(-f(2:nf+1)*uk(1:nf)+ce*yr(k+d:-1:k+d-

24、min(d,nc);yrk(1:nc-d)-g*y(k);yk(1:na-1)/f(1); %求控制量 %更新數(shù)據(jù) thetae_1=thetae(:,k); for i=d+nb:-1:2 uk(i)=uk(i-1); end uk(1)=u(k); for i=na:-1:2 yk(i)=yk(i-1); end yk(1)=y(k); for i=nc:-1:2 yrk(i)=yrk(i-1); xik(i)=xik(i-1); xiek(i)=xiek(i-1); end if nc>0 yrk(1)=yr(k); xik(1)=xi(k); xiek(1)=xie; enden

25、dfigure(1);subplot(2,1,1);plot(time,yr(1:L),'r:',time,y);xlabel('k'); ylabel('y_r(k)、y(k)');legend('y_r(k)','y(k)'); axis(0 L -50 50);subplot(2,1,2);plot(time,u);xlabel('k'); ylabel('u(k)'); axis(0 L -1000 1000);figure(2)subplot(211)plot(1:L,th

26、etae(1:na,:);xlabel('k'); ylabel('參數(shù)估計a');legend('a_1','a_2'); axis(0 L -2 2);subplot(212)plot(1:L,thetae(na+1:na+nb+1+nc,:);xlabel('k'); ylabel('參數(shù)估計b、c');legend('b_0','b_1','c_1'); axis(0 L 0 3)圖2-3-1系統(tǒng)自適應(yīng)控制效果由圖可知最后未知參數(shù)達(dá)到一穩(wěn)定值則

27、表明實現(xiàn)了最小方差控制器的功能，使系統(tǒng)達(dá)到了穩(wěn)定。圖2-3-2系統(tǒng)自適應(yīng)參數(shù)控制效果由圖可知穩(wěn)定性問題是系統(tǒng)設(shè)計中固有的問題，最后未知參數(shù)達(dá)到一穩(wěn)定值則表明實現(xiàn)了最小方差控制器的功能，使系統(tǒng)達(dá)到了穩(wěn)定。三、 簡述基于李雅普諾夫穩(wěn)定性理論，用被控對象輸入、輸出變量構(gòu)成的MRAS控制算法及特點。答：模型參考自適應(yīng)控制系統(tǒng)(Model Reference Adaptive System MRAS)是一類重要的自適應(yīng)控制系統(tǒng)，其關(guān)鍵問題是：如何設(shè)計自適應(yīng)參數(shù)調(diào)整律，在保證系統(tǒng)穩(wěn)定性的同時使誤差信號趨于零。模型參考自適應(yīng)法的主要思想是：將不含未知參數(shù)的方程作為參考模型，而將含有待估計參數(shù)的方程作為可調(diào)模

28、型（自適應(yīng)模型），兩個模型具有相同物理意義的輸出量，兩個模型同時工作，利用輸出量的誤差構(gòu)成合適的自適應(yīng)律以調(diào)節(jié)可調(diào)模型參數(shù)，以達(dá)到控制對象的輸出跟蹤參考模型的目的。然后根據(jù)穩(wěn)定性原理得到速度估計自適應(yīng)規(guī)律，系統(tǒng)和速度的漸進(jìn)收斂性由波波夫超穩(wěn)定性來保證。因此，模型參考自適應(yīng)控制系統(tǒng)的工作過程可以看成是參考模型與實際系統(tǒng)響應(yīng)之間誤差的調(diào)整過程。模型參考自適應(yīng)系統(tǒng)算法通常有三個基本設(shè)計方法：參數(shù)局部最優(yōu)化、穩(wěn)定性和超穩(wěn)定性設(shè)計方法。由于模型參考自適應(yīng)系統(tǒng)的時變非線性特點，所以穩(wěn)定性問題是系統(tǒng)設(shè)計中固有的問題，一個完整的模型參考自適應(yīng)控制系統(tǒng)設(shè)計必須包括穩(wěn)定性分析。以李雅普諾夫函數(shù)為基礎(chǔ)的設(shè)計方法能構(gòu)成功地用來設(shè)計穩(wěn)定的模型參考自適應(yīng)系統(tǒng)，但是用穩(wěn)定性的方法推導(dǎo)系統(tǒng)自適應(yīng)控制規(guī)律存在兩個問題：第一，由于沒有一個合適的方法可供尋求系統(tǒng)的李雅普諾夫函數(shù)所用，因此設(shè)計者的經(jīng)驗顯得特別重要。第二，由于李雅普諾夫函