自动控制理论常用MATLAB函数

1、自动控制理论中常用MATLAB函数第一部分 建立数学模型在MATLAB中，线性定常（linear time invariant, 简称为 LTI）系统可以用4种数学模型描述，即传递函数(TF)模型、零极点增益(ZPK)模型和状态空间(SS)模型以及SIMULINK结构图。前三种数学模型是用数学表达式表示的，且均有连续和离散两种类型，通常把它们统称为LTI模型。1传递函数模型（TF 模型）令单输入单输出线性定常连续和离散系统的传递函数分别为 (1-1)和 (1-2)在MATLAB中，连续系统和离散系统的传递函数都用分子/分母多项式系数构成的两个行向量num和den表示，即，系统的传递函数模型用M

2、ATLAB提供的函数tf( )建立。函数tf ( )不仅能用于建立系统传递函数模型，也能用于将系统的零极点增益模型和状态空间模型转换为传递函数模型。该函数的调用格式如下： G=tf(num,den) 返回连续系统的传递函数模型。 G=tf(num,den,Ts) 返回连续系统的传递函数模型。其中，Ts为采样周期，当Ts=-1或者Ts= 时，系统的采样周期未定义。 Gtf=tf(G) 可将任意的LTI模型转换为传递函数模型。例1-1 已知一个系统的传递函数为 建立传递函数模型。 在命令窗中运行下列命令 num=6;den=1 6 11 6;G=tf (num, den) 返回Transfer f

3、unction: 6-s3 + 6 s2 + 11 s + 62零极点增益模型（ZPK模型）系统的零极点增益模型是传递函数模型的一种特殊形式。令线性定常连续和离散系统的零极点形式的传递函数分别为 (1-3)和 (1-4)在MATLAB中，连续和离散系统的零点和极点都用行向量和表示，即，。系统的零极点增益模型用MATLAB提供的函数zpk ( )建立。函数zpk( )不仅能用来建立系统零极点增益模型，也能用于将系统的传递函数模型和状态空间模型转换为零极点增益模型。该函数的调用格式如下： G=zpk(z,p,k) 返回连续系统的传递函数模型。 G=zpk(z,p,k,Ts) 返回离散系统的零极点增

4、益模型。 Ts为采样周期，当Ts=-1或者Ts=时，系统的采样周期未定义。 Gzpk=zpk(G) 可将任意的LTI模型转换为零极点增益模型。例1-2 已知系统的传递函数为 建立系统的零极点增益模型。 在命令窗中运行下列命令 z= ;p=-1 -2 -3;k=6;G=zpk(z,p,k) 返回Zero/pole/gain:6-(s+1) (s+2) (s+3)注意：无零点时，设z为空。3状态空间模型(SS模型)令多输入多输出线性定常连续和离散系统的状态空间表达式分别为 （1-5）和 （1-6）在MATLAB中，连续系统和离散系统的状态空间模型都用MATLAB提供的函数ss ( )建立。函数ss

5、 ( )不仅能用于建立系统的状态空间模型，也能用于将系统的传递函数模型和零极点增益模型转换为状态空间模型。该函数的调用格式如下： G=ss(A,B,C,D) 返回连续系统的状态空间模型。 G=ss(A,B,C,D,Ts) 返回离散系统的状态空间模型。 Ts为采样周期，当Ts=-1或者Ts=时，系统的采样周期未定义。 Gss=ss(G) 可将任意的LTI模型转换为状态空间模型。例1-3 已知系统的状态空间表达式为 建立系统的状态空间模型。在命令窗中运行下列命令 A=0 1 0;0 0 1;-6 -11 -6;B=0;0;1;C=6 0 0;D=0;G=ss(A,B,C,D) 返回a = x1 x

6、2 x3 x1 0 1 0 x2 0 0 1 x3 -6 -11 -6b = u1 x1 0 x2 0 x3 1c = x1 x2 x3 y1 6 0 0d = u1 y1 0Continuous-time model.注意：D=0不能缺省。4模型转换上述三种LTI模型之间可以通过函数tf( )，zpk( )和ss( )相互转换。线性定常系统的传递函数模型和零极点增益模型是唯一的，但系统的状态空间模型是不唯一的。函数ss( )只能将传递函数模型和零极点增益模型转换为一种指定形式的状态空间模型。例1-4 已知系统的传递函数，建立TF模型，将其转换为ZPK模型和SS模型。再将转换得到的SS模型转换

