系统仿真综合实验报告.docx

1、系统仿真综合实验报告(此文档为word格式，下载后您可任意编辑修改！)合肥工业大学电气与自动化工程学院综合实验报告实 验 名 称:系统仿真综合实验姓 名:学 号:专 业 班 级:实 验 地 点:指 导 教 师:成 绩:日 期:2012年7月实验一 MATLAB基本操作实验目的1熟悉MATLAB实验环境，练习MATLAB命令、m文件、Simulink的基本操作。2利用MATLAB编写程序进行矩阵运算、图形绘制、数据处理等。3利用Simulink建立系统的数学模型并仿真求解。实验原理 MATLAB环境是一种为数值计算、数据分析和图形显示服务的交互式的环境。MATLAB有3种窗口，即：命令窗口（Th

2、e Command Window）、m-文件编辑窗口（The Edit Window）和图形窗口（The Figure Window），而Simulink另外又有Simulink模型编辑窗口。 1命令窗口（The Command Window） 当MATLAB启动后，出现的最大的窗口就是命令窗口。用户可以在提示符“”后面输入交互的命令，这些命令就立即被执行。在MATLAB中，一连串命令可以放置在一个文件中，不必把它们直接在命令窗口内输入。在命令窗口中输入该文件名，这一连串命令就被执行了。因为这样的文件都是以“.m”为后缀，所以称为m-文件。 2m-文件编辑窗口（The Edit Window）

3、 我们可以用m-文件编辑窗口来产生新的m-文件，或者编辑已经存在的m-文件。在MATLAB主界面上选择菜单“File/New/M-file”就打开了一个新的m-文件编辑窗口；选择菜单“File/Open”就可以打开一个已经存在的m-文件，并且可以在这个窗口中编辑这个m-文件。 3图形窗口（The Figure Window） 图形窗口用来显示MATLAB程序产生的图形。图形可以是2维的、3维的数据图形，也可以是照片等。 MATLAB中矩阵运算、绘图、数据处理等内容参见教材自动控制系统计算机仿真的相关章节。Simulink是MATLAB的一个部件，它为MATLAB用户提供了一种有效的对反馈控制系

4、统进行建模、仿真和分析的方式。有两种方式启动Simulink:1在Command window中，键入simulink，回车。2单击工具栏上Simulink图标。启动Simulink后，即打开了Simulink库浏览器（Simulink library browser）。在该浏览器的窗口中单击“Create a new model（创建新模型）”图标，这样就打开一个尚未命名的模型窗口。把Simulink库浏览器中的单元拖拽进入这个模型窗口，构造自己需要的模型。对各个单元部件的参数进行设定，可以双击该单元部件的图标，在弹出的对话框中设置参数。 实验内容1 用MATLAB可以识别的格式输入下面两个

5、矩阵 再求出它们的乘积矩阵C，并将C矩阵的右下角23子矩阵赋给D矩阵。赋值完成后，调用相应的命令查看MATLAB工作空间的占用情况。答案：1） 程序：A=1 2 3 3;2 3 5 7;1 3 5 7;3 2 3 9;1 8 9 4;B=1+4i 4 3 6 7 8;2 3 3 5 5 4+2i;2 6+7i 5 3 4 2;1 8 9 5 4 3;C=A*BD=C(4:5,4:6)结果：C = 1.0e+002 * Columns 1 through 4 0.1400 + 0.0400i 0.5200 + 0.2100i 0.5100 0.4000 0.2500 + 0.0800i 1.03

6、00 + 0.3500i 1.0300 0.7700 0.2400 + 0.0400i 0.9900 + 0.3500i 1.0000 0.7100 0.2200 + 0.1200i 1.0800 + 0.2100i 1.1100 0.8200 0.3900 + 0.0400i 1.1400 + 0.6300i 1.0800 0.9300 Columns 5 through 6 0.4100 0.3100 + 0.0400i 0.7700 0.5900 + 0.0600i 0.7000 0.5100 + 0.0600i 0.7900 0.6500 + 0.0400i 0.9900 0.7000

7、 + 0.1600iD = 82.0000 79.0000 65.0000 + 4.0000i 93.0000 99.0000 70.0000 +16.0000i工作空间:Name Value A B C D 82. 000 + 0. 000i,79. 000 + 0. 000i,65. 000 + 4. 000i;93. 000 + 0. 000i,99. 000 + 0. 000i,70. 000 + 16. 000i小结：matlab通过确认下标，可以对矩阵进行插入子块、提取子块和重排子块的操作。如果提取子块时，n或m是常数，则返回指定的行列；如果n或m是向量，则返回的是指定矩阵的子块。

