11第3章主要內容：

3.1時域分析的拉普拉斯變換法3.1.1連續(xù)時間函數的拉普拉斯變換3.1.2時域函數的拉氏反變換法3.1.3時域函數的部分分時展開法3.2時域分析的函數命令方法3.3作業(yè)與實驗2第3章主要內容：23.1時域分析的拉普拉斯變換法3

利用拉普拉斯變換法對線性系統進行時域分析的一般步驟：（1）對系統的傳遞函數模型進行部分分式展開，將其變?yōu)楹唵蝹鬟f函數之和；（2）利用拉普拉斯反變換，得到系統的輸出時間響應函數；（3）繪制系統的響應曲線;（4）通過改變系統的參數，觀察系統輸出響應的變化情況，對系統的時域特性（如性能指標等）進行分析。本節(jié)介紹上述過程的MATLAB實現方法。3.1時域分析的拉普拉斯變換法3利用拉普拉斯變2.1.1生成模型的常用函數命令格式4連續(xù)時間信號f(t)的單邊拉普拉斯變換F(s)的定義：記為：拉普拉斯反變換的定義：記為：2.1.1生成模型的常用函數命令格式4連續(xù)時間信號f(t)s=tf(‘s’)：生成以s為變量的遞函數s。此時，s既是傳遞函數也是指定變量。

sys=tf(num,den)：生成傳遞函數模型，num,den分別為模型的分子和分母多項式系數向量。sys=tf(num,den,‘td’,v)：生成延遲時間td=v的傳遞函數模型。sys=zpk(z,p,k)：生成零極點增益模型，z,p,k分別為零點、極點和增益向量。[num,den]=rmodel(n,P)：隨機生成一個n階連續(xù)的傳遞函數模型，該系統具有P個輸出。[num,den]=ord2(wn,z)：生成固有頻率為wn，阻尼系數為z的連續(xù)二階系統模型系統。[num,dent]=pade(L,n)：返回延遲環(huán)節(jié)G(s)=e-Ls近似為n階多項式傳遞函數的num和den。5s=tf(‘s’)：生成以s為變量的遞函數s。此時，s既是傳2.1.2有理多項式分式傳遞函數模型的建立6

控制理論中，系統傳遞函數的有理多項式分式模型通常表示為：

在MATLAB中，多項式是用其降冪系數的排列構成的數組來表示的，稱為多項式系數向量。

系統傳遞函數的分子和分母的多項式系數向量分別表示為：num=[bm,bm-1,…,b1,b0]den=[an,an-1,…,a1,a0]2.1.2有理多項式分式傳遞函數模型的建立67

在MATLAB語言中有專門對信號進行正、反拉氏變換的函數命令。常用拉氏變換的函數命令格式及其說明見表3-1。F=laplace(f)：f(t)的拉氏變換，結果為F(s)，默認變量為sF=laplace(f,v)：f(t)拉氏變換，并用v代替s，結果為F(v)F=laplace(f,v,u)：f(v)拉氏變換，并用u代替v，結果為F(u)f=ilaplace(F)

：F(s)的拉氏反變換，結果為f(t)，變量為tf=ilaplace(F,u)：F(s)拉氏反變換，用u代替s，結果為f(u)f=ilaplace(F,v,u)：

F(v)拉氏反變換，用u代v，結果為f(u)說明：默認情況下，MTLAB給出的是單邊拉氏變換；使用函數laplace()及ilaplace()之前，需要用syms命令對所有要用到的變量進行符號變量的定義。7在MATLAB語言中有專門【例3-1】試求函數

的拉氏變換式，并用拉氏反變換觀察變換結果。解：MATLAB程序如下：

>>clear;%清除所有變量symstAwbs

%定義符號變量t,A,w,b,sft=A*sin(w*t+b);%定義f(t)的符號函數ft的表達式Fs=laplace(ft)

%求ft的拉氏變換式Fs,即F(s)運行結果：Fs=A*(cos(b)*w/(s^2+w^2)+sin(b)*s/(s^2+w^2))8【例3-1】試求函數當取參數A=1、b=0、w=1時，帶入Fs>>A=1;b=0;w=1;%設定參數值Fs=A*(cos(b)*w/(s^2+w^2)+sin(b)*s/(s^2+w^2))

%計算給定參數下的F(s)運行結果：Fs=1/(s^2+1)即得到

F(s)=L[f(t)]=L[sin(t)]=1/(s2+1)9當取參數A=1、b=0、w=1時，帶入Fs9

可利用拉氏反變換對上述結果進行檢驗：>>ft=ilaplace(Fs)%求Fs的拉氏反變換式ft運行結果：ft=sin(t)即

