西建大《现代控制理论基础实验报告》.docx

1、西建大现代控制理论基础实验报告实 验 报 告课程 线性系统理论基础 实验日期 2013 年5月 日专业班级 姓名 学号 同组人 实验名称 MATLAB控制工具箱的应用及线性系统的运动分析 评分 批阅教师签字 一、实验目的1、学习掌握MATLAB控制工具箱中的基本命令的操作方法；2、掌握线性系统的运动分析方法。二、实验内容（1）自选控制对象模型，应用以下命令，并写出结果。1） step, damp, pzmap, rlocus, rlocfind, bode, margin, nyquist；2） tf2ss, ss2tf, tf2zp, zp2ss；3） ss2ss, jordan, cano

2、n, eig。（2）掌握线性系统的运动分析方法1)已知，求。(用三种方法求解)2) 利用MATLAB求解书上例2.8题，并画出状态响应和输出响应曲线，求解时域性能指标。(加图标题、坐标轴标注及图标)3) 利用MATLAB求解书上例2.12题，并画出状态响应和输出响应曲线。(加图标题、坐标轴标注及图标)4) P36 1.4（2） 1.5（3）；P56 2.3（3）三、实验环境大楼机房MATLAB6.X软件四、实验原理（或程序框图）及步骤1、学习掌握MATLAB控制工具箱中基本命令的操作设系统的模型如式(1-1)所示： (1-1)其中A为nn维系数矩阵；B为nm维输入矩阵；C为pn维输出矩阵；D为

3、pm维传递矩阵，一般情况下为0。系统的传递函数阵和状态空间表达式之间的关系如式(1-2)所示： (1-2)式(1-2)中，表示传递函数阵的分子阵，其维数是pm；表示传递函数阵的分母多项式，按s降幂排列的后，各项系数用向量表示。例1.1 已知SISO系统的状态空间表达式为(1-3)式，求系统的传递函数。 (1-3)程序:%首先给A、B、C阵赋值；A=0 1 0;0 0 1;-4 -3 -2;B=1;3;-6;C=1 0 0;D=0;%状态空间表达式转换成传递函数阵的格式为num,den=ss2tf(a,b,c,d,u)num,den=ss2tf(A,B,C,D,1) 程序运行结果:num = 0

4、 1.0000 5.0000 3.0000den = 1.0000 2.0000 3.0000 4.0000从程序运行结果得到系统的传递函数为： (1-4)例1.2 从系统的传递函数(1-4)式求状态空间表达式。程序:num =1 5 3;den =1 2 3 4;A,B,C,D=tf2ss(num,den)程序运行结果:A = B = -2 -3 -4 1 1 0 0 0 0 1 0 0C = D =1 5 3 0由于一个系统的状态空间表达式并不唯一, 例1.2程序运行结果虽然不等于式(1-3)中的A、B、C阵，但该结果与式(1-3)是等效的。不妨对上述结果进行验证。例1.3 对上述结果进行

5、验证编程。%将例1.2上述结果赋值给A、B、C、D阵；A =-2 -3 -4;1 0 0; 0 1 0;B =1;0;0;C =1 5 3;D=0；num,den=ss2tf(A，B，C，D,1)程序运行结果与例1.1完全相同。例1.4 给定系统，求系统的零极点增益模型和状态空间模型，并求其单位脉冲响应及单位阶跃响应。解：num=1 2 1 3;den=1 0.5 2 1;sys=tf(num,den) %系统的传递函数模型 Transfer function: s3 + 2 s2 + s + 3-s3 + 0.5 s2 + 2 s + 1sys1=tf2zp(num,den) %系统的零极点

6、增益模型sys1 = -2.1746 0.0873 + 1.1713i 0.0873 - 1.1713isys2=tf2ss(sys) %系统的状态空间模型模型；或用a,b,c,d=tf2ss(num,den)形式a = -0.5000 -2.0000 -1.0000 1.0000 0 0 0 1.0000 0b = 1 0 0c = 1.5000 -1.0000 2.0000d = 1impulse(sys2) %系统的单位脉冲响应 图1-1 系统的单位脉冲响应step(sys2) %系统的单位阶跃响应： 图1-2 系统的单位阶跃响应五、程序源代码1)已知，求。(用三种方法求解)状态转移矩阵