7、为TF模型。 在命令窗中运行下列命令num=5;den=1 4 5 2;Gtf=tf(num,den);Gzpk=zpk(Gtf) Gss=ss(Gtf)Gtf1=tf(Gss)返回Zero/pole/gain: 5-(s+2) (s+1)2a = x1 x2 x3 x1 -4 -2.5 -0.5 x2 2 0 0 x3 0 2 0b = u1 x1 1 x2 0 x3 0c = x1 x2 x3 y1 0 0 1.25d = u1 y1 0Continuous-time model.Transfer function: 5-s3 + 4 s2 + 5 s + 25. 函数tfdata( ),

8、 zpkdata( ), ssdata( )函数tfdata( ) 可得到传递函数模型的分子分母多项式系数。其调用格式为 num,den=tfdata(G,v) 其中G为系统LTI模型。和den分别为分子和分母多项式的系数向量。函数zpkdata( ) 可得到零极点增益模型的零点、极点和增益。其调用格式为 z,p,k=zpkdata(G,v) 其中G为系统LTI模型。z和p分别为零点和极点向量，k为增益。函数ssdata( ) 可得到状态空间模型的系数矩阵。其调用格式为 A,B,C,D=ssdata(G,v) 其中G为系统LTI模型。A,B,C,D为系数矩阵。6. 模型的连接1) 串联连接设线

9、性定常系统1和2的LTI模型分别为和。在MATLAB中, 两者的串联连接（参教材）可由下面命令实现。 G=G1*G2 其中G为整个系统的LTI模型。多个系统的串联连接可由下面命令实现。 G=G1*G2*Gn 2) 并联连接设线性定常系统1和2的LTI模型分别为和。在MATLAB中, 两者的并联连接（参教材）可由下面命令实现。 G=G1+G2 其中G为整个系统的LTI模型。 多个系统的并联连接可由下面命令实现。 G=G1+G2+Gn 3) 反馈连接设线性定常系统1和2的LTI模型分别为和。在MATLAB中, 两者的反馈连接（参教材）可由MATLAB提供的函数feedback( )实现。该函数的调

10、用格式为： G=feedback(G1,G2,Sign) 其中G为整个系统LTI模型。如果Sign=-1或省略Sign变量，则表示负反馈。如果Sign=1, 则表示正反馈。例1-5 已知单位负反馈系统的开环传递函数为 求闭环系统的传递函数。在命令窗中运行下列命令 num=1;den=1 1;G=tf(num,den); GB=feedback(G,1)返回Transfer function: 1-s + 2 7. 状态空间表达式的线性变换线性定常系统状态空间表达式为，假设存在一个非奇异矩阵T将原状态x变换为z=Tx， 则状态z对应的状态空间表达式为，其中，。MATLAB提供函数ss2ss( )

11、可完成状态空间模型的相似变换。该函数调用格式为 Gt=ss2ss(G,T) 其中G为原状态空间模型。T为变换矩阵。为经变换得到的状态空间模型。例1-6 考虑一个系统，它的状态空间表达式为 由于该系统的系统矩阵A为友矩阵，所以可由其特征值构造变换矩阵。令变换矩阵T为。在命令窗中运行下列命令 A=0 1;-2 -3;B=1 1;C=1 0;G=ss(A,B,C,0); T=inv(1 1;-1 -2); G1=ss2ss(G,T)返回a = x1 x2 x1 -1 0 x2 0 -2b = u1 x1 3 x2 -2c = x1 x2 y1 1 1d = u1 y1 0Continuous-tim

12、e model.8. 时间延迟系统的传递函数模型 带有延迟环节的系统不具有有理函数的标准形式，在MATLAB中建立这类系统的模型，要由一个属性设置函数set（）来实现。该函数的调用格式为： set(G,InputDelay,T);其中，T为延迟时间。由此修改后，模型G就已具有时间延迟特性。例1-12 已知系统的传递函数为建立TF模型。在命令窗中运行下列命令 num=1;den=1 1;G=tf(num,den); T=1; set(G,InputDelay,T); G返回 Transfer function: 1exp(-1*s) * - s + 19. Simulink建模 相关内容参阅有关