8、2 分别用for和while循环结构编写程序，求出答案：程序1：s=0;for k=0:63; s=s+2k;enddisp(The sum is),s结果：The sum iss = 1.8447e+019程序2： s=0;k=0;while k=63; s=s+2k; k=k+1;enddisp(The sum is),s结果：The sum iss = 1.8447e+019小结：在for循环语句，循环体内不能出现对循环控制变量的重新设置，否则会出错；while循环语句中，在语句内必须有可以修改循环控制变量的命令，否则将陷入死循环，除非循环语句中有控制退出循环的语句。3 选择合适的步距绘

9、制出下面的图形（1），其中（2），其中答案：（1）t=0.1:0.001:10;y=sin(t.(-1);plot(t,y)小结：t=0.1：0.001：10是代表一个步长为0.001的向量，t.(-1)是代表这个向量的点运算，而不是代表这个向量的向量运算，运算才不会出错。（2）t=linspace(-pi,pi,60);y=sin(tan(t)-tan(sin(t);plot(t,y)小结：linspace函数能够生成线性分度的向量。4 对下面给出的各个矩阵求取矩阵的行列式、秩、特征多项式、范数。， ， 答案：程序：A=7.5 3.5 0 0;8 33 4.1 0;0 9 103 -1.5;

10、0 0 3.7 19.3;B=5 7 6 5;7 10 8 7;6 8 10 9;5 7 9 10;C=1 2 3 4;5 6 7 8;9 10 11 12;13 14 15 16;D=3 -3 -2 4;5 -5 1 8;11 8 5 -7;5 -1 -3 -1;detA=det(A),rankA=rank(A),polyA=poly(A),normA=norm(A)detB=det(B),rankB=rank(B),polyB=poly(B),normB=norm(B)detC=det(C),rankC=rank(C),polyC=poly(C),normC=norm(C)detD=det

11、(D),rankD=rank(D),polyD=poly(D),normD=norm(D)结果：detA = 4.3222e+005rankA = 4polyA = 1.0e+005 * 0.0000 -0.0016 0.0713 -1.0659 4.3222normA = 103.7228detB = 1.0000rankB = 4polyB = 1.0000 -35.0000 146.0000 -100.0000 1.0000normB = 30.2887detC = 4.7332e-030rankC = 2polyC = 1.0000 -34.0000 -80.0000 -0.0000

12、-0.0000normC = 38.6227detD = 595.0000rankD = 4polyD = 1.0000 -2.0000 -32.0000 532.0000 595.0000normD = 16.6958小结：det、rank、poly、norm分别求出矩阵的行列式、秩、特征多项式以及范数5 求解下面的线性代数方程，并验证得出的解真正满足原方程。(a)， (b) 答案：(a).程序：A=7 2 1 -2;9 15 3 -2;-2 -2 11 5;1 3 2 13;B=4;7;-1;0;X=inv(A)*B结果：X = 0.4979 0.1445 0.0629 -0.0813 (

13、b)A=1 3 2 13;7 2 1 -2;7 15 3 -2; -2 -2 11 5;B=9 0;6 4;11 7;-2 -1;X=inv(A)*B结果X = 0.9247 0.4695 0.4140 0.2186 -0.1907 0.07910.5550 -0.0987小结：求矩阵A*B=C,已知A,C求B，则应用逆矩阵来求解即B=(A*(-1)*C6 假设有一组实测数据x0.10.20.30.40.50.60.70.80.91y2.32012.64702.97073.28853.60083.90904.21474.51914.82325.1275用最小二乘法拟合，求出相应的二次函数。答案

14、：程序x=0.1:0.1:1;y=2.3201 2.6470 2.9707 3.2885 3.6008 3.9090 4.2147 4.5191 4.8232 5.1275;p=polyfit(x,y,2)结果：p = -0.1563 3.2828 1.9967小结：polyfit是最小二乘拟合的函数，本题说明用最小二乘拟合求得的二次函数为y=-0.1563*x2+3.2828*x+1.99677 考虑线性微分方程 (1) 试用Simulink搭建起系统的仿真模型，并绘制出仿真结果曲线。(2) 将给定的微分方程转换成状态方程，并建立S函数，再利用Simulink进行仿真。答案：(1)参数设置：

15、Integrator：Initial condition分别设置成0.2，0.5，0.5，1Gain0,1,2, 3：gain 分别设置成3，3，4，5Sine Wave： Frequency(rad/sec)设置成 4 ；phase(rad)设置成 pi/3Transfer Fcn ：Numerator coefficients设置成 1 ； Denominator coefficients 设置成1 5Transfer Fcn1：Numerator coefficients设置成1 ； Numerator coefficients设置成1 3 Add ： List of signs 设置成

16、 +-仿真模型：（2）建立相应的S函数：functionsys,x0,str,ts= reds7_2 (t,x,u,flag,A,B,C,D)switch flag case 0 sys,x0,str,ts=mdlInitializeSizes(A,D); case 1 sys=mdlDerivatives(t,x,u,A,B); case 3 sys=mdlOutputs(t,x,u,C,D); case 2,4,9 sys=; otherwise error(Unhandled flag=,num2str(flag);endfunctionsys,x0,str,ts=mdlInitiali

17、zeSizes(A,D)sizes=simsizes;sizes.NumContStates=size(A,1);sizes.NumDiscStates=0;sizes.NumOutputs=size(A,1)+size(D,1);sizes.NumInputs=size(D,2);sizes.DirFeedthrough=1;sizes.NumSampleTimes=1;sys=simsizes(sizes);x0=1;0.5;0.5;0.2;str=;ts=-1 0;function sys=mdlDerivatives(t,x,u,A,B) sys=A*x+B*u;function sy

18、s=mdlOutputs(t,x,u,C,D) sys=C*x+D*u;x;由此建立仿真模型仿真曲线：小结：可以根据微分方程做出状态图并赋初值，得到系统的仿真模型以及仿真结果曲线；通过构建S函数，可以将状态方程封装在S函数，以简化运算。8 建立下图所示非线性系统的Simulink模型，并观察在单位阶跃信号输入下系统的输出曲线和误差曲线。答案：参数设置：Relay: Switch on point设置成0.5；Switch on point设置成0.5；Output when on 设置成2.5；Output when off 设置成0Relay1: Switch on point 设置成-0.

19、5；Switch on point设置成-0.5 ；Output when on 设置成2.5； Output when off 设置成0Transport Delay: Time delay设置成0.4Transfer Fcn: Numerator coefficients设置成30;Denominator coefficients设置成1 0.5 1Transfer Fcn1: Numerator coefficients设置成30;Denominator coefficients设置成1 6.5 1仿真模型：仿真曲线小结：1.用两个继电器并联得到具有死去特性的三位置继电特性 2.从图中可

20、以看出输出比较稳定实验二 经典控制理论实验目的以MATLAB及Simulink为工具，对控制系统进行时域、频域及根轨迹分析。实验原理1、 时域分析法是根据系统的微分方程（或传递函数），利用拉普拉斯变换直接解出动态方程，并依据过程曲线及表达式分析系统的性能。时域响应指标如图所示。延迟时间td，指响应曲线第一次达到其终值一半所需要的时间。上升时间tr，指响应曲线从终值10%上升到终值90%所需要的时间；对于有振荡的系统，也可定义为响应从零第一次上升到终值所需要的时间。上升时间是系统响应速度的一种度量。峰值时间tp，指响应超过终值达到第一个峰值所需要的时间。调节时间ts，指响应达到并保持在终值5%（

21、或2%）内所需要的时间。超调量%，指响应的最大偏离量h(tp)与终值h()之差的百分比，即：稳态误差，描述系统稳态性能的一种性能指标。2、 频域分析法通常从频率特性出发对系统进行研究。在工程分析和设计中，通常把频率特性画成一些曲线，从频率特性曲线出发进行研究。这些曲线包括幅频特性和相频特性曲线，幅相频率特性曲线，对数频率特性曲线以及对数幅相曲线等，其中以幅相频率特性曲线，对数频率特性曲线应用最广。对于最小相位系统，幅频特性和相频特性之间存在着唯一的对于关系，故根据对数幅频特性，可以唯一地确定相应的相频特性和传递函数。根据系统的开环频率特性去判断闭环系统的性能，并能较方便地分析系统参量对系统性能

22、的影响，从而指出改善系统性能的途径。3、 根轨迹是求解闭环系统特征根的图解方法。由于控制系统的动态性能是由系统闭环零极点共同决定，控制系统的稳定性由闭环系统极点唯一确定，利用根轨迹确定闭环系统的零极点在s平面的位置，分析控制系统的动态性能。实验内容*1控制系统数学模型的转换自动控制系统计算机仿真教材第4章中的所有例题【例4-1】已知系统的传递函数如下，利用MATLAB建立其相应的传递函数系统模型答案程序：clc,clearnum=5*2 0 3;den=conv(conv(conv(1 0 0,3 1),conv(1 2,1 2),5 0 3 8);printsys(num,den,s)G=t

23、f(num,den)结果：小结：G=tf(num,den)是用来求传递函数的系统模型的函数【例4-2】 已知系统传递函数如下 应用MATLAB函数将其转换为状态方程形式的模型答案：程序：clcclearnum=12 24 0 20;den=2 4 6 2 2;sys=tf(num,den)s=ss(sys) 结果：Transfer function: 12 s3 + 24 s2 + 202 s4 + 4 s3 + 6 s2 + 2 s + 2 a = x1 x2 x3 x4 x1 -2 -1.5 -0.5 -0.5 x2 2 0 0 0 x3 0 1 0 0 x4 0 0 1 0 b = u1

24、 x1 4 x2 0 x3 0 x4 0 c = x1 x2 x3 x4 y1 1.5 1.5 0 1.25 d = u1 y1 0 Continuous-time model.小结：tf(num,den)函数用来求传递函数的状态模型【例4-3】某线性定常系统的状态空间表达式如下，求该系统的传递函数。答案：clcclearA=0 1 1;0 0 1;-10 -17 -8;B=0;0;1;C=5 6 1;D=0;SS=ss(A,B,C,D)G1=tf(SS)结果：Transfer function: s2 + 11 s + 5s3 + 8 s2 + 27 s + 10小结：ss(A,B,C,D)

25、函数是用来求状态方程的传递函数【例4-4】对于例4-3题中的线性定常系统，将其转换为zpk形式 答案：程序：clcclearA=0 1 1;0 0 1;-10 -17 -8;B=0;0;1;C=5 6 1;D=0;z,p,k=ss2zp(A,B,C,D)G1=zpk(z,p,k)结果：z = -10.5249 -0.4751p = -0.4199 -3.7901 + 3.0745i -3.7901 - 3.0745ik = 1.0000 Zero/pole/gain: (s+10.52) (s+0.4751)(s+0.4199) (s2 + 7.58s + 23.82)小结：zpk(z,p,k

26、)是用来求状态方程的zpk形式的函数【例4-5】 某系统的传递函数如下，求它的部分分式形式。答案：程序：clcclearnum=20 10;den=1 15 74 120;R,P,H=residue(num,den)结果R = -55.0000 90.0000 -35.0000P = -6.0000 -5.0000 -4.0000H = 小结：R,P,H=residue(num,den)函数是用来求传递函数的部分分式函数【例4-6】 某连续系统的状态空间表达式如下，采用零阶保持器将其离散化，设采样周期为0.1秒，求离散化的系统方程。 答案：程序：A=0 1 0;0 0 1;-6 -11 -6;

27、B=1 0;2 -1;0 2;C=1 -1 0;2 1 -1;D=0;T=0.1;G=ss(A,B,C,D);Gn=c2d(G,T)结果a = x1 x2 x3 x1 0.9991 0.0984 0.004097 x2 -0.02458 0.9541 0.07382 x3 -0.4429 -0.8366 0.5112 b = u1 u2 x1 0.1099 -0.004672 x2 0.1959 -0.0902 x3 -0.1164 0.1936 c = x1 x2 x3 y1 1 -1 0 y2 2 1 -1 d = u1 u2 y1 0 0 y2 0 0 Sampling time: 0.

28、1Discrete-time model.小结：Gn=c2d(G,T)采用零阶保持器将系统模型离散化，求出采样周期为T离散化的系统方程的函数【例4-7】 已知系统如图所示，利用MATLAB求出系统的状态空间表达式。其中，sys1: sys2: 答案：程序：A1=-9 17;-1 3;B1=0 -1;-1 0;C1=-3 2;-13 18;D1=-1 0;-1 0;num=2;den=1 2;G1=ss(A1,B1,C1,D1);G2=tf(num,den);sys=feedback(G1,G)结果a = x1 x2 x3 x1 -9 17 1 x2 -1 3 0 x3 -26 36 -2 b

29、= u1 u2 x1 0 -1 x2 -1 0 x3 -2 0 c = x1 x2 x3 y1 -3 2 0 y2 -13 18 0 d = u1 u2 y1 -1 0 y2 -1 0 Continuous-time model.小结：sys=feedback(G1,G)是用来求反馈系统的函数【例4-8】 已知系统的状态空间表达式如下，求线性变换，将其变换成能控标准形。答案：程序：A=1 2 -1;0 2 1;1 -3 2;B=0;1;1;C=1 0 1;Qc=ctrb(A,B);A1=poly(A)if rank(Qc)=3 disp(The system is controllable ) Q=Qc*A1(3) A1(2) 1;A1(2)