7、的指数矩阵计算法a=0 1;-2 -3;syms t;eat1=expm(a*t)eat1 = -exp(-2*t)+2*exp(-t), exp(-t)-exp(-2*t) -2*exp(-t)+2*exp(-2*t), 2*exp(-2*t)-exp(-t)拉氏反变换计算法a=0 1;-2 -3;syms s t;G=inv(s*eye(size(a)-a)eat2=ilaplace(G)G = (s+3)/(s2+3*s+2), 1/(s2+3*s+2) -2/(s2+3*s+2), s/(s2+3*s+2)eat2 = -exp(-2*t)+2*exp(-t), exp(-t)-exp

8、(-2*t) -2*exp(-t)+2*exp(-2*t), 2*exp(-2*t)-exp(-t)非奇异变换法a=0 1;-2 -3;sym t;P,D=eig(a);Q=inv(P);eat3=P*expm(D*t)*Qeat3 = -exp(-2*t)+2*exp(-t), exp(-t)-exp(-2*t) -2*exp(-t)+2*exp(-2*t), 2*exp(-2*t)-exp(-t)例2.8a=-1 0;0 -2;b=1;1;c=1.5 0.5;d=0;G=ss(a,b,c,d);x0=2;3;syms s t;G0=inv(s*eye(size(a)-a);x1=ilapl

9、ace(G0)*x0G1=inv(s*eye(size(a)-a)*b;x2=ilaplace(G1/s)x=x1+x2y=c*xfor I=1:61; tt=0.1*(I-1); xt(:,I)=subs(x(:),t,tt); yt(I)=subs(y,t,tt);end;plot(0:60,xt;yt);x1 = 2*exp(-t) 3*exp(-2*t)x2 = 1-exp(-t) 1/2-1/2*exp(-2*t)x = exp(-t)+1 5/2*exp(-2*t)+1/2y =3/2*exp(-t)+7/4+5/4*exp(-2*t)输出响应输出响应输出响应例2.12a=0 1;

10、-0.16 -1;b=1;1;c=1 0x0=1;-1;syms z n k;thta=inv(z*eye(size(a)-a)*z;thtak=iztrans(thta,k)uz=z/(z-1);xk=iztrans(thta*x0+thta/z*b*uz)y=c*xkfor I=1:61; tt=0.1*(I-1); xt(:,I)=subs(x(:),t,tt); yt(I)=subs(y,t,tt);end;plot(0:60,xt;yt);thtak = 4/3*(-1/5)k-1/3*(-4/5)k, 5/3*(-1/5)k-5/3*(-4/5)k -4/15*(-1/5)k+4/

11、15*(-4/5)k, -1/3*(-1/5)k+4/3*(-4/5)kxk = -17/6*(-1/5)n+22/9*(-4/5)n+25/18 17/30*(-1/5)n-88/45*(-4/5)n+7/18y =-17/5*(-1/5)n+22/5*(-4/5)n+1P36 1.4（2） a=2 1 4;0 2 0;0 0 1;b=1 0;3 4;2 1;c=3 5 1;d=0 0;num,den=ss2tf(a,b,c,d,1)num = 0 20.0000 -29.0000 -13.0000den = 1 -5 8 -4num=1 4 2 2;0 3 1 1;den=1 2 3 0

12、2;A,B,C,D=tf2ss(num,den)A = -2 -3 0 -2 1 0 0 0 0 1 0 0 0 0 1 0B = 1 0 0 0C = 1 4 2 2 0 3 1 1D = 0 0P36 1.5（3）a=0 1;-6 -5;b=1;0;c=1 -1;d=0;G=ss(a,b,c,d);x0=1;1;syms z k;thta=inv(z*eye(size(a)-a)*z;thtak=iztrans(thta,k)uz=z/(z-1);xk=iztrans(thta*x0+thta/z*b*uz)y=c*xk thtak = 3*(-2)k-2*(-3)k, (-2)k-(-3

13、)k -6*(-2)k+6*(-3)k, -2*(-2)k+3*(-3)k xk = 3*(-2)n-5/2*(-3)n+1/2 -6*(-2)n+15/2*(-3)n-1/2y =9*(-2)n-10*(-3)n+1P56 2.3（3）a=0 1;-6 -5;b=1;0;c=1 -1;d=0;G=ss(a,b,c,d);x0=1;1;syms s t;G0=inv(s*eye(size(a)-a);x1=ilaplace(G0)*x0G1=inv(s*eye(size(a)-a)*b;x2=ilaplace(G1/s)x=x1+x2y=c*xfor I=1:61; tt=0.1*(I-1);

14、 xt(:,I)=subs(x(:),t,tt); yt(I)=subs(y,t,tt);end;plot(0:60,xt;yt);实 验 报 告课程 线性系统理论基础 实验日期 2013 年5月 日专业班级 姓名 学号 同组人 实验名称 系统的能控性、能观测性、稳定性分析及实现 评分 批阅教师签字 一、实验目的加深理解能观测性、能控性、稳定性、最小实现等观念。掌握如何使用MATLAB进行以下分析和实现。1、系统的能观测性、能控性分析；2、系统的稳定性分析；3、系统的最小实现。二、实验内容（1）能控性、能观测性及系统实现（a）了解以下命令的功能；自选对象模型，进行运算，并写出结果。gram,

15、ctrb, obsv, lyap, ctrbf, obsvf, mineral；（b）已知连续系统的传递函数模型，当a 分别取-1，0，1时，判别系统的能控性与能观测性；（c）已知系统矩阵为，判别系统的能控性与能观测性；（d）求系统的最小实现。（2）稳定性（a）代数法稳定性判据已知单位反馈系统的开环传递函数为：，试对系统闭环判别其稳定性（b）根轨迹法判断系统稳定性已知一个单位负反馈系统开环传递函数为，试在系统的闭环根轨迹图上选择一点，求出该点的增益及其系统的闭环极点位置，并判断在该点系统闭环的稳定性。（c）Bode 图法判断系统稳定性已知两个单位负反馈系统的开环传递函数分别为用Bode 图法判

16、断系统闭环的稳定性。（d）判断下列系统是否状态渐近稳定、是否BIBO稳定。三、实验环境大楼机房MATLAB6.X软件四、实验原理（或程序框图）及步骤1、系统能控性、能观性分析设系统的状态空间表达式如（1-1）所示。系统的能控性、能观测性分析是多变量系统设计的基础，包括能控性、能观测性的定义和判别。系统状态能控性定义的核心是：对于线性连续定常系统（1-1）,若存在一个分段连续的输入函数u(t)，在有限的时间（t1-t0）内，能把任一给定的初态x(t0)转移至预期的终端x(t1)，则称此状态是能控的。若系统所有的状态都是能控的，则称该系统是状态完全能控的。能控性判别分为状态能控性判别和输出能控性判

17、别。状态能控性分为一般判别和直接判别法，后者是针对系统的系数阵A是对角标准形或约当标准形的系统，状态能控性判别时不用计算，应用公式直接判断，是一种直接简易法；前者状态能控性分为一般判别是应用最广泛的一种判别法。输出能控性判别式为： （2-1）状态能控性判别式为： （2-2）系统状态能观测性的定义：对于线性连续定常系统（2-1）,如果对t0时刻存在ta，t0ta0Flagz=1;endenddisp(系统的零极点模型为);z,p,k系统的零极点模型为z = -2.7306 + 2.8531i -2.7306 - 2.8531i -1.5388 p = -4.0000 -3.0000 -2.000

18、0 -1.0000k = 1.0000 图2-1 系统的阶跃响应if Flagz=1disp(系统不稳定);else disp(系统是稳定的);end运行结果为:系统是稳定的step(A,B,C,D);五、程序源代码(1)能控性、能观测性及系统实现(b)已知连续系统的传递函数模型，当a 分别取-1，0，1时，判别系统的能控性与能观测性；num=1 -1; %a=-1den=1 10 27 18;a,b,c,d=tf2ss(num,den);Qc=ctrb(a,b)Qc = 1 -10 73 0 1 -10 0 0 1 rank(Qc)ans =3Qo=obsv(a,c) Qo = 0 1 -1

19、 1 -1 0 -11 -27 -18 rank(Qo)ans = 3num=1 0; %a=0den=1 10 27 18;a,b,c,d=tf2ss(num,den);Qc=ctrb(a,b)Qc = 1 -10 73 0 1 -10 0 0 1 rank(Qc)ans =3 Qo=obsv(a,c)Qo = 0 1 0 1 0 0 -10 -27 -18 rank(Qo)ans = 3num=1 1; %a=1den=1 10 27 18;a,b,c,d=tf2ss(num,den);Qc=ctrb(a,b)Qc = 1 -10 73 0 1 -10 0 0 1 rank(Qc)ans

20、= 3 Qo=obsv(a,c)Qo = 0 1 1 1 1 0 -9 -27 -18 rank(Qo)ans = 2（c）已知系统矩阵为，判别系统的能控性与能观测性a=6.666 -10.6667 -0.3333;1 0 1;0 1 2;b=0;1;1;c=1 0 2;Qc=ctrb(a,b)Qc = 0 -11.0000 -84.9926 1.0000 1.0000 -8.0000 1.0000 3.0000 7.0000 rank(Qc)ans = 3 Qo=obsv(a,c)Qo = 1.0000 0 2.0000 6.6660 -8.6667 3.6667 35.7689 -67.4

21、375 -3.5551 rank(Qo)ans =3(d)求系统的最小实现。 num=0 0 1 1;den=1 10 27 18;G=tf(num,den);Gs=ss(G);Gm=minreal(Gs);1 state removed. Am=Gm.aAm = 3.5391 -12.1540 5.1323 -12.5391 Bm=Gm.bBm = 0.0606 -0.242 Cm=Gm.cCm = 0.2500 0.0625 Dm=Gm.dDm = 0 （2）稳定性（a）已知单位反馈系统的开环传递函，试对系统闭环判别其稳定性Anum=100 200; den=1 21 20; A,B,C,

22、D=tf2ss(num,den)z,p,k=ss2zp(A,B,C,D,1)z = -2p = -20-1k = 100step(A,B,C,D)(b)已知一个单位负反馈系统开环，试在系统的闭环根轨迹图上选择一点，求出该点的增益及其系统的闭环极点位置，并判断在该点系统闭环的稳定性传递函数为num=1 3;den=1 13 54 83 63;W=tf(num,den)rlocus(W); title(系统根轨迹)(c)已知两个单位负反馈系统的开环传递函数分别为用Bode 图法判断系统闭环的稳定性。num=2.7; den=1 5 4 0; bode(num,den); title(系统的伯德图)

23、 num=2.7;den=1 5 -4 0;bode(num,den);title(系统的伯德图)实 验 报 告课程 线性系统理论基础 实验日期 2013 年5月 日专业班级 姓名 学号 同组人 实验名称 状态反馈极点配置方法的研究 评分 批阅教师签字 一、实验目的1掌握状态反馈系统的极点配置；2研究不同配置对系统动态特性的影响。二、实验内容原系统如图3-2所示。图中，X1和X2是可以测量的状态变量。图3-1 系统结构图试设计状态反馈矩阵,使系统加入状态反馈后其动态性能指标满足给定的要求: (1) 已知：K=10，T=1秒，要求加入状态反馈后系统的动态性能指标为： %20%，ts1秒。(12)

24、 已知：K=1，T=0.05秒，要求加入状态反馈后系统的动态性能指标为： %5%，ts0.5秒。 状态反馈后的系统，如图3-3所示：图3-2 状态反馈后系统结构图三、实验环境大楼机房MATLAB6.X软件四、实验原理（或程序框图）及步骤一个受控系统只要其状态是完全能控的，则闭环系统的极点可以任意配置。极点配置有两种方法：采用变换矩阵T，将状态方程转换成可控标准型，然后将期望的特征方程和加入状态反馈增益矩阵K后的特征方程比较，令对应项的系数相等，从而决定状态反馈增益矩阵K；基于Carlay-Hamilton理论，它指出矩阵状态矩阵A满足自身的特征方程，改变矩阵特征多项式的值，可以推出增益矩阵K，

25、这种方法推出增益矩阵K的方程式叫Ackermann公式。例4.1 某控制系统的状态方程描述如下：通过状态反馈使系统的闭环极点配置在P=-30，-1.2,-2.44i位置上,求出状态反馈阵K,并绘制出配置后系统的时间响应曲线。解： A=-10 -35 -50 -24;1 0 0 0;0 1 0 0;0 0 1 0;B=1;0;0;0;C=1 7 24 24;D=0;disp(原系统的极点为);p=eig(A)运算结果为：原极点的极点为p = -4.0000 -3.0000 -2.0000 -1.0000P=-30;-1.2;-2.4+sqrt(-16);-2.4-sqrt(-16);K=place(A,B,P)K = 26.0000 172.5200 801.7120 759.3600disp(配