13、书籍。第二部分 系统分析1. 系统响应的数值求解函数step( ) 可直接求取线性连续系统的单位阶跃响应。该函数的调用格式为： y=step(G,t);其中G为给定系统LTI模型，t为时间向量。通常取t=0:dt:t-end，其中t-end为终值时间，而dt为时间步长。y为系统输出。 y=step(G);这时时间向量t自动生成。 y,t,x=step(G,t); y,t,x=step(G);如果G为状态空间模型，则x为系统状态向量，否则x将返回空矩阵。 step(G); step(G,t);后两种格式不返回任何变量，而自动地绘制单位阶跃响应输出曲线。例2-1 已知系统为 或 初始状态为零，试绘

14、制输出响应曲线，并确定性能指标。 在命令窗中运行运行命令，可建立状态空间模型 A=0 1;-2 -2;B=0;1;C=1 1;D=0;Gss=ss(A,B,C,D);或运行命令，建立传递函数模型 num=1 1;den=1 2 2;Gtf=tf(Gss);也可将转换为传递函数，即在命令窗中输入再运行 step(Gtf) 或 step(Gss)返回响应曲线图2.1，通过点击鼠标右键选择从菜单中选择Characteristics，点击Peak Response（超调量），Settilng Time（调节时间），Rise Time（上升时间），Steady State（稳态），可得系统各性能指标。图

15、2.1 函数impulse( ) 可直接求取线性系统的单位脉冲响应。该函数的调用格式与函数step( )的调用格式相似。 函数lsim( ) 可直接求取线性系统在任意输入信号作用下的响应。该函数的调用格式为： y=lsim( G,u,t) 其中u为与时间向量t对应的输入向量。当然还可以和step( )函数一样有其它的调用格式。函数initial( ) 可求解系统的零输入响应。该函数的调用格式为： y=initial( G,x0) y,t,x=initial( G,x0) 其中G为状态空间模型，x0为初始状态。例2-2 已知系统为 初始状态为，试求u(t)为单位阶跃函数时系统状态响应和输出响应，

