系统仿真实验报告.docx
- 文档编号:10203510
- 上传时间:2023-02-09
- 格式:DOCX
- 页数:44
- 大小:240.85KB
系统仿真实验报告.docx
《系统仿真实验报告.docx》由会员分享,可在线阅读,更多相关《系统仿真实验报告.docx(44页珍藏版)》请在冰豆网上搜索。
系统仿真实验报告
中南大学系
统仿真技术实验报告
作者:
黄伟
学号:
0909122327
班级:
自动化1205
指导老师:
胡杨
时间:
2014-6-5
一、目的要求
1.了解MATLAB命令窗口和程序文件的调用。
2.熟悉如下MATLAB的基本运算:
①矩阵的产生、数据的输入、相关元素的显示;
②矩阵的加法、乘法、左除、右除;
③特殊矩阵:
单位矩阵、“1”矩阵、“0”矩阵、对角阵、随机矩阵的产生和运算;
④多项式的运算:
多项式求根、多项式之间的乘除。
3熟悉MATLAB基本绘图命令,掌握如下绘图方法:
①坐标系的选择、图形的绘制;
②图形注解(题目、标号、说明、分格线)的加入;
③图形线型、符号、颜色的选取。
4.熟悉MATLAB程序设计的方法和思路;
5.掌握循环、分支语句的编写,学会使用lookfor、help命令。
6.参照已知程序,改动程序中的参数和输入量,验证输出结果。
7.使用lookfor、help命令,验证输出结果
8.掌握MATLAB符号计算的特点和常用基本命令;
9.掌握SIMULINK的使用。
10.掌握MATLAB在控制系统时间响应分析中的应用;
11.掌握MATLAB在系统根轨迹分析中的应用;
12.掌握MATLAB控制系统频率分析中的应用;
13.掌握MATLAB在控制系统稳定性分析中的应用
14.理解欧拉法和龙格-库塔法的基本思想;
15.理解数值积分算法的计算精度、速度、稳定性与步长的关系;
二、仪器设备
计算机1台
Matlab软件1个
三、内容步骤及数据处理
I.MATLAB中矩阵与多项式的基本运算
1.eye(m)
eye(3)
ans=
100
010
001
2.one(n)、ones(m,n)
one
(1)
ans=
1
ones(2,3)
ans=
111
111
3.zeros(m,n)
zeros(3,4)
ans=
0000
0000
0000
4.rand(m,n)
5.diag(v)
p=rand(4,4);%3,4题%随即数列
>>disp(p)
0.64770.74470.36850.9294
0.45090.18900.62560.7757
0.54700.68680.78020.4868
0.29630.18350.08110.4359
>>diag(p)%提取对角线数据
ans=
0.6477
0.1890
0.7802
0.43596.A\B、A/B、inv(A)*B、B*inv(A)
A=rand(4,4);
>>B=rand(4,4);
>>A\B
ans=
-0.17642.15690.04780.2260
0.09261.04880.65841.2567
1.4410-2.3537-0.1777-1.2059
0.1421-1.88930.42390.0172
A/B
ans=
-0.44100.39650.9073-0.5299
-1.56360.94080.54770.9849
3.24340.5827-1.3008-1.3753
1.78290.5311-0.3391-0.8652
>>inv(A)*B
ans=
-0.17642.15690.04780.2260
0.09261.04880.65841.2567
1.4410-2.3537-0.1777-1.2059
0.1421-1.88930.42390.0172
>>B*inv(A)
ans=
-1.2220-0.2555-1.84303.3873
0.45350.76600.8263-0.7193
-0.8123-0.4002-2.58684.1538
-1.92150.1006-2.27683.7547
7.roots(p)
P=[121];%求多项式根
>>roots(P)
ans=
-1
-1
8.Poly
P=[1234;5678;9012;1111];%求特征多项式
>>Poly(P)
ans=
1.0000-9.0000-30.000030.0000-0.00009.conv、deconv
10.A*B与A.*B的区别
11.who与whos的使用
P=[1234;5678;9012;1111];
A=[1111;0000;1234;1222];
>>A*P%向量
ans=
1691215
0000
42182430
31162126
>>A.*P%对应项相乘
ans=
1234
0000
9038
1222
Who%显示参数
Yourvariablesare:
ABPansnp
Whos%显示参数及属性
NameSizeBytesClassAttributes
A4x4128double
B4x4128double
P4x4128double
ans4x4128double
n1x18double
p3x372double
12.disp、size(a)、length(a)的使用
P=[1234;5678];
>>disp(P)
1234
5678
>>size(P)
ans=
24
>>length(P)
ans=
4
II.MATLAB绘图命令
1.plot
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,’r’)%红色
2.loglog
t=[0123456789101112];
x=10.^t;
y=100.^t;1
>>loglog(x,y,’g’)%绿色
3.semilogx
4.Semilogy
x=0:
.1:
10;semilogx(10.^x,x)
x=0:
.1:
10;semilogy(x,10.^x)
5.polar
th=[pi/200:
pi/200:
2*pi]';
r=cos(2*th);
polar(th,r),grid
6.title
th=[pi/200:
pi/200:
2*pi]';
r=cos(2*th);
polar(th,r);grid;title('极坐标')
7.xlabel
x=0:
0.1:
10;semilogx(10.^x,x);xlabel('x轴对数变换')
8.Ylabel
x=0:
.1:
10;semilogy(x,10.^x);ylabel('y轴对数变换')
9.text
10.grid
x=0:
.1:
10;semilogy(x,10.^x);text(2,100,'(2,100)');grid%9,10题
11.bar
x=[456787654];
bar(x)%绘制直方图
12.Stairs
x=0:
0.1:
10;%阶梯
y=x;
stairs(x,y)
13.contour
x=1:
1:
10;%等高线
y=1:
1:
10;
[xx,yy]=meshgrid(y,x);
z=rand(10,10);
contour(xx,yy,z,15)
III.MATLAB程序设计
1.计算1~1000之内的斐波那契亚数列
f=[1,1];
i=1;
whilef(i)+f(i+1)<10000
f(i+2)=f(i)+f(i+1);
i=i+1;
end
>>f,i
f=
Columns1through13
1123581321345589144233
Columns14through20
3776109871597258441816765
i=
19
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=
1.00000.50000.33330.2500
0.50000.33330.25000.2000
0.33330.25000.20000.1667
请比较程序2与程序3的区别
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的值:
');
if(x==0)
disp('分母为零不合法,请重新输入!
');
else
y=(3^x)/x+cos(x);
end
y
请输入x的值:
5
y=
48.8837
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 输入完毕后在MATLAB的命令窗口0输入ex2_4(200),得到运行结果。 ans= 1123581321345589144 IV.MATLAB的符号计算与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*5^(1/2)+1/2 rv= 1.6180339887498948482045868343656 -.61803398874989484820458683436560 rv4= 1.618 -.6180 rv20= 1.6180339887498948482 -.61803398874989484820 3. >>a=sym('a');b=sym('b');c=sym('c');d=sym('d'); w=10;x=5;y=-8;z=11; A=[a,b;c,d] B=[w,x;y,z] det(A) det(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. >>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. symsa11a12a13a21a22a23a31a32a33b1b2b3; A=[a11,a12,a13;a21,a22,a23;a31,a32,a33]; b=[b1;b2;b3]; XX=A\b%用左除运算求解 (X=linsolve(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*a21*a12-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*a21*a12-a31*a13*a22+a31*a12*a23+a32*a21*a13-a11*a32*a23) (-a31*b1*a22-a11*a32*b2+a11*b3*a22-b3*a21*a12+a31*a12*b2+a32*a21*b1)/(a11*a33*a22-a33*a21*a12-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) 7.simulink仿真 V.MATLAB在控制系统分析中的应用 1.求下面系统的单位阶跃响应 %程序如下: num=[4];den=[1,1,4]; step(num,den) [y,x,t]=step(num,den); tp=spline(y,t,max(y))%计算峰值时间 max(y)%计算峰值 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函数功能 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)') 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= 1.6541 11.二阶系统为: 令wn=1,分别作出ξ=2,1,0.707,0.5时的nyquist曲线。 %程序如下: n=[1]; d1=[1,4,1];d2=[1,2,1];d3=[1,1.414,1];d4=[1,1,1]; nyquist(n,d1); holdon nyquist(n,d2);nyquist(n,d3);nyquist(n,d4); 12.已知系统的开环传递函数为 绘制系统的Nyqusit图,并讨论系统的稳定性 %程序如下: G=tf(1000,conv([1,3,2],[1,5])); nyquist(G);axis('square') 13.分别由w的自动变量和人工变量作下列系统的nyquist曲线: %程序如下: n=[1];d=[1,1,0];nyquist(n,d);%自动变量 n=[1];d=[1,1,0];w=[0.5: 0.1: 3];nyquist(n,d,w);%人工变量 14.一多环系统,其结构图如下,使用Nyquist频率曲线判断系统的稳定性。 %程序如下: k1=16.7/0.0125;z1=[0]; p1=[-1.25-4-16]; [num1,den1]=zp2tf(z1,p1,k1); [num,den]=cloop(num1,den1); [z,p,k]=tf2zp(num,den);p figure (1) nyquist(num,den) figure (2) [num2,den2]=cloop(num,den); impulse(num2,den2); p= -10.5969+36.2148i -10.5969-36.2148i -0.0562 15.已知系统为: 作该系统的nichols曲线。 %程序如下: n=[1];d=[1,1,0]; nichols(n,d); 16.已知系统的开环传递函数为: 当k=2时,分别作nichols曲线和波特图。 %程序如下: num=1; den=conv(conv([10],[11]),[0.51]); subplot(1,2,1); nichols(num,den);grid;%nichols曲线 su
- 配套讲稿:
如PPT文件的首页显示word图标,表示该PPT已包含配套word讲稿。双击word图标可打开word文档。
- 特殊限制:
部分文档作品中含有的国旗、国徽等图片,仅作为作品整体效果示例展示,禁止商用。设计者仅对作品中独创性部分享有著作权。
- 关 键 词:
- 系统 仿真 实验 报告