【正文】 三、附錄。分析系統(tǒng)的狀態(tài)穩(wěn)定和輸出穩(wěn)定是否一致。（4）令， ，初始狀態(tài)為零。用函數(shù)initial( )求出系統(tǒng)的零輸入響應(yīng)，并繪制相應(yīng)的狀態(tài)響應(yīng)曲線。與（1）的結(jié)果進(jìn)行比較。確定系統(tǒng)的穩(wěn)定性。二、 實(shí)驗(yàn)內(nèi)容1. 已知線性系統(tǒng)(a) (b) (c) (d)（1）用函數(shù)eig( )，pole( )和zpkdata( )求出系統(tǒng)的特征值和極點(diǎn)。Gm=minreal(G)返回2 states removed.a = x1 x1 1b = u1 x1 c = x1 y1 d = u1 y1 實(shí)驗(yàn)四 穩(wěn)定性一、 實(shí)驗(yàn)?zāi)康恼莆障到y(tǒng)穩(wěn)定性的概念。 C=[0 1 2]。1。0 1 3]。在命令窗中運(yùn)行下列命令 A=[0 0 1。 [Ak Bk Ck Tk]=kalmdec(A,B,C)返回Ac = 0 Bc = 0Cc = Tc = 4. 最小實(shí)現(xiàn)MATLAB 提供的函數(shù)minreal( )可直接得出系統(tǒng)的最小實(shí)現(xiàn)，其調(diào)用格式為 其中G為系統(tǒng)的LTI對象，為系統(tǒng)的一個(gè)最小實(shí)現(xiàn)。0]。 B=[1。1 0 3。Tk=rot90(T,2)。Bk=rot90(Bk1,2)。Ck1=C*inv(T)。Ak1=T*A*inv(T)。n2=rank(obsv(Ac(inc,inc),Cc(inc)))。T=To*Tc。To=[To2,zeros(m2,n1)。end[m1,n1]=size(To1)。if ncn, inc=1:nnc。 ic=nnc+1:n。nc=rank(ctrb(A,B))。 （335）可利用函數(shù)ctrbf( )、obsvf( )編制按能控性和能觀測性分解函數(shù)kalmdec( )，該函數(shù)調(diào)用格式為: 其中為給定系統(tǒng)的狀態(tài)空間模型，為分解以后系統(tǒng)的狀態(tài)空間模型，為線性變換矩陣。C=[0 1 2]。1。0 1 3]。end在命令窗中運(yùn)行下列命令 A=[0 0 1。 Co=C*inv(To)。 end end Ao=To*A*inv(To)。p]。 for i=0:m p(1:nm,m+1i:ni)=eye(nm,nm)。 end if rank(t1)==m break。 if rank(t2)==i t1=t2。 for i=2:n*m t2=[t1。n=size(A,1)。該函數(shù)的程序如下function [Ao, Bo, Co, To, k]=odes(A,B,C) %按能觀測性分解V=obsv(A,C)。為相應(yīng)的線性變換矩陣。我們編制的函數(shù)odes ( )能夠?qū)⑾到y(tǒng)分解為這種形式。 C=[0 1 2]。1。0 1 3]。例35 按能觀測分解例34中系統(tǒng)可在命令窗中運(yùn)行下列命令進(jìn)行分解 A=[0 0 1。MATLAB提供的函數(shù)obsvf( )可將系統(tǒng)(或狀態(tài))分解為如下形式: （333）該函數(shù)的調(diào)用格式為 其中為給定系統(tǒng)的狀態(tài)空間模型，為分解后系統(tǒng)的狀態(tài)空間模型。C=[0 1 2]。1。0 1 3]。end在命令窗中運(yùn)行下列命令 A=[0 0 1。 Cc=C*Tc。 end end Ac=inv(Tc)*A*Tc。]。 for i=0:m p(1:nm,m+1i:ni)=eye(nm,nm)。 end if rank(t1)==m break。 if rank(t2)==i t1=t2。if m~=n t1=U(:,1)。m=rank(U)。k表示能控狀態(tài)變量的數(shù)目。該函數(shù)的調(diào)用格式為 其中為給定系統(tǒng)的狀態(tài)空間模型，為分解后系統(tǒng)的狀態(tài)空間模型。 [Ac Bc Cc Tc Kc]=ctrbf(A,B,C)返回Ac = Bc = 0 0 Cc = Tc = 0Kc = 1 1 0另一種按能控性分解的形式為 （332）將由ctrbf( )函數(shù)得到的各系數(shù)矩陣均利用MATLAB提供的函數(shù)rot90（ ）旋轉(zhuǎn)就可得這種形式。0]。 B=[1。1 0 3。為相應(yīng)線性變換矩陣，返回系統(tǒng)中能控狀態(tài)變量的數(shù)目。3. 線性系統(tǒng)的結(jié)構(gòu)分解1) 按能控性分解如果線性系統(tǒng)的狀態(tài)不完全能控，則可通過非奇異線性變換, 將系統(tǒng)(或狀態(tài))分解為 能控和不能控兩部分。 rank(Vo) 返回ans= 3因?yàn)閞ank(Vo)=3=n，故系統(tǒng)狀態(tài)完全能觀測。0 1 0]。0 1 3]。例33 判斷例1中系統(tǒng)的能觀測性在命令窗中運(yùn)行下列命令 A=[1 3 2。能觀測性矩陣的秩即稱為能觀測性指數(shù)，表示系統(tǒng)能觀測狀態(tài)變量的數(shù)目。2. 能觀測性n階線性定常連續(xù)或離散系統(tǒng)狀態(tài)完全能觀測的充分必要條件是：能觀測性矩陣的秩為n。rank(Uy)返回ans= 1因?yàn)閞ank(Uy)=1=m，故系統(tǒng)是輸出能控的。Uc=ctrb(A,B)。C=[1 1]。B=[5。矩陣可以通過能控性矩陣得到，即例32 判斷系統(tǒng)的輸出能控性 在命令窗中運(yùn)行下列命令 A=[4 5。1 1]。B=[2 1。0 2 0。system is not pletely state controllable39。system is pletely state controllable39。n=length(A)。也可編制判斷系統(tǒng)狀態(tài)能控性的函數(shù)sctrb( )，其程序如下：function sctrb(A,B) %判斷系統(tǒng)狀態(tài)的能控性Uc=ctrb(A,B)。 Uc=ctrb(A,B)。 1 1。0 1 3]。例31 判斷系統(tǒng)的能控性 在命令窗中運(yùn)行下列命令 A=[1 3 2。能控性矩陣可用MATLAB提供的函數(shù)ctrb( )自動(dòng)產(chǎn)生，其調(diào)用格式為:其中A,B分別為系統(tǒng)矩陣和輸入矩陣，為能控性矩陣。判斷所得系統(tǒng)的能控性和能觀測性，驗(yàn)證其是否最小實(shí)現(xiàn)。它與（1）中的傳遞函數(shù)模型是否一致？為何？令初始狀態(tài)為零，用MATLAB計(jì)算系統(tǒng)的單位階躍輸出響應(yīng)，并繪制和記錄相應(yīng)的曲線。它與（1）中的傳遞函數(shù)模型是否一致？為何？令初始狀態(tài)為零，用MATLAB計(jì)算系統(tǒng)的單位階躍輸出響應(yīng)，并繪制和記錄相應(yīng)曲線。它與（1）中所得的傳遞函數(shù)模型是否一致？為何？令初始狀態(tài)為零，用MATLAB計(jì)算系統(tǒng)的單位階躍輸出響應(yīng)，并繪制和記錄相應(yīng)曲線。令初始狀態(tài)為零，用MATLAB計(jì)算系統(tǒng)的單位階躍輸出響應(yīng)，繪制和記錄相應(yīng)的曲線。Vo=obsv(A,C)。0 3]。Uc=ctrb(A,B)。B=[0。 Gm=minreal(Gss)1 state removed. a = x1 x2 x1 2 4 x2 0 3 b = u1 x1 0 x2 c = x1 x2 y1 0 d = u1 y1 0 Continuoustime model.最小實(shí)現(xiàn)能控性 A=[2 4。G=zpk(z,p,k)。p=[1 2 3]。 Vo=obsv(A,C)。0 1]。G=ss(A,B,C,D)。C=[1 1]。B=[4。rank(Uc)ans = 1能觀性： A=[3 4。0]。0 1]。G=ss(A,B,C,D)。C=[1 1]。B=[4。rank(Uy)ans = 1(3) a能控性： A=[3 4。Uc=ctrb(A,B)。C=[1 1]。B=[4。 rank(Vo) ans = 1輸出能控性： A=[3 4。 C=[1 1]。 rank(Uc)ans = 1能觀性： A=[3 4。1]。1 0]。觀察和記錄這些曲線。當(dāng)輸入改變時(shí), 每個(gè)狀態(tài)變量的響應(yīng)曲線是否隨著改變？能否根據(jù)這些曲線判斷系統(tǒng)狀態(tài)的能控性？（3） 將給定的狀態(tài)空間表達(dá)式變換為對角標(biāo)準(zhǔn)型，判斷系統(tǒng)的能控性和能觀測性，與（1）的結(jié)果是否一致？為何？（4）令（3）中系統(tǒng)的初始狀態(tài)為零， 輸入分別為單位階躍函數(shù)和單位脈沖函數(shù)。用MATLAB函數(shù)計(jì)算系統(tǒng)的狀態(tài)響應(yīng)和輸出響應(yīng)，并繪制相應(yīng)的響應(yīng)曲線。說明狀態(tài)能 控性和輸出能控性之間有無聯(lián)系。學(xué)會(huì)用MATLAB求最小實(shí)現(xiàn)。學(xué)會(huì)用MATLAB進(jìn)行結(jié)構(gòu)分解。學(xué)會(huì)用MATLAB判斷能控性和能觀測性。在命令窗中運(yùn)行該程序后得到系統(tǒng)的階躍響應(yīng)曲線如圖6。y= dstep(G,H,C,D,u,n)。u=1。C=[1 1]。H=[0。G=[0 1。w39。color39。pos39。 %計(jì)算phek=iztrans(phez) %計(jì)算在命令窗中運(yùn)行該程序后返回phek =[ 4/3*(1/5)^n1/3*(4/5)^n, 5/3*(1/5)^n5/3*(4/5)^n][ 4/15*(1/5)^n+4/15*(4/5)^n, 1/3*(1/5)^n+4/3*(4/5)^n]編制程序%ex252求出系統(tǒng)的單位階躍響應(yīng)并繪制階躍響應(yīng)曲線。syms z。該程序如下%ex251G=[0 1。 圖5 狀態(tài)和輸出響應(yīng)曲線 圖6 階躍響應(yīng)曲線3. 線性定常離散系統(tǒng)狀態(tài)方程的解如果線性定常離散系統(tǒng)的狀態(tài)空間表達(dá)式為 （327）初始狀態(tài)為，那么狀態(tài)方程的解為 （328）系統(tǒng)的狀態(tài)轉(zhuǎn)移矩陣為 （329）系統(tǒng)的輸出為 （3210）例25 已知系統(tǒng) 求狀態(tài)轉(zhuǎn)移矩陣，求初始狀態(tài)為零時(shí)系統(tǒng)的單位階躍響應(yīng)。FontSize39。x239。FontSize39。x139。FontSize39。y39。k39。:k39。,tt)。 yt(i)=subs(y,39。t39。 %計(jì)算輸出響應(yīng)的解析解for i=1:61tt=*(i1)。x=ilaplace(G*B*u)。phet=ilaplace(G)。syms s。0]。5 6]。)。,39。,[50 50 200 150],39。 %ex24figure(39。)可得狀態(tài)響應(yīng)和輸出響應(yīng)的數(shù)值解以及相應(yīng)的曲線，如圖4。,t,y,39。plot(t,x,39。u=exp(t)。G=ss(A,B,C,D)。C=[1 2]。B=[2。在命令窗中運(yùn)行下列命令A(yù)=[0 1。)返回圖3。,t,y,39。 plot(t,x,39。y=yo+yu。)返回圖2。,t,yu,39。plot(t,xu,39。)；u=ones(size(t))。,39。,[50 50 200 150],39。figure(39。)返回圖1。,t,yo,39。 plot(t,xo,39。1]。 t=0::10。 D=0。1]。2 3]。在命令窗中運(yùn)行下列命令，建立狀態(tài)空間模型，計(jì)算系統(tǒng)在初始狀態(tài)作用下的狀態(tài)響應(yīng)和輸出響應(yīng)，并繪制相應(yīng)的響應(yīng)曲線。該函數(shù)的調(diào)用格式為 或 其中G為狀態(tài)空間模型，x0為初始狀態(tài)。216。該函數(shù)的調(diào)用格式為：其中u為與時(shí)間向量t對應(yīng)的輸入向量。216。 函數(shù)impulse( ) 可直接求取線性系統(tǒng)的單位脈沖響應(yīng)。 或 此時(shí)不返回任何變量，而自動(dòng)地繪制單位階躍響應(yīng)輸出曲線。