16、并绘制状态响应曲线和输出响应曲线。在命令窗中运行下列命令，建立状态空间模型，计算系统在初始状态作用下的状态响应和输出响应，并绘制相应的响应曲线。 A=0 1;-2 -3; B=0;1; C=1 1; D=0; G=ss(A,B,C,D); t=0:0.5:10; x0=1;-1;yo,t,xo=initial(G,x0,t); plot(t,xo,:,t,yo,-)返回图2-2。 图2-2 图2-3 在命令窗中继续运行下列命令，计算系统在输入作用下的状态响应和输出响应，并绘制相应的响应曲线。u=ones(size(t);yu,t,xu=lsim(G,u,t);plot(t,xu,:,t,yu,

17、-)返回图2-3。再继续运行下列命令求系统总的状态响应和输出响应，并绘制相应的响应曲线。y=yo+yu; x=xo+xu; plot(t,x,:,t,y,-)返回图2-4。 图2-4 图2-5 例2-3 已知系统 或 求出系统在初始状态为零，且时系统的状态响应和输出响应。在命令窗中运行下列命令A=0 1;-5 -6;B=2;0;C=1 2;D=0;G=ss(A,B,C,D);t=0:0.5:20;u=exp(-t);y,t,x=lsim(G,u,t);plot(t,x,:k,t,y,-k)可得状态响应和输出响应的数值解以及相应的曲线，如图2-6。 图2-6 图2-7例2-4 已知系统 求初始状

18、态为零时系统的单位阶跃响应。编制程序%ex24求出系统的单位阶跃响应并绘制阶跃响应曲线。该程序如下%ex24G=0 1;-0.16 -1;H=0;1;C=1 -1;D=0;u=1;n=20;y= dstep(G,H,C,D,u,n);dstep(G,H,C,D,u,n);在命令窗中运行该程序后得到系统的阶跃响应曲线如图2-7。2. 系统的能控性和能观测性分析（1）线性定常系统状态能控性的判断n阶线性定常连续或离散系统状态完全能控的充分必要条件是：能控性矩阵的秩为n。能控性矩阵可用MATLAB提供的函数ctrb( )自动产生，其调用格式为： Uc=ctrb(A,B) 其中A,B分别为系统矩阵和输

19、入矩阵，为能控性矩阵。能控性矩阵的秩即称为能控性指数，表示系统能控状态变量的数目，可由MATLAB提供的函数rank( )求出。例2-5 判断系统的能控性。 在命令窗中运行下列命令 A=1 3 2;0 2 0;0 1 3; B=2 1; 1 1;-1 -1; Uc=ctrb(A,B); rank(Uc) 返回ans = 2因为rank(Uc)=2n，所以系统的状态不完全能控。（2）线性定常系统状态能观测性的判断n阶线性定常连续或离散系统状态完全能观测的充分必要条件是：能观测性矩阵的秩为n。能观测性矩阵可以用MATLAB提供的函数obsv( )自动产生，其调用格式为： Vo=obsv(A,C)

20、其中A, C分别为系统矩阵和输出矩阵，为能观测性矩阵。能观测性矩阵的秩即称为能观测性指数，表示系统能观测状态变量的数目。可由MATLAB提供的函数rank( )求出。例2-6 判断例2-5中系统的能观测性。在命令窗中运行下列命令 A=1 3 2;0 2 0;0 1 3; C=2 1 1;0 1 0; Vo=obsv(A,C); rank(Vo) 返回ans= 3因为rank(Vo)=3=n，故系统状态完全能观测。3. 系统的稳定性分析（1）根据系统的极点和特征值判定稳定性线性定常连续系统渐近稳定的充分必要条件是：系统传递函数的全部极点均位于S左半平面；或系统矩阵A的特征值均具有负实部。利用MA

21、TLAB提供的下列几种函数可以确定系统的极点和特征值。 函数eig( )的调用格式为： V=eig(A) 返回方阵A的特征值。函数roots( )的调用格式为： roots(den) 其中den为多项式的系数行向量。可计算多项式方程的解。函数pole ( )的调用格式为： pole(G) 其中G为系统的LTI对象。可计算系统传递函数的极点。函数pzmap( )的调用格式为： pzmapG) 其中G 为LTI对象。绘制系统的零点和极点。例2-7 已知系统 或 判断系统的稳定性。在命令窗中运行下列命令 A=0 1;-5 -6;B=2;0;C=1 2;D=0;G=ss(A,B,C,D); eig(A

22、) 返回ans = -1-5即系统的特征值为-1，-5，故系统是稳定的。或在命令窗中运行下列命令 num=2 -8;den=1 6 5;Gtf=tf(num,den); roots(den)返回ans = -5-1即系统的极点为-1，-5，故系统是稳定的。（2）用李氏第二法判定稳定性线性定常系统为。若对于任意给定的对称正定矩阵Q，存在唯一满足Lyapunov方程正定矩阵P，则系统大范围渐近稳定。可用MATLAB提供的函数lyap( )求解Lyapunov方程，该函数调用格式为： P=lyap(A,Q) 返回矩阵P 矩阵P的符号性质可用函数posdef( )判定。该函数的调用格式为： S=pos