f(t)=L-1[F(s)]=L-1[1/(s2+1)]=sin(t)說明:MATLAB中的拉氏變換的函數命令laplace()，在默認情況下是指拉氏右變換，其運行結果是單邊函數。如L[f(t)]=laplace(sin(t))，實際上是指對ft=sin(t)u(t)的雙邊拉氏變換，

其中u(t)為單位階躍函數。10可利用拉氏反變換對上述結果進行檢驗：10

利用MATLAB拉氏反變換可以求出傳遞函數的時間函數，或系統的時域輸出響應?！纠?-3】試求傳遞函數

的時間函數g(t)。解：求解G(s)的拉氏反變換程序如下：113.1.2時域函數的拉氏反變換法利用MATLAB拉氏反變換可以求出傳遞函數1>>clear;symss;Gs=(s^3+4*s^2+8*s+5)/(s^2+3*s+4);gt=ilaplace(Gs)運行結果：gt=Dirac(t)+Dirac(1,t)-1/7*exp(-3/2*t)*7^(1/2)*sin(1/2*7^(1/2)*t)+exp(-3/2*t)*cos(1/2*7^(1/2)*t)其中，即得：12>>clear;symss;12【例3-4】某質量-彈簧-阻尼(m-k-c)系統如圖3-1所示。其中，f(t)為外加力，x(t)為m的位移。試求：1）系統的MATLAB傳遞函數；2）輸入為單位階躍（f(t)=u(t)）時的時間響應函數x(t)；3）繪制m=3;k=2;c=1時的時域響應曲線。13【例3-4】某質量-彈簧-阻尼(m-k-c)系統如圖3-1所解：1）由牛頓第二定律可得到該系統模型的數學表示：上式的拉氏變換后的傳遞函數模型為：生成上述模型的MATLAB模型：>>symsmkcsGs=1/(m*s^2+c*s+k)%生成G(s)的符號傳遞函數Gs運行結果：Gs=1/(m*s^2+c*s+k)14解：1）由牛頓第二定律可得到該系統模型的數學表示：14152）f(t)=u(t)的拉氏變換為F(s)=1/s。所以輸出傳函為求上式的MATLAB拉氏反變換：>>Xs=Gs*(1/s)%建立輸出傳函X(s)xt=ilaplace(Xs)%求解輸出X(s)的時域響應

152）f(t)=u(t)的拉氏變換為F(s)=1/s。所以16運行結果：xt=-4*exp(-1/2*c/m*t)/(4*k*m-c^2)*m*cos(1/2*((4*k*m-c^2)/m^2)^(1/2)*t)+1/k*exp(-1/2*c/m*t)/(4*k*m-c^2)*c^2*cos(1/2*((4*k*m-c^2)/m^2)^(1/2)*t)-1/k*exp(-1/2*c/m*t)/(4*k*m-c^2)*m*((4*k*m-c^2)/m^2)^(1/2)*c*sin(1/2*((4*k*m-c^2)/m^2)^(1/2)*t)+1/k3）繪制m=3;k=2;c=1時的響應曲線。程序如下：>>symstm=3;k=2;c=1;%參數賦值

