系统仿真实验报告.docx
- 文档编号:11995847
- 上传时间:2023-06-03
- 格式:DOCX
- 页数:33
- 大小:181.95KB
系统仿真实验报告.docx
《系统仿真实验报告.docx》由会员分享,可在线阅读,更多相关《系统仿真实验报告.docx(33页珍藏版)》请在冰点文库上搜索。
系统仿真实验报告
(此文档为word格式,下载后您可任意编辑修改!
)
系统仿真实验报告
实验一MATLAB中矩阵与多项式的基本运算1
实验二MATLAB绘图命令6
实验三MATLAB程序设计8
实验四MATLAB的符号计算与SIMULINK的使用12
实验五MATLAB在控制系统分析中的应用15
实验六连续系统数字仿真的基本算法26
实验一MATLAB中矩阵与多项式的基本运算
一、实验任务
1.了解MATLAB命令窗口和程序文件的调用。
2.熟悉如下MATLAB的基本运算:
①矩阵的产生、数据的输入、相关元素的显示;
②矩阵的加法、乘法、左除、右除;
③特殊矩阵:
单位矩阵、“1”矩阵、“0”矩阵、对角阵、随机矩阵的产生和运算;
④多项式的运算:
多项式求根、多项式之间的乘除。
二、基本命令训练
1.eye(m)
>>eye(4)
ans=
1000
0100
0010
0001
2.one(n)、ones(m,n)
>>ones(3,4)
ans=
1111
1111
1111
3.zeros(m,n)
>>zeros(2,3)
ans=000
000
4.rand(m,n)
>>rand(3,4)
ans=
0.95010.48600.45650.4447
0.23110.89130.01850.6154
0.60680.76210.82140.7919
5.diag(v)
>>v=[111];
X=diag(v)
X=
100
010
001
>>X=diag(v,1)
X=
0100
0010
0001
0000
>>X=diag(v,2)
X=
00100
00010
00001
00000
00000
6.A\B、A/B、inv(A)*B、B*inv(A)
Matlab提供了两种除法运算:
左除(\)和右除(/)。
一般情况下,x=a\b是方程a*x=b的解,而x=b/a是方程x*a=b的解。
a=[1 2 3;4 2 6;7 4 9];b=[4;1;2];x=a\b
x=
-1.5000
2.0000
0.5000
>>a=[1127;385;436];b=flipud(a)
b=
436
385
1127
>>a/b
ans=
001
010
100
如果a为非奇异矩阵,则a\b和b/a可通过a的逆矩阵与b阵得到:
a\b=inv(a)*b
b/a=b*inv(a)
7.roots(p)
求解多项式x^3-6x^2-72x-27=0的根
>>p=[1-6-72-27];roots(p)
ans=
12.1229
-5.7345
-0.3884
8.Poly
>>r=[123];poly(r)
ans=
1-611-6
创建以之为根的方程为:
x^3-6*x^2+11*x-6=0
>>A=[123;456;780];
P=poly(A)
P=
1.0000-6.0000-72.0000-27.0000
创建特征方程为λ^3-6*λ^2-72*λ-27=0
9.conv、deconv
>>a=[1234];b=[3-5823-7];conv(a,b)
ans=
31736438771-28
>>deconv(a,b)
ans=
0
>>deconv(b,a)
ans=
3-11
>>a=[10-1];b=[1-1];deconv(a,b)
ans=
11
>>conv(a,b)
ans=
1-1-11
10.A*B与A.*B的区别
>>a=[12;34];b=[45;67];c=a*b
c=
1619
3643
>>a.*b
ans=
410
1828
11.who与whos的使用
>>who
Yourvariablesare:
APaansbcr
>>whos
NameSizeBytesClass
A3x372doublearray
P1x432doublearray
a2x232doublearray
ans2x232doublearray
b2x232doublearray
c2x232doublearray
r1x324doublearray
Grandtotalis32elementsusing256bytes
12.disp、size(a)、length(a)的使用
>>d=[12;34];disp(d)
12
34
>>disp('d')
d
>>size(d)
ans=
22
>>length(d)
ans=
2
实验二MATLAB绘图命令
一、实验任务
熟悉MATLAB基本绘图命令,掌握如下绘图方法:
1.坐标系的选择、图形的绘制;
2.图形注解(题目、标号、说明、分格线)的加入;
3.图形线型、符号、颜色的选取。
二、基本命令训练
1.plot2.loglog3.semilogx4.semilogy
5.polar6.title7.xlabel8.ylabel
9.text10.grid11.bar12.stairs
13.contour
三、实验举例
1.t=[0:
pi/360:
2*pi*22/3];
x=93*cos(t)+36*cos(t*4.15);
y=93*sin(t)+36*sin(t*4.15);
plot(y,x),grid;
2.t=0:
0.05:
100;
x=t;y=2*t;z=sin(2*t);
plot3(x,y,z,'b:
')
3.t=0:
pi/20:
2*pi;
y=sin(x);
stairs(x,y)
4.th=[pi/200:
pi/200:
2*pi]';
r=cos(2*th);
polar(th,r),grid
5.th=[0:
pi/10:
2*pi];
x=exp(j*th);
plot(real(x),imag(x),'r*');
grid;
实验三MATLAB程序设计
一、实验任务
1.熟悉MATLAB程序设计的方法和思路;
2.掌握循环、分支语句的编写,学会使用lookfor、help命令。
二、程序举例
1.计算1~1000之内的斐波那契亚数列
f=[1,1];
i=1;
whilef(i)+f(i+1)<1000
f(i+2)=f(i)+f(i+1);
i=i+1;
end
f,i
f=
Columns1through14
1123581321345589144233377
Columns15through16
610987
i=15
2.m=3;
n=4;
fori=1:
m
forj=1:
n
a(i,j)=1/(i+j-1);
end
end
formatrat
a
a=
11/21/31/4
1/21/31/41/5
1/31/41/51/6
3.m=3;
n=4;
fori=1:
m
forj=1:
n
a(i,j)=1/(i+j-1);
end
end
a
a=
11/21/31/4
1/21/31/41/5
1/31/41/51/6
4.x=input('请输入x的值:
');
ifx==10
y=cos(x+1)+sqrt(x*x+1);
else
y=x*sqrt(x+sqrt(x));
end
y
x=input('请输入x的值:
')
请输入x的值:
10
x=10
>>ifx==10
y=cos(x+1)+sqrt(x*x+1);
else
y=x*sqrt(x+sqrt(x));
end
y
y=10.0543
>>x=input('请输入x的值:
')
请输入x的值:
8
x=8
>>ifx==10
y=cos(x+1)+sqrt(x*x+1);
else
y=x*sqrt(x+sqrt(x));
end
y
y=26.3253
5.去掉多项式或数列开头的零项
p=[0001302009];
fori=1:
length(p),ifp
(1)==0,p=p(2:
length(p));
end;
end;
p
p=
1302009
6.建立MATLAB的函数文件,程序代码如下,以文件名ex2_4.m存盘
functionf=ffibno(n)
%ffibno计算斐波那契亚数列的函数文件
%n可取任意自然数
%程序如下
f=[1,1];
i=1;
whilef(i)+f(i+1) f(i+2)=f(i)+f(i+1); i=i+1; end >>exe2_4(200) ans= 1123581321345589144 >>lookforffibno exe2_4.m: %ffibno计算斐波那契亚数列的函数文件 >>helpexe2_4 ffibno计算斐波那契亚数列的函数文件 n可取任意自然数 三、程序设计题 用一个MATLAB语言编写一个程序: 输入一个自然数,判断它是否是素数,如果是,输出“Itisoneprime”,如果不是,输出“Itisnotoneprime.”。 要求通过调用子函数实现。 最好能具有如下功能: ①设计较好的人机对话界面,程序中含有提示性的输入输出语句。 ②能实现循环操作,由操作者输入相关命令来控制是否继续进行素数的判断。 如果操作者希望停止这种判断,则可以退出程序。 ③如果所输入的自然数是一个合数,除了给出其不是素数的结论外,还应给出至少一种其因数分解形式。 例: 输入6,因为6不是素数。 则程序中除了有“Itisnotoneprime”的结论外,还应有: “6=2*3”的说明。 >> close all; c=input('是否进行素数运算 1为是 0为否: ') while c==1 a=input('请输入一个自然数: '); if a==1 disp('1既不是质数也不是合数'); else if factor(a)==a disp('It IS one prime'); else disp('It is not one prime'); b=factor(a); fprintf('%3d=',a) for j=1: (length(b)-1) fprintf('%3d *',b(j)) end fprintf('%3d \n *',b(length(b))) end end c=input('是否进行素数运算 1为是 0为否: ') end 是否进行素数运算 1为是 0为否: 1 c = 1 请输入一个自然数: 1 1既不是质数也不是合数 是否进行素数运算 1为是 0为否: 1 c = 1 请输入一个自然数: 6 It is not one prime 6= 2 * 3 实验四MATLAB的符号计算与SIMULINK的使用 一、实验任务 1.掌握MATLAB符号计算的特点和常用基本命令; 2.掌握SIMULINK的使用。 二、程序举例 1.求矩阵对应的行列式和特征根 a=sym('[a11a12;a21a22]'); da=det(a) ea=eig(a) da= a11*a22-a12*a21 ea= 1/2*a11+1/2*a22+1/2*(a11^2-2*a11*a22+a22^2+4*a12*a21)^(1/2) 1/2*a11+1/2*a22-1/2*(a11^2-2*a11*a22+a22^2+4*a12*a21)^(1/2) 2.求方程的解(包括精确解和一定精度的解) r1=solve('x^2-x-1') rv=vpa(r1) rv4=vpa(r1,4) rv20=vpa(r1,20) r1= 1/2*5^(1/2)+1/2 1/2-1/2*5^(1/2) rv= 1. -. rv4= 1.618 -.6180 rv20= 1. -. 3.a=sym('a');b=sym('b');c=sym('c');d=sym('d');%定义4个符号变量 w=10;x=5;y=-8;z=11;%定义4个数值变量 A=[a,b;c,d]%建立符号矩阵A B=[w,x;y,z]%建立数值矩阵B det(A)%计算符号矩阵A的行列式 det(B)%计算数值矩阵B的行列式 A= [a,b] [c,d] B= 105 -811 ans= a*d-b*c ans= 150 4.symsxy; s=(-7*x^2-8*y^2)*(-x^2+3*y^2); expand(s)%对s展开 collect(s,x)%对s按变量x合并同类项(无同类项) factor(ans)%对ans分解因式 ans= 7*x^4-13*x^2*y^2-24*y^4 ans= 7*x^4-13*x^2*y^2-24*y^4 ans= (8*y^2+7*x^2)*(x^2-3*y^2) 5.对方程AX=b求解 A=[34,8,4;3,34,3;3,6,8]; b=[4;6;2]; X=linsolve(A,b)%调用linsolve函数求解 A\b%用另一种方法求解 X= 0.0675 0.1614 0.1037 ans= 0.0675 0.1614 0.1037 6.对方程组求解 a11*x1+a12*x2+a13*x3=b1 a21*x1+a22*x2+a23*x3=b2 a31*x1+a32*x2+a33*x3=b3 symsa11a12a13a21a22a23a31a32a33b1b2b3; A=[a11,a12,a13;a21,a22,a23;a31,a32,a33]; b=[b1;b2;b3]; X=linsolve(A,b)%调用linsolve函数求的解 XX=A\b%用左除运算求解 XX=(a12*a23*b3-a12*b2*a33+a13*a32*b2-a13*b3*a22+b1*a33*a22-b1*a32*a23)/(a11*a33*a22-a33*a12*a21-a31*a13*a22+a31*a12*a23+a32*a21*a13-a11*a32*a23) -(a11*a23*b3-a11*b2*a33-a21*a13*b3+b2*a31*a13-a23*a31*b1+a21*b1*a33)/(a11*a33*a22-a33*a12*a21-a31*a13*a22+a31*a12*a23+a32*a21*a13-a11*a32*a23) (-a31*b1*a22-a11*a32*b2+a11*b3*a22-b3*a12*a21+a31*a12*b2+a32*a21*b1)/(a11*a33*a22-a33*a12*a21-a31*a13*a22+a31*a12*a23+a32*a21*a13-a11*a32*a23) 7.symsabtxyz; f=sqrt(1+exp(x)); diff(f)%未指定求导变量和阶数,按缺省规则处理 f=x*cos(x); diff(f,x,2)%求f对x的二阶导数 diff(f,x,3)%求f对x的三阶导数 f1=a*cos(t);f2=b*sin(t); diff(f2)/diff(f1)%按参数方程求导公式求y对x的导数 ans= 1/2/(1+exp(x))^(1/2)*exp(x) ans= -2*sin(x)-x*cos(x) ans= -3*cos(x)+x*sin(x) ans= -b*cos(t)/a/sin(t) 三、SIMULINK的使用 1.在命令窗口中输入: simulink(回车)得到如下simulink模块: 2.双击打开各模块,选择合适子模块构造控制系统,例如: 3.双击各子模块可修改其参数,选择Simulation菜单下的start命令运行仿真,在示波器(Scope)中观察结果。 实验五MATLAB在控制系统分析中的应用 一、实验任务 1.掌握MATLAB在控制系统时间响应分析中的应用; 2.掌握MATLAB在系统根轨迹分析中的应用; 3.掌握MATLAB控制系统频率分析中的应用; 4.掌握MATLAB在控制系统稳定性分析中的应用 二、基本命令 1.step2.impulse3.initial4.lsim5.rlocfind 6.bode7.margin8.nyquist9.Nichols10.cloop 三、程序举例 1.求下面系统的单位阶跃响应 %程序如下: num=[4];den=[1,1,4]; step(num,den) [y,x,t]=step(num,den); tp=spline(y,t,max(y))%计算峰值时间 max(y)%计算峰值 tp= 1.6062 ans= 1.4441 2.求如下系统的单位阶跃响应 %程序如下: a=[0,1;-6,-5];b=[0;1];c=[1,0];d=0; [y,x]=step(a,b,c,d); plot(y) 3.求下面系统的单位脉冲响应: %程序如下: num=[4];den=[1,1,4]; impulse(num,den) 4.已知二阶系统的状态方程为: 求系统的零输入响应和脉冲响应。 %程序如下: a=[0,1;-10,-2];b=[0;1]; c=[1,0];d=[0]; x0=[1,0]; subplot(1,2,1);initial(a,b,c,d,x0) subplot(1,2,2);impulse(a,b,c,d) 5: 系统传递函数为: 输入正弦信号时,观察输出信号的相位差。 %程序如下: num=[1];den=[1,1]; t=0: 0.01: 10; u=sin(2*t); holdon plot(t,u,‘r: ’) lsim(num,den,u,t) 6.有一二阶系统,求出周期为4秒的方波的输出响应 %程序如下: num=[251]; den=[123]; t=(0: .1: 10); period=4; u=(rem(t,period)>=period./2);%看rem函数功能 u=double(u) lsim(num,den,u,t) 7.已知开环系统传递函数,绘制系统的根轨迹,并分析其稳定性 %程序如下: num=[12]; den1=[143]; den=conv(den1,den1); figure (1) rlocus(num,den) [k,p]=rlocfind(num,den) figure (2) k=55; num1=k*[12]; den=[143]; den1=conv(den,den); [num,den]=cloop(num1,den1,-1); impulse(num,den) title('impulseresponse(k=55)') figure(3) k=56; num1=k*[12]; den=[143]; den1=conv(den,den); [num,den]=cloop(num1,den1,-1); impulse(num,den) title('impulseresponse(k=56)') Selectapointinthegraphicswindow selected_point= -1.5687+0.6708i k= 2.4232 p= -3.7489 -2.3284 -0.9613+0.8137i -0.9613-0.8137i 8.作如下系统的bode图 %程序如下: n=[1,1];d=[1,4,11,7]; bode(n,d) 9.系统传函如下 求有理传函的频率响应,然后在同一张图上绘出以四阶伯德近似表示的系统频率响应 %程序如下: num=[1];den=conv([12],conv([12],[12])); w=logspace(-1,2);t=0.5; [m1,p1]=bode(num,den,2); p1=p1-t*w'*180/pi; [n2,d2]=pade(t,4); numt=conv(n2,num);dent=(conv(den,d2)); [m2,p2]=bode(numt,dent,w); subplot(2,1,1);semilogx(w,20*log10(m1),w,20*log10(m2),'g--'); gridon;title('bodeplot');xlabel('frequency');ylabel('gain'); subplot(2,1,2);semilogx(w,p1,w,p2,'g--');gridon; xlabel('frequency');ylabel('phase'); 10.已知系统模型为 求它的幅值裕度和相角裕度 %程序如下: n=[3.5];d=[1232]; [Gm,Pm,Wcg,Wcp]=margin(n,d) Gm= 1.1433 Pm= 7.1688 Wcg= 1.7323 Wcp
- 配套讲稿:
如PPT文件的首页显示word图标,表示该PPT已包含配套word讲稿。双击word图标可打开word文档。
- 特殊限制:
部分文档作品中含有的国旗、国徽等图片,仅作为作品整体效果示例展示,禁止商用。设计者仅对作品中独创性部分享有著作权。
- 关 键 词:
- 系统 仿真 实验 报告