23、def(P) 该函数的程序如下。function s=posdef(P) %判定矩阵的正定性r=length(P);for i=1:rpp(i)=det(P(1:i,1:i);endk=find(pp A=0 1;-5 -6;Q=eye(2);P=lyap(A,Q) s=posdef(P)返回P = 1.0000 -0.5000 -0.5000 0.5000s =matrix is positive definite matrix4. 绘制系统的根轨迹利用MATLAB提供的函数rlocus()可以绘制系统的根轨迹。该函数的调用格式有以下几种： rlocus(num,den) rlocus(nu

24、m,den,K) rlocus(G) rlocus(G,K)其中，G为开环系统的对象模型，K为用户自己选择的增益向量。如果用户不给出K向量，则该命令函数会自动选择K向量。 MATLAB中还有一个函数rlocfind()，该函数允许用户求取根轨迹上指定点处的开环增益值，并将该增益下所有的闭环极点显示出来。这个函数的调用函数为： K,P=rlocus(G,K)例2-9 已知系统的开环传递函数模型为绘制系统的根轨迹。在命令窗中运行下列命令 num=1;den=conv(1,1,1,2),0;G=tf(num,den); rlocus(G); grid返回根轨迹图2-8，继续运行下面语句， K,P=r

25、locfind(G)并在根轨迹图上点击一点，返回Select a point in the graphics windowselected_point = 0.2275 + 1.8758iK = 12.3349P = -3.4546 0.2273 + 1.8759i 0.2273 - 1.8759i图2-85. 绘制系统的奈奎斯特图利用MATLAB提供的函数nyquist()可以绘制系统的根轨迹。当命令不包含左端返回变量时，函数nyquist()仅产生奈奎斯特图，调用格式有以下几种： nyquist(num,den) nyquist(num,den,w) nyquist(G) nyquist(

26、G,w)其中，w为用户给出的频率向量，在这些频率点上，将对系统的频率响应进行计算。若没有指定的w向量，则该函数自动选择频率向量进行计算。当命令包含左端返回变量时，函数nyquist()不产生奈奎斯特图，而是将计算结果返回到矩阵re、im和w中。矩阵re和im分别表示频率响应的实部和虚部，它们都是由向量w中指定的频率点计算得到的。调用格式有以下几种： re,im,w=nyquist(G) re,im,w=nyquist(G,w)例2-10 给定二阶典型环节试绘制奈奎斯特图。在命令窗中运行下列命令 num=1;den=1 0.8 1;nyquist(num,den)返回奈奎斯特图2-9。图2-96

27、. 绘制系统的伯德图利用MATLAB提供的函数bode( )可以绘制系统的伯德图。当命令不包含左端返回变量时，函数bode( )直接产生伯德图。当命令包含左端返回变量时，函数bode( )计算出的幅值和相角将返回到相应的矩阵中，这时不产生伯德图。调用格式有以下几种： mag,phase,w=bode(num,den) mag,phase,w=bode(num,den,w) mag,phase,w=bode(G) mag,phase,w=bode(G,w)其中，magphase包含系统频率响应的幅值和相角，这些幅值和相角是在用户指定的频率点上计算得到的。如果用户不指定频率w，MATLAB会自动产

28、生w向量。相角的单位是度，幅值的单位不是分贝。例2-11 已知系统的开环传递函数模型为绘制系统的伯德图。在命令窗中运行下列命令 num=10*1 1;den=1 7 0;bode(num,den) grid返回伯德图2-10。图2-10利用MATLAB提供的函数margin( )可以不仅绘制系统的伯德图，而且可以求取系统的幅值裕度和相角裕度。调用格式如下： Gm,Pm,Wcg,Wcp=margin(G)其中，Gm为幅值裕度的值；Wcg为相角穿越频率；Pm为相角裕度的值；Wcp为幅值穿越频率。例2-12 对例2-11中系统求相应的幅值裕度和相角裕度。在命令窗中运行下列命令 num=10*1 1;

29、den=1 7 0;G=tf(num,den); margin(G) grid 返回伯德图2-11。图2-11在命令窗中运行下列命令 Gm,Pm,Wcg,Wcp=margin(G); Gm,Pm,Wcg,Wcp返回ans = Inf 126.0768 NaN 7.2725第三部分 系统的综合1. 状态反馈和输出反馈设线性定常系统的状态空间表达式为 如果采用状态反馈控制规律，其中r是参考输入，则状态反馈闭环系统的状态空间表达式为 其传递函数为如果采用输出反馈控制律，其中r是参考输入，则输出反馈闭环系统的状态空间表达式为 其传递函数为。2. 极点配置如果SISO线性定常系统完全能控，则可通过适当的状态反馈, 将闭环系统极点配置到任意期望的位置。假定期望的闭环极点为，,，则闭环系统的特征方程为 原系统的特征方程为 状态反馈矩阵K可由下面两种方法得出。1) 其中 2) （Ackermann公式） 其中 MATLAB提供的函数acker( )是用Ackermann公式求解状态反馈阵K。该函数的调用格式为： K=acker(A,B,P)其中A和B分别为系统矩阵和输入矩阵。P是期望极点构成的向量。 MATLAB提供的函数place( )也可求出状态反馈阵K。该函数的调用格式为： K=place