16運行結果：xt=-4*exp(-1/2*c/m*t)/(4*k*m-c^2)*m*cos(1/2*((4*k*m-c^2)/m^2)^(1/2)*t)+1/k*exp(-1/2*c/m*t)/(4*k*m-c^2)*c^2*cos(1/2*((4*k*m-c^2)/m^2)^(1/2)*t)-1/k*exp(-1/2*c/m*t)/(4*k*m-c^2)*m*((4*k*m-c^2)/m^2)^(1/2)*c*sin(1/2*((4*k*m-c^2)/m^2)^(1/2)*t)+1/k運行結果：xt=-1/2*exp(-1/6*t)*cos(1/6*23^(1/2)*t)-1/46*exp(-1/6*t)*23^(1/2)*sin(1/6*23^(1/2)*t)+1/217xt=-4*exp(-1/2*c/m*t)/(4*k*m-c>>t=[0:0.1:10]*pi;%取繪圖區(qū)間采樣點xt=-1/2*exp(-1/6*t).*cos(1/6*23^(1/2)*t)-...1/46*exp(-6*t).*23^(1/2).*sin(1/6*23^(1/2)*t)+1/2;%計算所有采樣點上的函數值plot(t,xt);grid%繪圖，加柵格結果如圖3-2所示。18>>t=[0:0.1:10]*pi;1919

對于復雜的分式多項式傳函模型，首先，采用部分分式法將其分解成如下部分分式展開的形式：其中，pi為極點，

ri為有理函數Y(s)的留數，即k為余項；其次，對其求拉氏反變換；最后，由拉氏變換簡表得到時間響應函數。203.1.3時域函數的部分分時展開法對于復雜的分式多項式傳函模型，首先，采用部20321通常有下列幾中情況及其相應的處理方法：1）單極點時，即pi為不相同極點。則Y(s)的拉普拉斯逆變換為：2）有m重極點時。假設有m個重極點

，Y(s)的拉普拉斯逆變換為：3）有共軛極點時，逆變換較為復雜，需具體分析處理（參考教材），

21通常有下列幾中情況及其相應的處理方法：，22

在MATLAB中，對于上述三種情況均可利用如下函數命令格式[r,p,k]=residue(num,den)對分式多項式進行部分分式展開，然后再根據極點的不同情況進行相應的處理。其中，num,den為Y(s)的分子和分母多項式系數向量，[r,p,k]為返回的留數，極點，

和余項。22在MATLAB中，對于上述三種情況均可利用如下【例3-5】試求輸出傳遞函數的時間響應。解：MATLAB程序如下：>>clear;num=[1-2];den=[13310];%系數向量[r,p,k]=residue(num,den)%求部分分式參數運行結果：r=2.00002.00003.0000-2.000023【例3-5】試求輸出傳遞函數23p=-1.0000-1.0000-1.00000k=[]可見，系統有3個重根p1=p2=p3=-1和p4=0，則由式（3-3）和式（3-4）得到由式（3-4）知，m=3,對應的逆變換為y(t)=L-1[Y(s)]?？梢灾苯永美献儞Q表得出結果。也可以通過MATLAB求解，程序如下：24p=24>>symssft1=ilaplace(3/(s+1)^3),ft2=ilaplace(2/(s+1)^2),ft3=ilaplace(2/(s+1)),ft4=ilaplace(-2/s)%求各分式逆變換得結果：ft1=3/2*t^2*exp(-t)ft2=2*t*exp(-t)ft3=2*exp(-t)ft4=-2>>ft=ft1+ft2+ft3+ft4%求f(t)的結果：ft=3/2*t^2*exp(-t)+2*t*exp(-t)+2*exp(-t)-2即25>>symss25【例3-7】設一個RLC網絡結構如圖所示。其中ui為輸入信號，uc為輸出信號，網絡各變量的初始狀態(tài)均為零。試求：1）建立RLC網絡的傳遞函數

；2）當ui=Eu(t)，用拉氏逆變換函數命令求uc；3）繪制不同參數時的單位階躍響應uc(t)。解：1）建立RLC網絡的傳遞函數

(1)建立時域模型。根據R、L、C電壓定律及回路電壓方程

，可得RLC網絡的微分方程模型：其中，

26【例3-7】設一個RLC網絡結構如圖所示。26（2）建立傳函模型

對上述微分模型的兩邊取拉氏變換，經整理得到RLC網絡拉氏變換的標準二階系統傳遞函數：2）當ui=Eu(t)時，由于L[ui]=L[E]=E/s，所以輸出的拉氏變換為：利用MATLAB拉氏逆變換函數命令求uc。程序如下：27（2）建立傳函模型27>>clear

;symstsEwnetwd%定義變量t,s,E,ωn,ζUc=E*wn^2/((s^2+2*et*wn*s+wn^2)*s);%建立函數uc=ilaplace(Uc,s,t)%拉氏逆變換求輸出uc(t)simple(uc)%簡化結果結果：ans=E*wn^2*(1/wn^2+1/(4*et^2*wn^2-4*wn^2)^(1/2)*(1/(-et*wn+1/2*(4*et^2*wn^2-4*wn^2)^(1/2))*exp((-et*wn+1/2*(4*et^2*wn^2-4*wn^2)^(1/2))*t)-1/(-et*wn-1/2*(4*et^2*wn^2-4*wn^2)^(1/2))*exp((-et*wn-1/2*(4*et^2*wn^2-4*wn^2)^(1/2))*t)))28>>clear

;283）繪制不同參數時的單位階躍響應uc(t)。>>clear

;E=1;wn=3;t=0:0.02:6;%設置參數E和wn,時間采樣點tet=[0.1,0.2,0.5,0.95];%設置4個不同的阻尼向量etfori=1:4;%對et循環(huán)計算uc=E*wn^2*(1/wn^2+1/(4*et(i)^2*wn^2-4*wn^2)^(1/2)*(1/(-et(i)*wn+1/2*(4*et(i)^2*wn^2-4*wn^2)^(1/2))*exp((-et(i)*wn+1/2*(4*et(i)^2*wn^2-4*wn^2)^(1/2)).*t)-1/(-et(i)*wn-1/2*(4*et(i)^2*wn^2-4*wn^2)^(1/2)).*exp((-et(i)*wn-1/2*(4*et(i)^2*wn^2-4*wn^2)^(1/2)).*t)));plot(t,uc,'k-');holdon;%繪制響應曲線并保持圖形end;293）繪制不同參數時的單位階躍響應uc(t)。29title('uc(t)的二階欠阻尼響應曲線')gtext('et=0.1')%在圖中對應的響應曲線的合適位置%點擊鼠標左鍵添加阻尼參數標示'et=0.1'gtext('et=0.2')gtext('et=0.5')gtext('et=0.95')xlabel('t');ylabel('uc(t)');%添加x和y軸標示't'和'uc(t)'grid;%添加網格運行得到響應曲線見圖3-4。30title('uc(t)的二階欠阻尼響應曲線')3031313.2時域分析的函數命令方法

上述介紹的MATLAB拉氏變換法，不但可以用來對通常情況下的微分方程、分式多項式、傳遞函數及系統結構等函數進行數學上的處理與分析，還可用于控制系統的時域分析。MATLAB還提供了特殊的函數命令，專用于對系統的分析。本節(jié)介紹常用MATLAB線性系統時域分析的函數命令格式及其使用方法。表3-2列出了常用時域響應函數命令格式及說明。323.2時域分析的函數命令方法上述介紹的MA【例3-8】已知系統傳遞函數為

試求其單位脈沖響應。解：使用專用函數命令來實現，程序如下：>>clear;n=[3];d=[1,1,3];%傳遞函數的系數向量Gs=tf(n,d);%建立模型G(s)impulse(Gs)%G(s)的單位脈沖響應運行結果如圖3-4所示。33【例3-8】已知系統傳遞函數為333434【例3-9】已知系統傳遞函數為求其單位階躍響應。解：程序如下：>>clear;n=[1,7];d=[1,1.5,8,7];Gs=tf(n,d);step(Gs)運行結果如圖3-5所示。35【例3-9】已知系統傳遞函數為3536【例3-14】已知某二階系統的超調量Mp=16.3%。試根據下列條件確定二階系統的開環(huán)傳遞函數G(s)：1）上升時間tr=0.1s；2）峰值時間tp=0.114s；3）調節(jié)時間ts=2s。

解：由二階系統模型結構可知，其傳遞函數由兩個參數決定：阻尼比ζ和自然振蕩頻率ωn?？筛鶕A系統的性能指標中的任意兩個進行求解即可確定模型參數。所以可以利用MATLAB的方程求解函數命令實現。具體格式：

X=solve('eq1','eq2',...,'eqn',val1,val2,...,valn)其中，eqi和vali，i=1,2,…,n,分別為方程和待求變量。36【例3-14】已知某二階系統的超調量Mp=16.3%。371）求解同時滿足超調量Mp=16.3%和上升時間tr=0.1s時的阻尼比ζ和自然振蕩頻率ωn。即求解下列方程組中的兩個變量ζ和ωn。設ζ=x，ωn=y，則上述方程改寫為371）求解同時滿足超調量Mp=16.3%和上升時間tr=038計算程序如下：>>clear;symsxyMptrx=solve('Mp-exp(-x*pi/sqrt(1-x^2))=0',x);%求x=ζx=subs(x,Mp,0.163)y=solve('tr-(pi-atan(sqrt(1-x^2)/x))/(y*sqrt(1-x^2))=0',y);%求y=ωny=subs(y,tr,0.1);y=subs(y,x)得到的結果為：x=0.5000；y=24.1852。即得二階系統模型為：38計算程序如下：392）求解同時滿足超調量Mp=16.3%和峰值時間tp=0.114s時的阻尼比ζ和自然振蕩頻率ωn。即求解下列方程組計算程序如下：>>clear;symsxyMptpx=solve('Mp-exp(-x*pi/sqrt(1-x^2))=0',x);%求阻尼比x=subs(x,Mp,0.163)y=solve('tp-pi/(y*sqrt(1-x^2))=0',y);%求y=ωny=subs(y,tp,0.114);y=subs(y,x)所得結果：x=ζ=0.5000；y=ωn=31.8219。將參數ζ和ωn代入G(s)即可。392）求解同時滿足超調量Mp=16.3%和峰值時間tp=040【例2-12】已知傳遞函已知傳遞函數G1和G2分別為試以G1作為前向通道，以G2作為反饋通道時的負反饋等效閉環(huán)模型G。解：實現反饋的等效傳函程序如下：方法1：利用多項式系數求解。>>clear;num1=[1,11];den1=[1,3,21];num2=[2,1];den2=[3,5,