控制系统的Matlab仿真与设计课后答案

合集下载

控制系统的MATLAB仿真与设计课后答案

控制系统的MATLAB仿真与设计课后答案

控制系统的MATLAB 仿真与设计课后答案第二章1>>x=[15 22 33 94 85 77 60]>>x(6)>>x([1 3 5])>>x(4:end)>>x(find(x>70))2>>T=[1 -2 3 -4 2 -3] ;>>n=length(T);>>TT=T';>>for k=n-1:-1:0>>B(:,n-k)=TT.^k;>>end>>B>>test=vander(T)3>>A=zeros(2,5);>>A(:)=-4:5>>L=abs(A)>3>>islogical(L)>>X=A(L)4>>A=[4,15,-45,10,6;56,0,17,-45,0] >>find(A>=10&A<=20)5>>p1=conv([1,0,2],conv([1,4],[1,1]));>>p2=[1 0 1 1];>>[q,r]=deconv(p1,p2);>>cq='商多项式为 '; cr='余多项式为 ';>>disp([cq,poly2str(q,'s')]),disp([cr,poly2str(r,'s')]) 6>>A=[11 12 13;14 15 16;17 18 19];>>PA=poly(A)>>PPA=poly2str(PA,'s')第三章1>>n=(-10:10)';>>y=abs(n);>>plot(n,y,'r.','MarkerSize',20)>>axis equal>>grid on>>xlabel('n')2>>x=0:pi/100:2*pi;>>y=2*exp(-0.5*x).*sin(2*pi*x);>>plot(x,y),grid on;3>>t=0:pi/50:2*pi;>>x=8*cos(t);>>y=4*sqrt(2)*sin(t);>>z=-4*sqrt(2)*sin(t);>>plot3(x,y,z,'p');>>title('Line in 3-D Space');>>text(0,0,0,'origin');>>xlabel('X'),ylable('Y'),zlable('Z');grid;4>>theta=0:0.01:2*pi;>>rho=sin(2*theta).*cos(2*theta); >>polar(theta,rho,'k');5>>[x,y,z]=sphere(20);>>z1=z;>>z1(:,1:4)=NaN;>>c1=ones(size(z1));>>surf(3*x,3*y,3*z1,c1);>>hold on>>z2=z;>>c2=2*ones(size(z2));>>c2(:,1:4)=3*ones(size(c2(:,1:4))); >>surf(1.5*x,1.5*y,1.5*z2,c2);>>colormap([0,1,0;0.5,0,0;1,0,0]); >>grid on>>hold off第四章1>>for m=100:999m1=fix(m/100);m2=rem(fix(m/10),10);m3=rem(m,10);if m==m1*m1*m1+m2*m2*m2+m3*m3*m3 disp(m)endend2M文件:function[s,p]=fcircle(r)s=pi*r*r;p=2*pi*r;主程序:[s,p]=fcircle(10)3>>y=0;n=100;for i=1:ny=y+1/i/i;end>>y4 M文件:function f=factor(n)if n<=1f=1;elsef=factor(n-1)*n; end主程序:>>s=0;for i=1:5s=s+factor(i);end>>s5>>sum=0;i=1;while sum<sum=sum+i;i=i+1;end;>>n=i-26for循环M文件:function k=jcsum(n) k=0;for i=0:nk=k+2^i;end主程序:>>jcsum(63)While循环M文件:function k=jcsum1(n)k=0;i=0;while i<=nk=k+2^i;i=i+1;end主程序:>>jcsum1(63)第五章1>>A=[2,1,-5,1;1,-5,0,7;0,2,1,-1;1,6,-1,-4]; >>b=[13,-9,6,0]';>>x=A\b2M文件:function f=fxyz(u)x=u(1);y=u(2);z=u(3);f=x+y.^2./x/4+z.^2./y+2./z;主程序:[U,fmin]=fminsearch('fxyz',[0.5,0.5,0.5]) 3>>X=linspace(0,2*pi,50);。

控制系统仿真答案

控制系统仿真答案

控制系统仿真答案一.选择题二.名词解释及简答题1.系统的三个属性是什么?请解释其具体含义。

答:系统的三个属性:整体性、相关性、隶属性。

整体性:各部分(子系统)不能随意分割。

相关性:各部分(子系统) 以一定的规律或方式相联系,由此决定了其特有的性能。

隶属性:不能清楚的分出系统“内部”与“外部”,常常需要根据研究的问题来确定哪些属于系统的内部因素,哪些属于外部环境,其界限也是随不同的研究目的而变化,将这一特性称之为隶属性。

2.试利用图形表示仿真的基本内容以及其相互关系。

答:3.简述传统的设计流程以及基于虚拟样机的设计流程的区别(可用图形表示)。

答:三.判断题,正确的在括号内打“√”,错误的打“╳”,并改正错误结论重新阐述。

1.(√)为了限制所研究问题涉及的范围,一般用系统边界把被研究的系统与系统环境区分开来。

2.(╳)模型按数学模型的形式分为:物理仿真、数学仿真、数学-物理混合仿真或半实物仿真。

改1:按模型的性质分为:物理仿真、数学仿真、数学-物理混合仿真或半实物仿真。

改2:按数学模型的形式分为:连续系统仿真、离散事件系统仿真、离散—连续系统仿真。

3.(√)series函数可以将两个系统按串联方式连接,它即适合于连续时间系统,也适合于离散时间系统。

4.(√)S-函数为Simulink的“系统”函数,它是能够响应Simulink求解器命令的函数,采用非图形化的方法实现一个动态系统。

5.(╳)Matlab在执行运算符的优先级的时候,逻辑运算符“<”的优先级要高于矩阵乘法运算“*”的优先级。

改:Matlab在执行运算符的优先级的时候,逻辑运算符“<”的优先级要低于矩阵乘法运算“*”的优先级。

四、2、main(){float x,money;Scanf(“%f/n”,&x);If(x<50) then money=2.5*xElse if(x>=50 and x<100)Else if (x>=100)Then money=2.5*x*(1-0.2) ;Printf(“money=%f”,money);}五.综合题1.编程实现以下图形绘制。

Matlab与控制系统仿真部分习题答案

Matlab与控制系统仿真部分习题答案

【4.2】程序:num=[5,0];den=conv([1,1],conv([1,2],[1,3])); [numc,denc]=cloop(num,den);[z,p,k]=tf2zp(numc,denc);[A,B,C,D]=tf2ss(numc,denc);g_zp=zpk(z,p,k)g_tf=tf(numc,denc)g_ss=ss(A,B,C,D)运行结果:Zero/pole/gain:5 s----------------------------------(s+0.4432) (s^2 + 5.557s + 13.54)Transfer function:5 s----------------------s^3 + 6 s^2 + 16 s + 6a =x1 x2 x3x1 -6 -16 -6x2 1 0 0x3 0 1 0b =u1x1 1x2 0x3 0c =x1 x2 x3y1 0 5 0d =u1y1 0【4.3】程序:A=[0 0 0 -1;1 0 0 -2;0 1 0 -3;0 0 1 -4]; B=[0;0;0;1];C=[1 0 0 0];g_ss=ss(A,B,C,D)[num,den]=ss2tf(A,B,C,D);g_tf=tf(num,den)[z,p,k]=ss2zp(A,B,C,D);g_zpk=zpk(z,p,k)运行结果:a =x1 x2 x3 x4x1 0 0 0 -1x2 1 0 0 -2x3 0 1 0 -3x4 0 0 1 -4b =u1x1 0x2 0x3 0x4 1c =x1 x2 x3 x4y1 1 0 0 0d =u1y1 0Continuous-time model.Transfer function:-3.109e-015 s^3 - s^2 - 3.331e-015 s - 4.441e-016 -------------------------------------------------s^4 + 4 s^3 + 3 s^2 + 2 s + 1Zero/pole/gain:- s^2----------------------------------------------(s+0.6724) (s+3.234) (s^2 + 0.0936s + 0.4599)【5.1】(1)程序num=[0,10];den=conv([1,0],[1,7,17]); [numc,denc]=cloop(num,den,-1); G=tf(numc,denc)[y,t]=step(G);plot(t,y,'b-')C=dcgain(G);n=1;while y(n)<0.1*Cn=n+1;endm=1;while y(m)<0.9*Cm=m+1;endrisetime=t(m)-t(n)[Y,k]=max(y); percentovershoot=100*(Y-C)/Ci=length(t);while(y(i)>0.98*C)&(y(i)<1.02*C) i=i-1;endsettlingtime=t(i)运行结果:Transfer function:10-----------------------s^3 + 7 s^2 + 17 s + 10risetime =2.7312percentovershoot =-0.4399settlingtime =5.1372图:0123456700.10.20.30.40.50.60.70.80.91(2)程序k=[10,100,1000];t=linspace(1,20,200);num=1;den=conv([1,0],[1,7,17]);for j=1:3;s1=tf(num*k(j),den);sys=feedback(s1,1)y(:,j)=step(sys,t);endplot(t,y(:,1),'r',t,y(:,2),'b',t,y(:,3),'g')gtext('k=10');gtext('k=100');gtext('k=1000') 运行结果:Transfer function:10-----------------------s^3 + 7 s^2 + 17 s + 10Transfer function:100------------------------s^3 + 7 s^2 + 17 s + 100Transfer function:1000-------------------------s^3 + 7 s^2 + 17 s + 1000图:024681012141618200.20.40.60.811.21.41.61.8图:02468101214161820-3-2-1123422【6.1】程序:(1)num1=[1,1];den1=conv([1,0,0],conv([1,2],[1,4]));sys1=tf(num1,den1)rlocus(sys1)运行结果:-12-10-8-6-4-2024-8-6-4-202468Root LocusReal Axis I m a g i n a r y A x i s(2)num2=[1,1];den2=conv([1,0],conv([1,-1],[1,4,16]));sys2=tf(num2,den2)rlocus(sys2)运行结果:-10-8-6-4-2024-8-6-4-202468Root LocusReal Axis I m a g i n a r y A x i s(3)num3=[1,8];den3=conv([1,0,0],conv([1,3],conv([1,5],conv([1,7],[1,15])))); sys3=tf(num3,den3)rlocus(sys3)运行结果:-30-25-20-15-10-5051015-20-15-10-505101520Root LocusReal Axis I m a g i n a r y A x i s【6.3】程序:num=[1,2];den=conv([1,0],conv([1,4],conv([1,8],[1,2,5])));sys=tf(num,den)rlocus(sys)[k,poles]=rlocfind(sys)运行结果:Transfer function:s + 2---------------------------------------s^5 + 14 s^4 + 61 s^3 + 124 s^2 + 160 sSelect a point in the graphics windowselected_point =0.0296 + 2.2826i k =135.8815poles =-7.3248-5.41040.0145 + 2.3021i0.0145 - 2.3021i -1.2939图:-20-15-10-5051015-15-10-551015Root LocusReal Axis I m a g i n a r y A x i s【7.3】程序(1)画波特图num=[50];den=conv([1,0],conv([1,10],[3,1]));sys=tf(num,den)sys1=feedback(sys,1)bode(sys)grid图(1)-150-100-50050100M a g n i t u d e (d B)10-210-1100101102103-270-225-180-135-90P h a s e (d e g )Bode DiagramFrequency (rad/sec)程序(2)画奈奎斯特图num=[50];den=conv([1,0],conv([1,10],[3,1]));sys=tf(num,den)sys1=feedback(sys,1)nyquist(sys)grid图(2)-16-14-12-10-8-6-4-20-300-200-100100200300Nyquist DiagramReal Axis I m a g i n a r y A x i s程序(3)画零极点图num=[50];den=conv([1,0],conv([1,10],[3,1]));sys=tf(num,den)sys1=feedback(sys,1)pzmap(sys1)gird图(3)P ole-Zero MapReal Axis I m a g i n a r y A x i s -12-10-8-6-4-20-1.5-1-0.50.511.5程序(4)计算相角裕量和幅值裕量num=[50];den=conv([1,0],conv([1,10],[3,1]));sys=tf(num,den)sys1=feedback(sys,1)[gm,pm,wcg,wcp]=margin(sys)运行结果Transfer function:50---------------------3 s^3 + 31 s^2 + 10 sTransfer function:50--------------------------3 s^3 + 31 s^2 + 10 s + 50gm =2.0667pm =7.5615wcg =1.8257wcp =1.2645程序(5)绘制阶跃响应曲线num=[50];den=conv([1,0],conv([1,10],[3,1])); sys=tf(num,den)sys1=feedback(sys,1)step(sys1)图(5)00.20.40.60.811.21.41.61.82Step ResponseTime (sec)A m p l i t u d e【7.4】程序如下:num=[300];den=conv([1,0,0],conv([0.2,1],[0.02,1]));sys=tf(num,den)margin(sys)grid波特图如下:-150-100-50050100M a g n i t u d e (d B )10-1100101102103-360-315-270-225-180P h a s e (d e g )Bode DiagramGm = Inf , P m = -78 deg (at 11 rad/sec)Frequency (rad/sec)【9.3】程序:A=[-2 2 -1;0 -2 0;1 -4 0];B=[0;0;1];C=[1,0,0];D=0;M=ctrb(A,B)m=rank(M)if m==3;disp('系统可控')elsedisp('系统不可控')endN=obsv(A,C)n=rank(N)if n==3;disp('系统可观')elsedisp('系统不可观') endsys=ss(A,B,C,D) [num,den]=ss2tf(A,B,C,D) sys1=tf(num,den)[z,p,k]=ss2zp(A,B,C,D)运行结果:M =0 -1 20 0 01 0 -1m =2系统不可控N =1 0 0-2 2 -13 -4 2n =2系统不可观a =x1 x2 x3x1 -2 2 -1x2 0 -2 0x3 1 -4 0b =u1x1 0x2 0x3 1c =x1 x2 x3y1 1 0 0d =u1y1 0 Continuous-time model.0 0 -1 -2den =1 4 5 2Transfer function:-s - 2---------------------s^3 + 4 s^2 + 5 s + 2z =-2p =-1-1-2k = -1【10.1】(1)程序:A=[0,1,0,0;0,5,0,0;0,0,-7,0;0,0,0,-8]; B=[0;1;0;1];C=[1,2,3,4];D=zeros(1,1);G_ss=ss(A,B,C,D)运行结果:a =x1 x2 x3 x4x1 0 1 0 0x2 0 5 0 0x3 0 0 -7 0x4 0 0 0 -8u1x1 0x2 1x3 0x4 1c =x1 x2 x3 x4y1 1 2 3 4d =u1y1 0(2):程序:[num1,den1]=ss2tf(A,B,C,D); p=roots(den1)i=0;for k=1:1:length(p)if real(p(k))>0i=i+1;endendif i>0disp('系统不稳定');elsedisp('系统稳定');end运行结果:p =5.0000-8.0000-7.0000系统不稳定(3)(4)程序:AA=[0,1,0;0,5,0;0,0,-8];BB=[0;1;1];P=[-1,-2,-8];K=acker(AA,BB,P);i=4;K(4)=0;Kpp=eig(A-B*K)sys1=tf(num1,den1);[y1,t]=step(sys1);plot(t,y1)hold onA_feedback=A-B*K;[num2,den2]=ss2tf(A_feedback,B,C,D); sys2=tf(num2,den2);[y2,t]=step(sys2);plot(t,y2,'r')gridgtext('反馈前')gtext('反馈后')运行结果:K =2 8 0 0pp =-8-2-1-7图形:01234560123456【13.1】程序:A=[0,1;0,0];B=[0;1];C=[1,0];D=zeros(1,1);G_ss=ss(A,B,C,D)M=ctrb(A,B);if rank(M)==2disp('系统完全能控'); elsedisp('系统不完全能控'); endS=[1,0];N=obsv(A,S);if rank(N)==2disp('(A,S)可观测'); elsedisp('(A,S)不可观测'); endR=1;Q=[1,0;0,0];[K,P,E]=Lqr(A,B,Q,R)A_new=A-B*K;G_new=ss(A_new,B,C,D);t=linspace(0,5,100)';y1=step(G_ss,t);y2=step(G_new,t);plot(t,y1,'r:',t,y2,'b-')gridgtext('反馈前')gtext('反馈后')运行结果:a =x1 x2x1 0 1x2 0 0b =u1x1 0x2 1c =x1 x2y1 1 0d =u1y1 0Continuous-time model. 系统完全能控(A,S)可观测K =1.0000 1.4142P =1.4142 1.00001.0000 1.4142E =-0.7071 + 0.7071i-0.7071 - 0.7071i图形:00.51 1.52 2.53 3.54 4.5502468101214。

MATLAB语言:Simulink系统仿真习题与答案

MATLAB语言:Simulink系统仿真习题与答案

一、单选题1、将模块连接好之后,如果要分出一根连线,操作方法是()。

A.把鼠标指针移到分支点的位置,按住鼠标左键拖曳到目标模块的输入端B.双击分支点的位置,按住鼠标左键拖曳到目标模块的输入端C.把鼠标指针移到分支点的位置,按下Ctrl键并按住鼠标拖曳到目标模块的输入端D.把鼠标指针移到分支点的位置,按下Shift键并按住鼠标拖曳到目标模块的输入端正确答案:C2、在一个模型窗口上按住一个模块并同时按Shift键移动到另一个模型窗口,则()。

A.在两个模型窗口都有这个模块B.在后一个窗口有这个模块C.在前一个窗口有这个模块D.在两个窗口都有模块并添加连线正确答案:A3、为子系统定制参数设置对话框和图标,使子系统本身有一个独立的操作界面,这种操作称为子系统的()。

A.包装B.封装C.集成D.组合正确答案:B4、使用S函数时,要在模型编辑窗口添加()。

A.Sine Wave模块B.S-Program模块C.Subsystem模块D.S-Function模块正确答案:D二、多选题1、启动Simulink的方法有()。

A.在命令行窗口中输入simulink命令B.在“主页”选项卡中单击SIMULINK命令组中的“Simulink”命令按钮C.在“主页”选项卡中单击“文件”命令组中的“新建”命令按钮D.在“主页”选项卡中单击“文件”命令组中的“新建脚本”命令按钮正确答案:A、B、C2、根据控制信号的控制方式不同,条件执行子系统分为()。

A.事件驱动子系统B.使能子系统C.触发子系统D.使能加触发子系统正确答案:B、C、D3、以下关于S函数的描述中,正确的有()。

A.利用S函数可以对Simulink模块库进行扩充B.S函数只能用MATLAB语言编写C.S函数有现成的模板程序D.S函数模块能够被封装正确答案:A、C、D三、判断题1、建立系统仿真模型是在Simulink模型编辑窗口中进行的。

正确答案:√2、利用触发子系统能够将锯齿波转换为方波。

MATLAB语言与控制系统仿真_参考题答案_第3章

MATLAB语言与控制系统仿真_参考题答案_第3章

3.5 MATLAB 绘图实训3.5.1 实训目的1.学会MATLAB 绘图的基本知识;2.掌握MATLAB 子图绘制、图形注释、图形编辑等基本方法;3.学会通过MATLAB 绘图解决一些实际问题;4.练习二维、三维绘图的多种绘图方式,了解图形的修饰方法;5.学会制作简单的MATLAB 动画。

图3-46 炮弹发射示意图3.5.2 实训内容1. 炮弹发射问题〔1炮弹发射的基础知识炮弹以角度α射出的行程是时间的函数,可以分解为水平距离)(t x 和垂直距离)(t y 。

)cos()(0αtv t x = %水平方向的行程; 205.0)sin()(gt tv t y -=α %垂直方向的行程;其中,0v 是初速度;g 是重力加速度,为9.82m/s ;t 是时间。

〔2炮弹发射程序举例:分析以下程序以及图3-47各个图形的实际意义。

a=pi/4; v0=300; g=9.8;t=0:0.01:50; x=t*v0*cos<a>;y=t*v0*sin<a>-0.5*g*t.^2;subplot<221>;plot<t,x>;grid;title<‘时间-水平位移曲线'>; subplot<222>;plot<t,y>;grid;title<‘时间-垂直位移曲线'>; subplot<223>;plot<x,y>;grid;title<‘水平位移-垂直位移曲线'>; subplot<224>;plot<y,x>;grid;title<‘垂直位移-水平位移曲线'>; 图3-4745角发射曲线 〔3编程解决炮弹发射问题①假设在水平地面上以垂直于水平面的角度向上发射炮弹,即发射角90=α,假设初速度分别为[310,290,270]m/s,试绘制时间-垂直位移曲线,编程求取最高射程;绘图要求:◆ 标题设为"炮弹垂直发射问题";◆ 在图上通过添加文本的方式表明初速度; ◆ 在x 轴标注"时间";◆ 在y 轴上标注"垂直距离"; ◆ 添加网格线;◆ 将310m/s 的曲线改为线粗为2的红色实线; ◆ 将290m/s 的曲线改为线粗为3的绿色点划线;◆ 将270m/s 的曲线改为线粗为2的蓝色长点划线;a=pi/2; v1=310; g=9.8;t=0:0.01:50; x1=t*v1*cos<a>;y1=t*v1*sin<a>-0.5*g*t.^2;plot<t,y1>;grid; title<'炮弹垂直发射问题'>; xlabel<'时间'>; ylabel<'垂直距离'>; hold on; v2=290;x2=t*v2*cos<a>;y2=t*v2*sin<a>-0.5*g*t.^2; plot<t,y2>; v3=270;x3=t*v3*cos<a>;y3=t*v3*sin<a>-0.5*g*t.^2; plot<t,y3>;zgsc=[max<y1>; max<y2>; max<y3>] %三次发射的最高射程 运行结果如下: zgsc =1.0e+003 * 4.9031 4.29083.7194最高射程分别为:4903.1米,4290.8米,3719.4米。

MATLAB与控制系统仿真练习题(含图)及答案

MATLAB与控制系统仿真练习题(含图)及答案

MATLAB 与控制系统仿真练习题(含图)1、已知函数x x e x f x sin cos )(-=,作出函数的大致图像。

>> syms x>> y=exp(x)*cos(x)-sin(x); >> ezplot(y)2、求下列极限:(1)30sin lim xx x x -→ >> syms x>> y=(x-sin(x))/(x^3);>> limit(y,x,0)ans =1/6(2) xx x ⎪⎭⎫ ⎝⎛+∞→11lim >> y=(1+1/x)^x;>> limit(y,x,inf)ans =exp(1)3、求下列函数的导数:(1)x e y x sin =>> syms x>> y=exp(x)*sin(x);>> diff(y,x)ans =exp(x)*sin(x)+exp(x)*cos(x)(2) x e x x y 22sin +=>> syms x>> y=sin(x)+x^2*exp(2*x);>> diff(y,x)ans =cos(x)+2*x*exp(2*x)+2*x^2*exp(2*x)4、求.)1(532⎰-dx x x 和.sin ⎰xdx e x(1).)1(532⎰-dx x xsyms x>> int(x^2*(1-x^3)^5)ans =-1/18*x^18+1/3*x^15-5/6*x^12+10/9*x^9-5/6*x^6+1/3*x^3(2).sin ⎰xdx e x>> int(exp(x)*sin(x))ans =-1/2*exp(x)*cos(x)+1/2*exp(x)*sin(x)5、求.)(102⎰-dx x x 和.1102⎰-dx x x (1) .)(102⎰-dx x x>> syms x>> int(x-x^2,0,1)ans =1/6(2) .1102⎰-dx x x>> syms x>> int(x*(1-x^2)^0.5,0,1)ans =1/36、已知二元函数),(cos )sin(2xy xy z +=试求y x z x z y z x z ∂∂∂∂∂∂∂∂∂222,,,。

控制系统计算机辅助设计:MATLAB语言与应用(第2版)薛定宇-课后习题答案

第1章控制系统计算机辅助设计概述第2章 MATLAB语言程序设计基础第3章线性控制系统的数学模型第4章线性控制系统的计算机辅助分析第5章 Simulink在系统仿真中的应用第6章控制系统计算机辅助设计第1章控制系统计算机辅助设计概述【1】已阅,略【2】已阅,略【3】已经掌握help命令和Help菜单的使用方法【4】区别:MATLAB语言实现矩阵的运算非常简单迅速,且效率很高,而用其他通用语言则不然,很多通用语言所实现的矩阵运算都是对矩阵维数具有一点限制的,即使限制稍小的,但凡维数过大,就会造成运算上的溢出出错或者运算出错,甚至无法处理数据的负面结果【5】【8】(1)输入激励为正弦信号(2)输入激励为脉冲模拟信号(3)输入激励为时钟信号(4) 输入激励为随机信号(5) 输入激励为阶跃信号δ=0.3δ=0.05δ=0.7结论:随着非线性环节的死区增大,阶跃响应曲线的范围逐渐被压缩,可以想象当死区δ足够大时,将不再会有任何响应产生。

所以可以得到结论,在该非线性系统中,死区的大小可以改变阶跃响应的幅值和超调量。

死区越大,幅值、超调量将越小,而调整时间几乎不受其影响第2章 MATLAB语言程序设计基础【1】>> A=[1 2 3 4;4 3 2 1;2 3 4 1;3 2 4 1]A =1 2 3 44 3 2 12 3 4 13 24 1>> B=[1+4i,2+3i,3+2i,4+i;4+i,3+2i,2+3i,1+4i;2+3i,3+2i,4+i,1+4i;3+2i,2+3i,4+i,1+4i]B =1.0000 + 4.0000i2.0000 +3.0000i 3.0000 + 2.0000i4.0000 + 1.0000i4.0000 + 1.0000i 3.0000 + 2.0000i 2.0000 + 3.0000i 1.0000 + 4.0000i2.0000 +3.0000i 3.0000 + 2.0000i4.0000 + 1.0000i 1.0000 + 4.0000i3.0000 + 2.0000i 2.0000 + 3.0000i4.0000 + 1.0000i 1.0000 + 4.0000i>> A(5,6)=5A =1 2 3 4 0 04 3 2 1 0 02 3 4 1 0 03 24 1 0 00 0 0 0 0 5∴若给出命令A(5,6)=5则矩阵A的第5行6列将会赋值为5,且其余空出部分均补上0作为新的矩阵A,此时其阶数为5×6【2】相应的MATLAB命令:B=A(2:2:end,:)>> A=magic(8)A =64 2 3 61 60 6 7 579 55 54 12 13 51 50 1617 47 46 20 21 43 42 2440 26 27 37 36 30 31 3332 34 35 29 28 38 39 2541 23 22 44 45 19 18 4849 15 14 52 53 11 10 568 58 59 5 4 62 63 1>> B=A(2:2:end,:)B =9 55 54 12 13 51 50 1640 26 27 37 36 30 31 3341 23 22 44 45 19 18 488 58 59 5 4 62 63 1∴从上面的运行结果可以看出,该命令的结果是正确的【3】>> syms x s; f=x^5+3*x^4+4*x^3+2*x^2+3*x+6f =x^5 + 3*x^4 + 4*x^3 + 2*x^2 + 3*x + 6>> [f1,m]=simple(subs(f,x,(s-1)/(s+1)))f1 =19 - (72*s^4 + 120*s^3 + 136*s^2 + 72*s + 16)/(s + 1)^5m =simplify(100)【4】>> i=0:63; s=sum(2.^sym(i))s =0615【5】>> for i=1:120if(i==1|i==2) a(i)=1;else a(i)=a(i-1)+a(i-2);endif(i==120) a=sym(a); disp(a); endend[ 1, 1, 2, 3, 5, 8, 13, 21, 34, 55, 89, 144, 233, 377, 610, 987, 1597, 2584, 4181, 6765, 10946, 17711, 28657, 46368, 75025, 121393, 196418, 317811, 514229, 832040, 1346269, 2178309, 3524578, 5702887, 9227465, , , , , 5, 1, 6, 7, 3, 70, 03, 73, 76, 49, , 074, 099, 173, 272, 2445, 3717, 6162, 9879, 6041, 55920, 81961, 37881, 19842, 106, 177565, 035288, 212853, 248141, 0460994, , 1170129, 1879264, 8065, , , , 00884757, , 0, 5, 6, 1, 0, 88, , 673, 58, 931, , 120, , 029, 4, 2, 9905, 3072, 2977, 46049, 69026, 15075, 40, 99176, 083277, 082453, 165730, 248183, 7576, 62096, , 4738105, 5814114, 9, 186333, , 284885, 9, 3488322, 9, 0, 0]【6】>>k=1;for i=2:1000for j=2:iif rem(i,j)==0if j<i, break;endif j==i, A(k)=i; k=k+1; break; endendendenddisp(A);Columns 1 through 132 3 5 7 11 13 17 19 23 29 31 37 41 Columns 14 through 2643 47 53 59 61 67 71 73 79 83 89 97 101 Columns 27 through 39103 107 109 113 127 131 137 139 149 151 157 163 167 Columns 40 through 52173 179 181 191 193 197 199 211 223 227 229 233 239 Columns 53 through 65241 251 257 263 269 271 277 281 283 293 307 311 313 Columns 66 through 78317 331 337 347 349 353 359 367 373 379 383 389 397 Columns 79 through 91401 409 419 421 431 433 439 443 449 457 461 463 467 Columns 92 through 104479 487 491 499 503 509 521 523 541 547 557 563 569 Columns 105 through 117571 577 587 593 599 601 607 613 617 619 631 641 643 Columns 118 through 130647 653 659 661 673 677 683 691 701 709 719 727 733 Columns 131 through 143739 743 751 757 761 769 773 787 797 809 811 821 823 Columns 144 through 156827 829 839 853 857 859 863 877 881 883 887 907 911 Columns 157 through 168919 929 937 941 947 953 967 971 977 983 991 997【7】说明:h和D在MATLAB中均应赋值,否则将无法实现相应的分段函数功能syms x; h=input(‘h=’); D=input(‘D=’);y=h.*(x>D)+(h.*x/D).*(abs(x)<=D)-h.*(x<-D)【10】function y=fib(k)if nargin~=1,error('出错:输入变量个数过多,输入变量个数只允许为1!');endif nargout>1,error('出错:输出变量个数过多!');endif k<=0,error('出错:输入序列应为正整数!');endif k==1|k==2,y=1;else y=fib(k-1)+fib(k-2);endend【13】【14】>> t=[-1:0.001:-0.2,-0.1999:0.0001:0.1999,0.2:0.001:1];y=sin(1./t);plot(t,y);grid on;【15】(1) >> t=-2*pi:0.01:2*pi;r=1.0013*t.^2;polar(t,r);axis('square')(2) >> t=-2*pi:0.001:2*pi;r=cos(7*t/2);polar(t,r);axis('square')(3) >> t=-2*pi:0.001:2*pi;r=sin(t)./t;polar(t,r);axis('square')(4) >> t=-2*pi:0.001:2*pi;r=1-cos(7*t).^3;polar(t,r);axis('square')【17】(1)z=xy>> [x,y]=meshgrid(-3:0.01:3,-3:0.01:3);z=x.*y;mesh(x,y,z);>> contour3(x,y,z,50);(1)z=sin(xy)>> [x,y]=meshgrid(-3:0.01:3,-3:0.01:3);z=sin(x.*y);mesh(x,y,z);>> contour3(x,y,z,50);第3章线性控制系统的数学模型【1】(1) >> s=tf('s');G=(s^2+5*s+6)/(((s+1)^2+1)*(s+2)*(s+4))Transfer function:s^2 + 5 s + 6--------------------------------s^4 + 8 s^3 + 22 s^2 + 28 s + 16(2) >> z=tf('z',0.1);H=5*(z-0.2)^2/(z*(z-0.4)*(z-1)*(z-0.9)+0.6)Transfer function:5 z^2 - 2 z + 0.2---------------------------------------z^4 - 2.3 z^3 + 1.66 z^2 - 0.36 z + 0.6Sampling time (seconds): 0.1【2】(1)该方程的数学模型>> num=[6 4 2 2];den=[1 10 32 32];G=tf(num,den)Transfer function:6 s^3 + 4 s^2 + 2 s + 2------------------------s^3 + 10 s^2 + 32 s + 32(2)该模型的零极点模型>> G=zpk(G)Zero/pole/gain:6 (s+0.7839) (s^2 - 0.1172s + 0.4252)-------------------------------------(s+4)^2 (s+2)(3)由微分方程模型可以直接写出系统的传递函数模型【5】(1) >> P=[0;0;-5;-6;-i;i];Z=[-1+i;-1-i];G=zpk(Z,P,8)Zero/pole/gain:8 (s^2 + 2s + 2)-------------------------s^2 (s+5) (s+6) (s^2 + 1)(2) P=[0;0;0;0;0;8.2];Z=[-3.2;-2.6];H=zpk(Z,P,1,'Ts',0.05,'Variable','q')Zero/pole/gain:(q+3.2) (q+2.6)---------------q^5 (q-8.2)Sampling time (seconds): 0.05【8】(1)闭环系统的传递函数模型>> s=tf('s');G=10/(s+1)^3;Gpid=0.48*(1+1/(1.814*s)+0.4353*s/(1+0.4353*s));G1=feedback(Gpid*G,1)Transfer function:7.58 s^2 + 10.8 s + 4.8--------------------------------------------------------------0.7896 s^5 + 4.183 s^4 + 7.811 s^3 + 13.81 s^2 + 12.61 s + 4.8 (2)状态方程的标准型实现>> G1=ss(G1)a =x1 x2 x3 x4 x5 x1 -5.297 -2.473 -2.186 -0.9981 -0.7598x2 4 0 0 0 0 x3 0 2 0 0 0 x4 0 0 2 0 0x5 0 0 0 0.5 0b =u1x1 2x2 0x3 0x4 0x5 0c =x1 x2 x3 x4 x5y1 0 0 0.6 0.4273 0.3799d =u1y1 0Continuous-time state-space model.(3)零极点模型>> G1=zpk(G1)Zero/pole/gain:9.6 (s^2 + 1.424s + 0.6332)--------------------------------------------------------(s+3.591) (s^2 + 1.398s + 0.6254) (s^2 + 0.309s + 2.707)【11】>> Ga=feedback(s/(s^2+2)*1/(s+1),(4*s+2)/(s+1)^2);Gb=feedback(1/s^2,50);G=3*feedback(Gb*Ga,(s^2+2)/(s^3+14))Transfer function:3 s^6 + 6 s^5 + 3 s^4 + 42 s^3 + 84 s^2 + 42 s---------------------------------------------------------------------------s^10 + 3 s^9 + 55 s^8 + 175 s^7 + 300 s^6 + 1323 s^5 + 2656 s^4 + 3715 s^3+ 7732 s^2 + 5602 s + 1400【13】c1=feedback(G5*G4,H3)=G5*G4/(1+G5*G4*H3)c2=feedback(G3,H4*G4)=G3/(1+G3*H4*G4)c3=feedback(c2*G2,H2)=c2*G2/(1+c2*G2*H2)=G3*G2/(1+G3*H4*G4+G3*G2*H1)G=feedback(G6*c1*c3*G1,H1)=G6*c1*c3*G1/(1+ G6*c1*c3*G1*H1)=G6*G5*G4*G3*G2*G1/(1+G3*H4*G4+G3*G2*H1+G5*G4*H3+G5*G4*H3*G3*H4*G4+G5*G4* H3*G3*G2*H1+G6*G5*G4*G3*G2*G1*H1)【14】>> s=tf('s');c1=feedback(0.21/(1+0.15*s),0.212*130/s);c2=feedback(c1*70/(1+0.0067*s)*(1+0.15*s)/(0.051*s),0.1/(1+0.01*s));G=(1/(1+0.01*s))*feedback(130/s*c2*1/(1+0.01*s)*(1+0.17*s)/(0.085*s),0.0044/(1+0.01*s)) Transfer function:0.004873 s^5 + 1.036 s^4 + 61.15 s^3 + 649.7 s^2 + 1911 s---------------------------------------------------------------------------4.357e-014 s^10 + 2.422e-011 s^9 +5.376e-009 s^8 +6.188e-007 s^7+ 4.008e-005 s^6 + 0.001496 s^5 + 0.03256 s^4 + 0.4191 s^3+ 2.859 s^2 + 8.408 s 第4章线性控制系统的计算机辅助分析【1】(1) >> num=[1];den=[3 2 1 2];G=tf(num,den);eig(G)ans =-1.00000.1667 + 0.7993i0.1667 - 0.7993i分析:由以上信息可知,系统的极点有2个是在s域的右半平面的,因此系统是不稳定的(2) >> num=[1];den=[6 3 2 1 1];G=tf(num,den);eig(G)ans =-0.4949 + 0.4356i-0.4949 - 0.4356i0.2449 + 0.5688i0.2449 - 0.5688i分析:由以上信息可知,系统的极点有2个是在s域的右半平面的,因此系统是不稳定的(3) >> num=[1];den=[1 1 -3 -1 2];G=tf(num,den);eig(G)ans =-2.0000-1.00001.00001.0000分析:由以上信息可知,系统的极点有2个是在s域的右半平面的,因此系统是不稳定的(4) >> num=[3 1];den=[300 600 50 3 1];G=tf(num,den);eig(G)ans =-1.9152-0.14140.0283 + 0.1073i0.0283 - 0.1073i分析:由以上信息可知,系统的极点有2个是在s域的右半平面的,因此系统是不稳定的(5) >> s=tf('s');G=0.2*(s+2)/(s*(s+0.5)*(s+0.8)*(s+3)+0.2*(s+2));eig(G)ans =-3.0121-1.0000-0.1440 + 0.3348i-0.1440 - 0.3348i分析:由以上信息可知,系统的所有极点都在s域的左半平面,因此系统是稳定的【2】(1) >> num=[-3 2];den=[1 -0.2 -0.25 0.05];H=tf(num,den,'Ts',0.5);abs(eig(H)')ans =0.5000 0.5000 0.2000分析:由以上信息可知,所有特征根的模均小于1,因此该系统是稳定的(2) >> num=[3 -0.39 -0.09];den=[1 -1.7 1.04 0.268 0.024];H=tf(num,den,'Ts',0.5);abs(eig(H)')ans =1.1939 1.1939 0.1298 0.1298分析:由以上信息可知,由于前两个特征根的模均大于1,因此该系统是不稳定的(3) >> num=[1 3 -0.13];den=[1 1.352 0.4481 0.0153 -0.01109 -0.001043];H=tf(num,den,'Ts',0.5);abs(eig(H)')ans =0.8743 0.1520 0.2723 0.2344 0.1230分析:由以上信息可知,所有特征根的模均小于1,因此该系统是稳定的(4) >> num=[2.12 11.76 15.91];den=[1 -7.368 -20.15 102.4 80.39 -340];H=tf(num,den,'Ts',0.5,'Variable','q');abs((eig(H))')ans =8.2349 3.2115 2.3415 2.3432 2.3432分析:由以上信息可知,所有特征根的模均大于1,因此该系统是不稳定的【3】(1) >> A=[-0.2,0.5,0,0,0;0,-0.5,1.6,0,0;0,0,-14.3,85.8,0;0,0,0,-33.3,100;0,0,0,0,-10];eig(A)ans =-0.2000-0.5000-14.3000-33.3000-10.0000分析:由以上信息可知,该连续线性系统的A矩阵的所有特征根的实部均为负数,因此该系统是稳定的(2)>>F=[17,24.54,1,8,15;23.54,5,7,14,16;4,6,13.75,20,22.5589;10.8689,1.2900,19.099,…-4-3.5-3-2.5-2-1.5-1-0.50x 10-6P ole-Zero Map Real Axis (seconds -1)I m a g i n a r y A x i s (s e c o n d s -1)21.896,3;11,18.0898,25,2.356,9];abs(eig(F)') ans =63.7207 23.5393 12.4366 13.3231 19.7275分析:由以上信息可知,该离散系统的F 矩阵的所有特征根的模均大于1,因此该系统是不稳定的 【4】>> A=[-3 1 2 1;0 -4 -2 -1;1 2 -1 1;-1 -1 1 -2]; B=[1 0;0 2;0 3;1 1];C=[1 2 2 -1;2 1 -1 2];D=[0 0;0 0];G=ss(A,B,C,D); tzero(G)pzmap(G)ans =-3.6124-1.2765结论:∴可以得到该系统的 零点为-3.6124、-1.2765 分析:由以上信息可知,【5】>> s=tf('s');Gc=sscanform(G,'ctrl') Go=sscanform(G,'obsv') a =x1 x2 x3 x4 x1 0 1 0 0 x2 0 0 1 0 x3 0 0 0 1 x4 -0.4 -1.4 -4.3 -4.3 b =u1 x1 0 x2 0 x3 0 x4 1 c =x1 x2 x3 x4 y1 0.4 0.2 0 0 d =u1 y1 0Continuous-time state-space model. a =x1 x2 x3 x4x1 0 0 0 -0.4x2 1 0 0 -1.4x3 0 1 0 -4.3x4 0 0 1 -4.3b =u1x1 0.4x2 0.2x3 0x4 0c =x1 x2 x3 x4y1 0 0 0 1d =u1y1 0Continuous-time state-space model.【9】(1)>> num=[18 514 5982 36380 122664 222088 185760 40320];den=[1 36 546 4536 22449 67284 118124 109584 40320];[R1,P1,K1]=residue(num,[den 0]);[R1,P1]ans =-1.2032 -8.0000-1.0472 -7.00000.2000 -6.00000.7361 -5.0000-2.8889 -4.00002.2250 -3.0000-2.0222 -2.00003.0004 -1.00001.0000 0>> [n,d]=rat(R1);sym([n./d]')ans =[ -379/315, -377/360, 1/5, 53/72, -26/9, 89/40, -91/45, 7561/2520, 1][阶跃响应的解析解]y(t)=(-379/315)*e-8t+(-377/360)*e-7t+(1/5)*e-6t+(53/72)*e-5t+(-26/9)*e-4t+(89/40)*e-3t+ (-90/45)*e-2t+(7561/2520)*e-t+1(2) >> num=[18 514 5982 36380 122664 222088 185760 40320];den=[1 36 546 4536 22449 67284 118124 109584 40320];[R2,P2,K2]=residue(num,den);[R2,P2]ans =9.6254 -8.00007.3306 -7.0000-1.2000 -6.0000-3.6806 -5.000011.5556 -4.0000-6.6750 -3.00004.0444 -2.0000-3.0004 -1.0000>> [n,d]=rat(R2);sym([n./d]')ans =[ 3032/315, 887/121, -6/5, -265/72, 104/9, -267/40, 182/45, -7561/2520][脉冲响应的解析解]y(t)=(3032/315)*e-8t+(887/121)*e-7t+(-6/5)*e-6t+(-265/72)*e-5t+(104/9)*e-4t+(-267/40)*e-3t+(182/45)*e-2t+(-7561/2520)*e-t(3) >> syms t;u=sin(3*t+5);Us=laplace(u)Us =(3*cos(5) + s*sin(5))/(s^2 + 9)>> s=tf('s');Us=(3*cos(5)+s*sin(5))/(s^2+9);num=[18 514 5982 36380 122664 222088 185760 40320];den=[1 36 546 4536 22449 67284 118124 109584 40320];G=tf(num,den); Y=Us*G;num=Y.num{1}; den=Y.den{1};[R3,P3,K3]=residue(num,den); [R3,P3]ans =1.1237 -8.00000.9559 -7.0000-0.1761 -6.0000-0.6111 -5.00002.1663 -4.0000-1.1973 - 0.0010i 0.0000 + 3.0000i-1.1973 + 0.0010i 0.0000 - 3.0000i-1.3824 -3.00000.8614 -2.0000-0.5430 -1.0000>> [n,d]=rat(R3);sym([n./d]')ans =[109/97, 282/295, -59/335, -965/1579, 951/439, - 449/375 + (18*i)/17981, - 449/375 - (18*i)/17981, -1663/1203, 317/368, -82/151]Linear Simulation Results Time (seconds)A m p l i t u d e [正弦信号时域响应的解析解]y(t)=(109/97)*e -8t +(282/295)*e -7t +(-59/335)*e -6t +(-965/1579)*e -5t +(-449/375)*e -4t +(-1663/1203)*e -3t +(317/368)*e -2t +(-82/151)*e -t -2.3947sin(3t)[输出波形]>> num=[18 514 5982 36380 122664 222088 185760 40320];den=[1 36 546 4536 22449 67284 118124 109584 40320]; G=tf(num,den); t=[1:.1:20]';u=sin(3*t+5); lsim(G,u,t);分析:由解析解可知,输出信号的稳态部分是振荡的,并且其幅值与相位始终 在到达稳态的时候保持不变,因此 右图所示的输出波形与解析解所得的结论是一致的【10】(1)因为PI 或PID 控制器均含有Ki/s 节,则当Kp →∞,即|e(t)|一环节后,如果要求|e(t)|→0(2)不稳定系统能用PI 或PID 在积分控制中,控制器的输出与输入误差信号的积分成正比关系。

MATLAB语言及控制系统仿真_参考答案解析_第6章

6.4 控制系统频域分析MATLAB 仿真实训6.4.1实训目的 1. 学会利用MATLAB 绘制开环系统的伯德图; 2. 学会利用MATLAB 绘制开环系统的极坐标图;3. 掌握通过编程或相关命令求取系统稳定裕度的方法;4.通过仿真进一步理解掌握系统频域分析的有关知识。

6.4.2实训内容1.已知单位负反馈系统的开环传递函数为()()()10.5s 10.2s 1s 1.0ks G +++=)(要求编程绘制50=k 时的极坐标图,确定曲线与负实轴的交点坐标及频率值;n=50;d=conv([0.1,1],conv([0.2,1],[0.5,1])); sys=tf(n,d); nyquist(sys)曲线与负实轴的交点坐标为-3.76; 曲线与负实轴的交点频率值9.2;2.绘制下列系统的伯德图,并要求在图上显示出幅值裕度、相角裕度等信息。

(1)()()()18s 12s 2.6s G ++=>> n=2.6;>> d=conv([2,1],[8,1]); >> sys=tf(n,d)Transfer function: 2.6----------------------16 s^2 + 10 s + 1>> margin(sys)从图上信息可知,幅值裕度为无穷,相角裕度为88.2度。

(2)()()()110s 1s 10s G ++=s>> n=10;>> d=conv([1,0],conv([1,1],[10,1])); >> sys=tf(n,d)Transfer function: 10------------------------ 10 s^3 + 11 s^2 + s >> margin(sys)由图上信息可知,幅值裕度为-19.2dB,相角裕度为-34.3度。

(3)()()()10.5s 1s 1)8(10s s G +++=s从图上信息可知,幅值裕度为无穷,相角裕度为13.2度。

控制系统的Matlab仿真与设计课后答案

控制系统的Matlab仿真与设计课后答案第一篇:控制系统的Matlab仿真与设计课后答案MATLAB课后习题答案 2.1 x=[15 22 33 94 85 77 60] x(6)x([1 3 5])x(4:end)x(find(x>70))2.3 A=zeros(2,5);A(:)=-4:5L=abs(A)>3 islogical(L)X=A(L)2.4 A=[4,15,-45,10,6;56,0,17,-45,0] find(A>=10&A<=20)2.5 p1=conv([1,0,2],conv([1,4],[1,1]));p2=[1 0 1 1];[q,r]=deconv(p1,p2);cq='商多项式为';cr='余多项式为';disp([cq,poly2str(q,'s')]),disp([cr,poly2str(r,'s')])2.6 A=[11 12 13;14 15 16;17 18 19];PA=poly(A)PPA=poly2str(PA,'s')3.1 n=(-10:10)';y=abs(n);plot(n,y,'r.','MarkerSize',20)axis equal grid on xlabel('n')3.2 x=0:pi/100:2*pi;y=2*exp(-0.5*x).*sin(2*pi*x);plot(x,y),grid on;3.3 t=0:pi/50:2*pi;x=8*cos(t);y=4*sqrt(2)*sin(t);z=-4*sqrt(2)*sin(t);plot3(x,y,z,'p');title('Line in 3-D Space');text(0,0,0,'origin');xlabel('X'),ylable('Y'),zlable('Z');grid;3.4theta=0:0.01:2*pi;rho=sin(2*theta).*cos(2*theta);polar(theta,rho,'k');3.5[x,y,z]=sphere(20);z1=z;z1(:,1:4)=NaN;c1=ones(size(z1));surf(3*x,3*y,3*z1,c1);hold on z2=z;c2=2*ones(size(z2));c2(:,1:4)=3*ones(size(c2(:,1:4)));surf(1.5*x,1.5*y,1.5*z2,c2);col ormap([0,1,0;0.5,0,0;1,0,0]);grid on hold off 第四章function f=factor(n)if n<=1 f=1;elsef=factor(n-1)*n;endfunction[s,p]=fcircle(r)s=pi*r*r;p=2*pi*r;function k=jcsum1(n)k=0;i=0;while i<=n k=k+2^i;i=i+1;end function k=jcsum(n)k=0;for i=0:n k=k+2^i;end4.1for m=100:999m1=fix(m/100);m2=rem(fix(m/10),10);m3=rem(m,10);ifm==m1*m1*m1+m2*m2*m2+m3*m3*m3 disp(m)end end4.2[s,p]=fcircle(10)4.3y=0;n=100;for i=1:ny=y+1/i/i;end y4.4s=0;for i=1:5s=s+factor(i);end s4.5sum=0;i=1;while sum<2000 sum=sum+i;i=i+1;end;n=i-2 4.6jcsum(63)jcsum1(63)4.1 for m=100:999m1=fix(m/100);m2=rem(fix(m/10),10);m3=rem(m,10);if m==m1*m1*m1+m2*m2*m2+m3*m3*m3disp(m)end end 4.3 y=0;n=100;for i=1:ny=y+1/i/i;end y 4.4 s=0;for i=1:5s=s+factor(i);end s 4.5sum=0;i=1;while sum<2000sum=sum+i;i=i+1;end;n=i-2 4.6i=0;k=0;while i<=63k=k+2^i;i=i+1;end k ii=0;k=0;for i=0:63k=k+2^i;end i k第五章functionf=fxyz(u)x=u(1);y=u(2);z=u(3);f=x+y.^2./x/4+z.^2./y+2./z;5.1A=[2,1,-5,1;1,-5,0,7;0,2,1,-1;1,6,-1,-4];b=[13,-9,6,0]';x=Ab5.2[U,fmin]=fminsearch('fxyz',[0.5,0.5,0.5])5.3X=linspace(0,2*pi,50);Y=sin(X);P=polyfit(X,Y,3)AX=linspace(0,2*pi,50);Y=sin(X);Y1=polyval( P,X)plot(X,Y,':O',X,Y1,'-*')5.4x=0:2.5:10;h=[0:30:60]';T=[95,14,0,0,0;88,48,32,12,6;67,64,54,48,41];xi=[0:0.5:10];hi= [0:10:60]';temps=interp2(x,h,T,xi,hi,'cubic');mesh(xi,hi,temps);5.1 A=[2,1,-5,1;1,-5,0,7;0,2,1,-1;1,6,-1,-4];b=[13,-9,6,0]';x=Ab 5.3X=linspace(0,2*pi,50);Y=sin(X);P=polyfit(X,Y,3)AX=linspace(0,2*p i,50);Y=sin(X);Y1=polyval(P,X)plot(X,Y,':O',X,Y1,'-*')6.1syms x y=finverse(1/tan(x))6.2syms x yf=1/(1+x^2);g=sin(y);fg=compose(f,g)6.3syms xg=(exp(x)+x*sin(x))^(1/2);dg=diff(g)6.4F=int(int('x*exp(-x*y)','x'),'y')6.5syms xF=ztrans(x*exp(-x*10))6.6a=[0 1;-2-3];syms sinv(s*eye(2)-a);6.7 f=solve('a*x^2+b*x+c')6.8 f=solve('x+y+z=1','x-y+z=2','2*x-y-z=1')6.9y=dsolve('D2y+2*Dy+2*y=0','y(0)=1','Dy(0)=0')ezplot(y), grid on6.10a=maple('simplify(sin(x)^2+cos(x)^2);')6.11f=maple('laplace(exp(-3*t)*sin(t),t,s);')6.12 syms t xF=sin(x*t+2*t);L=laplace(F)第七章function[sys,x0,str,ts]=ww(t,x,u,flag)%¶¨ÒåÁ¬ÐøÏµÍ³µÄSº¯Êý A=[0,1;-0.4,-0.2];B=[0;0.2];C=[1,0];D=0;switch flag, case 0,[sys,x0,str,ts]=mdlInitializeSizes(A,B,C,D);case 1,sys=mdlDerivatives(t,x,u,A,B,C,D);case 2,sys=mdlUpdate(t,x,u);case 3,sys=mdlOutputs(t,x,u,A,B,C,D);case 4,sys=mdlGetTimeOfNextVarHit(t,x,u);case 9,sys=mdlTerminate(t,x,u);otherwiseerror(['Unhandled flag = ',num2str(flag)]);end%=============================== function [sys,x0,str,ts]=mdlInitializeSizes(A,B,C,D)sizes = simsizes;sizes.NumContStates = 2;sizes.NumDiscStates = 0;sizes.NumOutputs = 1;sizes.NumInputs = 1;sizes.DirFeedthrough = 1;sizes.NumSampleTimes = 1;sys = simsizes(sizes);x0 = [0;0];str = [];ts = [0 0];%=============================== function sys=mdlDerivatives(t,x,u,A,B,C,D)sys = A*x+B*u;%===============================function sys=mdlUpdate(t,x,u)sys = [];%=============================== functionsys=mdlOutputs(t,x,u,A,B,C,D)sys = C*x+D*u;%=============================== function sys=mdlGetTimeOfNextVarHit(t,x,u)sampleTime = 1;sys = t + sampleTime;%===============================function sys=mdlT erminate(t,x,u)sys = [];7.17.27.37.47.57.67.7第八章8.1num=[5];den=[1,2,2];sys=tf(num,den)8.1.2s = tf('s');H = [5/(s^2+2*s+2)];H.inputdelay =28.1.3h=tf([0.5,0],[1,-0.5,0.5],0.1)8.2num=2*[1,0.5];den=[1,0.2,1.01];sys=tf(num,den)[z,p,k]=tf2zp(num,den);zpk(z,p,k)[A,B,C,D]=tf2ss(num,den);ss(A,B,C,D)8.3 num=[1,5];den=[1,6,5,1];ts=0.1;sysc=tf(num,den);sysd=c2d(sysc,ts,'tustin')8.4.08.4.1 %¶Ôϵͳ·½¿òͼÿ¸ö»·½Ú½øÐбàºÅ,ÓÐ8¸öͨµÀ,ÁÐдÿ¸öͨµÀ´«µÝº¯Êý r1=1;r2=2;c1=3;c2=4;G1=r1;G2=tf(1,[c1,0]);G3=1;%ÊÇ·ÖÀëµãºÍ»ãºÏµãµÄÁ¬Ïß,²»Äܺϲ¢,´«º¯Îª1 G4=-1;G5=1/r2;G6=tf(1,[c2,0]);G7=-1;G8=-1;%½¨Á¢ÎÞÁ¬½ÓµÄ״̬¿Õ¼äÄ£ÐÍG=append(G1,G2,G3,G4,G5,G6,G7,G8)%д³öϵͳµÄÁ¬½Ó¾ØÕóQ=[1 4 0 %ͨµÀ1µÄÊäÈëÊÇͨµÀ4 2 1 7 %ͨµÀ2µÄÊäÈëÊÇͨµÀ1,7 3 2 02 0 53 8 6 5 0 7 5 060];%¸ººÅÔÚ´«º¯ÖÐÌåÏÖ %ÁгöϵͳµÄ×ܵÄÊäÈëºÍÊä³ö¶ËµÄ±àºÅinputs=1;outputs=6;%Éú³É×éºÏºóϵͳµÄ״̬¿Õ¼äÄ£ÐÍsys=connect(G,Q,inputs,outputs)8.4.2r1=1;r2=2;c1=3;c2=4;[A,B,C,D]=linmod('x84');[num,den ]=ss2tf(A,B,C,D);sys=tf(num,den)8.5A=[1,1,0;0,1,0;0,0,2];B=[0,0;1,0;0,-2];n=size(A)Tc=ctrb(A,B);if n==rank(T c)disp('ϵͳÍêÈ«ÄÜ¿Ø');elsedisp('ϵͳ²»ÍêÈ«ÄÜ¿Ø');end 第九章function [rtab,info]=routh(den)info=[];vec1=den(1:2:length(den));nrT=length(vec1);vec2=den(2:2:length(den)-1);rtab=[vec1;vec2,zeros(1,nrT-length(vec2))];for k=1:length(den)-2, alpha(k)=vec1(1)/vec2(1);for i=1:length(vec2),a3(i)=rtab(k,i+1)-alpha(k)*rtab(k+1,i+1);endif sum(abs(a3))==0 a3=polyder(vec2);info=[info,'All elements in row ',...int2str(k+2)' are zeros;'];elseif abs(a3(1))info=[info,'Replaced first element;'];endrtab=[rtab;a3, zeros(1,nrT-length(a3))];vec1=vec2;vec2=a3;end9.1num=[2,5,1];den=[1,2,3];bode(num,den);grid on;figure;nyquist(num,den);9.2num=5*[1,5,6];den=[1,6,10,8];step(num,den);gridon;figure;impulse(num,den);grid on;9.3kosi=0.7;wn=6;num=wn^2;den=[1,2*kosi*wn,wn^2];step(num,den);grid on;figure;impulse(num,den);gridon;9.4den=[1,2,8,12,20,16,16];[rtab,info]=routh(den)a=rtab(:,1) if all(a>0)disp('ϵͳÊÇÎȶ¨µÄ');elsedisp('ϵͳÊDz»Îȶ¨µÄ');end9.5num=7*[1,5];den=conv([1,0,0],conv([1,10],[1,1]));[gm,pm,wg,wc]=margin(num,den)9.1 >> sys=tf([2,5,1],[1,2,3])Transfer function: 2 s^2 + 5 s + 1---------------s^2 + 2 s + 3>> rlocus(sys)>> nyquist(sys)>> bode(sys)9.2>> G=tf(conv([5],[1,5,6]),[1,6,10,8]);>> step(G)>> impulse(G) sys=tf([5,25,30],[1,6,10,8]);>> step(sys)>> impulse(sys)9.4>> GH=tf(conv([7],[1,5]),conv([1,0,0],conv([1,10],[1,1])));>> [Gm,Pm,Wcg,Wcp]=margin(GH)Gm =0 Pm =-47.2870 Wcg =0 Wcp = 1.4354 >>GH=tf(conv([7],[1,5]),conv([1,0,0],conv([1,10],[1,1])));>>[Gm,Pm,Wcg,Wcp]=margin(GH)GH_close=feedback(GH,1)step( GH_close),grid on Gm =0 Pm =-47.2870 Wcg =0 Wcp =1.4354 第十章function s=bpts2s(bp,ts,delta)kosi=sqrt(1-1./(1+((1./pi).*log(1./bp)).^2));wn=-log(delta.*sqrt(1-kosi.^2))/(kosi.*ts);s=-kosi.*wn+j.*wn.*sqrt(1-kosi.^2);function[ngc,dgc]=fa_lead(ng0,dg0,Pm,wc,w)ngv=polyval(ng0,j*wc);dgv=polyval(dg0,j*wc);g=ngv/dgv;thetag=angle(g);mg=abs(g);thetar=Pm*pi/180;tz=(1+mg*cos(thetar-thetag))/(-wc*mg*sin(thetar-thetag));tp=(cos(thetar-thetag)+mg)/(wc*sin(thetar-thetag));ngc=[tz,1];dgc=[tp,1];function[ngc,dgc]=fg_lag_pm(ng0,dg0,w,Pm)[mu,pu]=bode(ng0,dg0,w);wgc=spline(pu,w,Pm+5-180);%²åÖµÇóÈ¡Âú×ãÏà½ÇÔ£¶ÈµÄ½ÇƵÂÊ×÷ΪÆÚÍûµÄ¼ôÇÐÆµÂÊngv=polyval(ng0,j*wgc);dgv=polyval(dg0,j*wgc);g=ngv/dgv;alph=abs(1/g);T=10/alph*wgc, ngc=[alph*T,1];dgc=[T,1];function[ngc,dgc]=fg_lag_wc(ng0,dg0,w,wc)ngv=polyval(ng0,j*wc);dgv=polyval(dg0,j*wc);g=ngv/dgv;alph=abs(1/g);T=10/(alph*wc);ngc=[alph*T,1];dgc=[T,1];function[ngc,dgc]=fg_lead_pd(ng0,dg0,wc)ngv=polyval(ng0,j*wc);dg v=polyval(dg0,j*wc);g=ngv/dgv;mg0=abs(g);t=sqrt(((1/mg0)^2-1)/(wc^2));%·ùÖµÏà¼ÓΪÁãngc=[t,1];dgc=[1];%Gc(s)=1+Tsfunction[ngc,dgc]=fg_lead_pm(ng0,dg0,Pm,w)[mu,pu]=bode(ng0,dg0,w);%¼ÆËãÔ-ϵͳµÄ¶ÔÊýƵÂÊÏìÓ¦Êý¾Ý[gm,pm,wcg,wcp]=margin(mu,pu,w);%ÇóÈ¡Ô-ϵͳµÄÏà½ÇÔ£¶ÈºÍ¼ôÇÐÆµÂÊalf=ceil(Pm-pm+5);%¼ÆËã¿ØÖÆÆ÷ÌṩµÄ×î´ó³¬Ç°½Ç¶È£¬phi=(alf)*pi/180;%½«×î´ó³¬Ç°½Çת»»Îª»¡¶Èµ¥Î»a=(1+sin(phi))/(1-sin(phi));%¼ÆËãaÖµdbmu=20*log10(mu);%ϵͳµÄ¶ÔÊý·ùÖµmm=-10*log10(a);%wm´¦µÄ¿ØÖÆÆ÷¶ÔÊý·ùÖµwgc=spline(dbm u,w,mm);%²îÖµÇóÈ¡wm£¬ÈÏΪwm£½wcT=1/(wgc*sqrt(a));%¼ÆËãTngc=[a*T,1];dgc=[T,1];function[ngc,dgc]=fg_lead_pm_wc(ng0,dg0,Pm,wc,w)[mu,pu]=bode(ng0,dg0,w);ngv=polyval(ng0,j*wc);dgv=poly val(dg0,j*wc);g=ngv/dgv;%ÇóÔ-ϵͳÔÚÆÚÍûµÄ¼ôÇÐÆµÂÊ´¦µÄƵÂÊÏìÓ¦Êý¾ÝG0(jwc) theta=180*angle(g)/pi;%Ô-ϵͳÔÚÆÚÍûµÄ¼ôÇÐÆµÂÊ´¦µÄÏà½ÇÔ£¶È£¬µ¥Î»Îª¶Èalf=ceil(Pm-(theta+180)+5);%×î´ó³¬Ç°½Ç phi=(alf)*pi/180;a=(1+sin(phi))/(1-sin(phi));dbmu=20*log10(mu);mm=-10*log10(a);wgc=spline(dbmu,w,mm);T=1/(wgc*sqrt(a));KK=128;s1=-2+i*2*sqrt(3);a=2ng0=[10];dg0=conv([1,0],conv([1,2],[1,8]));g0=tf(ng0,dg0);[ngc,dgc,k]=rg_lag(ng0,dg0,KK,s1,a);gc=tf(ngc,dgc)function s=kw2s(kosi,wn)s=-kosi.*wn+j*wn.*sqrt(1-kosi.^2);10.1ng0=[1];dg0=10000*[1 0-1.1772];g0=tf(ng0,dg0);%Âú×㿪»·ÔöÒæµÄΪУÕýϵͳµÄ´«µÝº¯Êýs=kw2s(0.7,0.5)%ÆÚÍûµÄ±Õ»·Ö÷µ¼¼«µãngc=rg_lead(ng0,dg0,s);gc=tf(ngc,1)g0c=tf(g0*gc);rlocus(g0 ,g0c);b1=feedback(g0,1);%δУÕýϵͳµÄ±Õ»·´«µÝº¯Êýb2=feedback(g0c,1);%УÕýºóϵͳµÄ±Õ»·´«µÝº¯Êýfigure,step(b1,'r--',b2,'b');gridon %»æÖÆÐ£ÕýǰºóϵͳµÄµ¥Î»½×Ô¾KK=20;s1=-2+i*sqrt(6);a=1ng0=[10];dg0=conv([1,0],[1,4]);g0=tf(ng0,dg0);[ngc,dgc,k]=rg_lag(ng0,dg0,KK,s1,a);gc=tf(ngc,dgc)g0c=tf(K K*g0*gc);b1=feedback(k*g0,1);b2=feedback(g0c,1);step(b1,'r--',b2,'b');grid ong0c=tf(KK*g0*gc);rlocus(g0,g0c);b1=feedback(k*g0,1);b2=feedback(g0c,1);figure,step(b1,'r--',b2,'b');grid onng0=[1];dg0=conv([1,0,0],[1,5]);g0=tf(ng0,dg0);w=logspace (-3,3);KK=1;Pm=50;[ngc,dgc]=lead4(ng0,dg0,KK,Pm,w);gc=tf(ngc,dgc);g0c=tf(KK*g0*gc);bode(KK*g0,w);holdon,bode(g0c,w);grid on,hold off [gm,pm,wcg,wcp]=margin(g0c)Kg=20*log10(gm)g1=feedback( g0c,1);bode(g1),grid on, [mag,phase,w]=bode(g1);a=find(mag<=0.707*mag(1));wb=w(a( 1))max(mag)b=find(mag==max(mag))wr=w(b)KK=40;Pm=50;ng0= KK *[1];dg0=conv([1,0],conv([1,1],[1,4]));g0=tf(ng0,dg0);w=logspace(-2,4);[ngc,dgc]=fg_lead_pm(ng0,dg0,Pm,w)gc=tf(ngc,dgc),g0c=tf(g0*gc);b1=feedback(g0,1);b2=feedback(g0c,1);step(b1,'r--', b2,'b');grid on figure, bode(g0,'r--',g0c,'b',w), grid on,[gm,pm,wcg,wcp]=margin(g0c), Km=20*log10(gm) KK=200;bp=0.3;ts=0.7;delta=0.05;ng0=[1];dg0=conv([1,0],conv([0.1,1],conv([0.021],conv([0.01,1],[0.005 1]))));g0=tf(ng0,dg0);w=logspace(-4,3);t=[0:0.1:3];[mag,phase]=bode(KK*g0,w);[gm0,pm0,wg0,wc0] =margin(mag,phase,w),gm0=20*log10(gm0)%gm0 =-15.6769 %2¡¢È·¶¨ÆÚÍûµÄ¿ª»·´«µÝº¯Êý mr=0.6+2.5*bp;wc=ceil((2+1.5*(mr-1)+2.5*(mr-1)^2)*pi/ts), h=(mr+1)/(mr-1)w1=2*wc/(h+1), w2=h*w1 w1=wc/10;w2=25;ng1=[1/w1,1];dg1=conv([1/w2,1],conv([1,0],[1 ,0]));g1=tf(ng1,dg1);g=polyval(ng1,j*wc)/polyval(dg1,j*wc);K=abs(1/g);%¼ôÇÐÆµÂÊ´¦·ùֵΪ1£¬ÇóKÖµ g1=tf(K*g1)%3¡¢È·¶¨·´À¡»·½Ú´«µÝº¯Êýh=tf(dg1,ng1);Kh=1/K;h=tf(Kh*h)%ÆÚÍûƵÂÊÌØÐԵĵ¹ÌØÐÔ%4¡¢ÑéËãÐÔÄÜÖ¸±êg2=feedback(KK*g0,h);%УÕýºó£¬ÏµÍ³µÄ¿ª»·´«µÝº¯Êýb1=feedback(KK*g0,1);b2=feedback(g2,1);bode(KK*g0,'r--',g2,'b',h,'g',w);grid onfigure,step(b1, 'r--',b2, 'b',t);grid on,[pos,tr,ts,tp]=stepchar(b2,delta)function[ngc,dgc]=lag2(ng0,dg0,w,KK,Pm)[mu,pu]=bode(KK*ng0,dg 0,w);wgc=spline(pu,w,Pm+5-180),ngv=polyval(KK*ng0,j*wgc);dgv=polyval(dg0,j*wgc);g=ngv/dgv;alph=abs(1/g), T=10/alph*wgc, ngc=[alph*T,1];dgc=[T,1];function[ngc,dgc]=lead4(ng0,dg0,KK,Pm,w)[mu,pu]=bode(KK*ng0,d g0,w);[gm,pm,wcg,wcp]=margin(mu,pu,w);alf=ceil(Pm-pm+5);phi=(alf)*pi/180;a=(1+sin(phi))/(1-sin(phi)), dbmu=20*log10(mu);mm=-10*log10(a);wgc=spline(dbmu,w,mm), T=1/(wgc*sqrt(a)),ngc=[a*T,1];dgc=[T,1];function[ngc,dgc]=ra_lead(ng0,dg0,s1)ngv=polyval(ng0,s1);dgv=pol yval(dg0,s1);g=ngv/dgv;thetag=angle(g);mg=abs(g);thetas=ang le(s1);ms=abs(s1);tz=(sin(thetas)-mg*sin(thetag-thetas))/(mg*ms*sin(thetag));tp=-(mg*sin(thetas)+sin(thetag+thetas))/(ms*sin(thetag));ngc=[tz,1]; dgc=[tp,1];function[ngc,dgc,k]=rg_lag(ng0,dg0,KK,s1,a)ngv=polyval(ng0,s1);dgv=polyval(dg0,s1);g=dgv/ngv;k=abs(g);%ÆÚÍûÖ÷µ¼¼«µã´¦µÄ¸ù¹ì¼£ÔöÒæ beta=k/KK;[kosi1,wn1]=s2kw(s1);zc=-wn1*sin(a*pi/180)/sin(pi-atan(sqrt(1-kosi1^2)/kosi1)-(a*pi/180));%ÀûÓÃÕýÏÒ¶¨Àí pc=beta*zc;ngc=beta*[1,-zc];dgc=[1,-pc];functionvarargout=rg_lead(ng0,dg0,s1)if nargout==1ngv=polyval(ng0,s1);dgv=polyval(dg0,s1);g=ngv/dgv;thetal=pi-angle(g);zc=real(s1)-imag(s1)/tan(thetal);t=-1/zc;varargout{1}=[t,1];elseif nargout==2 ngv=polyval(ng0,s1);dgv=polyval(dg0,s1);g=ngv/dgv;theta=angle(g);phi=angle(s1);if theta>0 phi_c=pi-theta;endif theta<0;phi_c=-thetaendtheta_z=(phi+phi_c)/2;theta_p=(phi-phi_c)/2;z_c=real(s1)-imag(s1)/tan(theta_z);p_c=real(s1)-imag(s1)/tan(theta_p);nk=[1-z_c];varargout{2}=[1-p_c];kc=abs(p_c/z_c);if theta<0 kc=-kcendvarargout{1}=kc*nk;elseerror('Êä³ö±äÁ¿ÊýÄ¿²»ÕýÈ·£¡');endfunction [bp,ts]=s2bpts(s,delta)[kosi,wn]=s2kw(s);bp=exp(-kosi.*pi./sqrt(1-kosi.^2));ts=-1./(kosi.*wn)*log(delta.*sqrt(1-kosi.^2));function[kosi,wn]=s2kw(s)kosi=1./sqrt(1+(imag(s)/real(s)).^2);wn=-real(s)./kosi;%Èç¹ûwnΪ¸ºÖµ£¬ÔòwnÈ¡Õý£¬²¢ÇÒkosiÈ¡·´iwn=(wn<0);wn(iwn)=-wn(iwn);kosi(iwn)=-kosi(iwn);function[pos,tr,ts,tp]=stepchar(g0,delta)[y,t]=step(g0);[mp,ind]=max(y);dimt=length(t);yss=y(dimt);pos=100*(mp-yss)/yss;tp=t(ind);for i=1:dimtif y(i)>=1 tr=t(i);break;end end;for i=1:length(y)ify(i)<=(1-delta)*yss|y(i)>=(1+delta)*yss ts=t(i);end end第十一章11.1a=[0 1 0;0 0 1;-1-5-6];b=[0 0 1]';p=[-2+4j;-2-4j;-10];K=acker(a,b,p)eig(a-b*K)11.2a=[0 1 0;0 0 1;-6-11-6];b=[1,0,0]';p=[-2+2*sqrt(3)*j;-2-2*sqrt(3)*j;-10];K=acker(a,b,p)eig(a-b*K)11.6A=[-1 0 0;0-2-3;0 0-3];B=[1 0;2 3;-3-3];C=[1 0 0;1 1 1 ];[G,K,L]=decoupling(A,B,C)11.8A=[0 20.6;1 0];b=[0 1]';c=[0 1];d=0;G=ss(A,b,c,d);Q=diag([1,0,0,0,0]);R=1;p=[-1.8+2.4j;-1.8-2.4j];[k,P]=lqr(A,b,Q,R);l=(acker(A',c',p))' Gc=-reg(G,k,l);zpk(Gc), eig(Gc.a), t=0:0.05:2;G_1=feedback(G*Gc,1);a1=eig(G_1.a), y_1=step(G_1,t);第十二章function[t,xx]=diffstate(G,H,x0,u0,N,T)xk=x0;u=u0;t=0 for k=1:N xk=G*xk+H*u;x(:,k)=xk;t=[t,k*T];end;xx=[x0,x];12.1function sys=M601(t,x)u=1;sys=[x(2);x(3);-800*x(1)-80*x(2)-24*x(3)+u];function[t,y]=ode4(A,B,C,D,x0,h,r,v,t0,tf)Ab=A-B*v*C;B=B;C=C;x=x0';y=0;t=t0;N=round((tf-t0)/h);for i=1:N k1=Ab*x+B*r;k2=Ab*(x+h*k1/2)+B*r;k3=Ab*(x+h*k2/2)+B*r;k4=Ab*(x+h *k3)+B*r;x=x+h*(k1+2*k2+2*k3+k4)/6;y=[y,C*x];t=[t,t(i)+h];end 12.1tspan=[0,10];x0=[0,0,0]';[t,y]=ode45('M601',tspan,x0);y1=800*y(:,1);plot(t,y1);12.2 num=10;den=conv([1,0],conv([1,2],[1,3]));[A,B,C,D]=tf2ss(num,den);x0=[0,0,0];v=1;t0=0;tf=10;h=0.01; r=1;[t,y]=ode4(A,B,C,D,x0,h,r,v,t0,tf);plot(t,y),grid12.3 12.4 g=[-2.8-1.4 0 0;1.4 0 0 0;-1.8-0.3-1.4-0.6;0 0 0.6 0];h=[1 0 1 0]';c=[0 0 0 1];d=0;x0=[0 0 0 0]';u=1;N=30;T=0.1;[t,xx]=diffstate(g,h,x0,u,N,T);plot(t,xx);y=c*xx;figurestairs(t,y)grid on12.6 第十四章14.1clear all;load optcar.mat;t=signals(1,:);p=signals(2,:);v=signals(3,:);a=signals(4,:);theta =signals(5,:);subplot(4,1,1);plot(t,p);gridon;ylabel('λÖÃ(m)');subplot(4,1,2);plot(t,v);gridon;ylabel('ËÙ¶È(m/s)');subplot(4,1,3);plot(t,a);gridon;ylabel('¼ÓËÙ¶È(m/s2)');subplot(4,1,4);plot(t,theta);gridon;ylabel('½Ç¶È(¶È)');14.1clear all loadcar.mat %½«µ¼Èëµ½car.matÖеķÂÕæÊµÑéÊý¾Ý¶Á³ö t=signals(1,:);x=signals(2,:);theta=signals(3,:);x1=signals(4,:);thet a1=signals(5,:);plot(t,x,t,x1);ylabel('С³µÎ»ÖÃ(m)'),grid on;%»æÖÆ¿ØÖÆÁ¦×÷ÓÃϽüËÆÄ£Ðͺ;«È·Ä£ÐÍxµÄµ¥Î»½×Ô¾ÏìÓ¦ÇúÏß figure % »æÖÆ¿ØÖÆÁ¦×÷ÓÃϽüËÆÄ£Ðͺ;«È·Ä£ÐÍthetaµÄµ¥Î»½×Ô¾ÏìÓ¦ÇúÏßplot(t,theta,t,theta1);ylabel('°Ú½ÇÖµ(rad)'),grid on;第二篇:控制系统的MATLAB仿真与设计课后答案第二章1>>x=[15 22 33 94 85 77 60] >>x(6)>>x([1 3 5])>>x(4:end)>>x(find(x>70))2>>T=[1-2 3-4 2-3];>>n=length(T);>>TT=T';>>for k=n-1:-1:0 >>B(:,n-k)=TT.^k;>>end >>B >>test=vander(T)3>>A=zeros(2,5);>>A(:) =-4:5 >>L=abs(A)>3 >>islogical(L)>>X=A(L)4>>A=[4,15,-45,10,6;56,0,17,-45,0] >>find(A>=10&A<=20)5>>p1=conv([1,0,2],conv([1,4],[1,1 ]));>>p2=[1 0 1 1];>>[q,r]=deconv(p1,p2);>>cq='商多项式为';cr='余多项式为';>>disp([cq,poly2str(q,'s')]),disp([cr,poly2str(r,'s')])6>>A=[11 12 13;14 15 16;17 18 19];>>PA=poly(A)>>PPA=poly2str(PA,'s')第三章1>>n=(-10:10)';>>y=abs(n);>>plot(n,y,'r.','MarkerSize',20)>>axisequal >>grid on>>xlabel('n')2>>x=0:pi/100:2*pi;>>y=2*exp(-0.5*x).*sin(2*pi*x);>>plot(x,y),gridon;3>>t=0:pi/50:2*pi;>>x=8*cos(t);>>y=4*sqrt(2)*sin(t);>>z=-4*sqrt(2)*sin(t);>>plot3(x,y,z,'p');>>title('Line in 3-D Space');>>text(0,0,0,'origin');>>xlabel('X'),ylable('Y'),zlable('Z');gr id;4>>theta=0:0.01:2*pi;>>rho=sin(2*theta).*cos(2*theta);>>po lar(theta,rho,'k');5>>[x,y,z]=sphere(20);>>z1=z;>>z1(:,1:4)=NaN ;>>c1=ones(size(z1));>>surf(3*x,3*y,3*z1,c1);>>holdon >>z2=z;>>c2=2*ones(size(z2));>>c2(:,1:4)=3*ones(size(c2(:,1 :4)));>>surf(1.5*x,1.5*y,1.5*z2,c2);>>colormap([0,1,0;0.5,0,0;1,0,0 ]);>>grid on >>hold off第四章1>>for m=100:999 m1=fix(m/100);m2=rem(fix(m/10),10);m3=rem(m,10);ifm==m1*m1*m1+m2*m2*m2+m3*m3*m3 disp(m)end end M文件:function[s,p]=fcircle(r)s=pi*r*r;p=2*pi*r;主程序:[s,p]=fcircle(10)3>>y=0;n=100;for i=1:n y=y+1/i/i;end >>y M文件:function f=factor(n)if n<=1 f=1;elsef=factor(n-1)*n;end主程序:>>s=0;for i=1:5 s=s+factor(i);end >>s5>>sum=0;i=1;while sum<2000 sum=sum+i;i=i+1;end;>>n=i-26 for循环M文件:function k=jcsum(n)k=0;for i=0:n k=k+2^i;end主程序: >>jcsum(63)While循环M文件: function k=jcsum1(n)k=0;i=0;while i<=n k=k+2^i;i=i+1;end主程序:>>jcsum1(63)第五章1>>A=[2,1,-5,1;1,-5,0,7;0,2,1,-1;1,6,-1,-4];>>b=[13,-9,6,0]';>>x=Ab M文件:functionf=fxyz(u)x=u(1);y=u(2);z=u(3);f=x+y.^2./x/4+z.^2./y+2./z;主程序:[U,fmin]=fminsearch('fxyz',[0.5,0.5,0.5])3>>X=linspace(0,2*p i,50);>>Y=sin(X);>>P=polyfit(X,Y,3)>>AX=linspace(0,2*pi,50);>>Y=sin(X);>>Y1=polyval(P,X)>>plot(X,Y,':O',X,Y1,'-*')4>>x=0:2.5:10;>>h=[0:30:60]';>>T=[95,14,0,0,0;88,48,32,12,6; 67,64,54,48,41];>>xi=[0:0.5:10];>>hi=[0:10:60]';>>temps=interp 2(x,h,T,xi,hi,'cubic');>>mesh(xi,hi,temps);第六章1>>syms x>>y=finverse(1/tan(x))2>>syms x y>>f=1/(1+x^2);g=sin(y);>>fg=compose(f,g)3>>syms x>>g=(exp(x)+x*sin(x))^(1/2);>>dg=diff(g)4>>F=int(int('x*e xp(-x*y)','x'),'y')5>>syms x>>F=ztrans(x*exp(-x*10))6>>a=[0 1;-2-3];>>syms s>>inv(s*eye(2)-a);7>>f=solve('a*x^2+b*x+c')8>>f=solve('x+y+z=1','x-y+z=2','2*x-y-z=1')9>>y=dsolve('D2y+2*Dy+2*y=0','y(0)=1','Dy(0)=0')>>ezpl ot(y),grid on10>>a=maple('simplify(sin(x)^2+cos(x)^2);')11>>f=maple(' laplace(exp(-3*t)*sin(t),t,s);')12>>syms t x>>F=sin(x*t+2*t);>>L=laplace(F)第七章第八章1-1>>h=tf([5,0],[1,2,2])1-2>>s = tf('s');>>H = [5/(s^2+2*s+2)];>>H.inputdelay =2 1-3>>h=tf([0.5,0],[1,-0.5,0.5],0.1)2>>num=2*[1,0.5];den=[1,0.2,1.01];>>sys=tf(num,d en)>>[z,p,k]=tf2zp(num,den);>>zpk(z,p,k)>>[A,B,C,D]=tf2ss(nu m,den);>>ss(A,B,C,D)3 >>num=[1,5];den=[1,6,5,1];ts=0.1;>>sys c=tf(num,den);>>sysd=c2d(sysc,ts,'tustin')>>r1=1;r2=2;c1=3;c2=4;>>[A,B,C,D]=linmod('x84');>>[num ,den]=ss2tf(A,B,C,D);>>sys=tf(num,den)5>>A=[1,1,0;0,1,0;0,0,2]; B=[0,0;1,0;0,-2];>>n=size(A)>>Tc=ctrb(A,B);if n==rank(Tc)disp('系统完全能控');elsedisp('系统不完全能控');end第九章1>>num=[2,5,1];den=[1,2,3];>>bode(num,den);grid on;>>figure;>>nyquist(num,den);2>>num=5*[1,5,6];den=[1,6,1 0,8];>>step(num,den);gridon;>>figure;>>impulse(num,den);gridon;3>>kosi=0.7;wn=6;>>num=wn^2;den=[1,2*kosi*wn,wn^2]; >>step(num,den);grid on;>>figure;>>impulse(num,den);grid on;4 M文件:function[rtab,info]=routh(den)info=[];vec1=den(1:2:length(den));nrT=len gth(vec1);vec2=den(2:2:length(den)-1);rtab=[vec1;vec2,zeros(1,nrT-length(vec2))];for k=1:length(den)-2, alpha(k)=vec1(1)/vec2(1);for i=1:length(vec2), a3(i)=rtab(k,i+1)-alpha(k)*rtab(k+1,i+1);endif sum(abs(a3))==0 a3=polyder(vec2);info=[info,'All elements in row ',...int2str(k+2)' are zeros;'];elseif abs(a3(1)) rtab=[rtab;a3, zeros(1,nrT-length(a3))];vec1=vec2;vec2=a3;end主程序:>>den=[1,2,8,12,20,16,16];>>[rtab,info]=routh(den)>>a=rt ab(:,1)if all(a>0)disp('系统是稳定的');elsedisp('系统是不稳定的');end5>>num=7*[1,5];den=conv([1,0,0],conv([1,10],[1,1]));>>[gm ,pm,wg,wc]=margin(num,den)第十章 M文件:。

matlab语言与控制系统仿真参考答案第1章

1.6 MATLAB操作基础-实训实训目的1.熟悉MATLAB语言环境,识别MATLAB桌面和命令窗口,命令历史窗口,工作空间窗口等;2.练习设置变量精度或变量显示方式;3.练习通过MATLAB编程解决一些实际问题;4.通过作图总结自变量增量设置对作图结果的影响;5.学会求解方程、方程组的基本方法;6.练习M文件的建立与执行;7.学会进入工具箱的演示系统,以便于进一步了解和学习感兴趣的知识,为以后的自主学习奠定基础。

1.根据表1-2中要求,先使用命令format改变变量的精度或显示方式,然后键入表达式,并将运行结果填入表1-2中,并练习使用clc 清除命令窗口中内容。

>> format bank;>> aa = 3.142.720.611.4164.0081.00>> format short>> aa =3.14162.71830.60651.414264.000081.0000>> format short e>> aa =3.1416e+0002.7183e+0006.0653e-001 1.4142e+000 6.4000e+001 8.1000e+001 >> format rat >> a a =355/113 1457/536 743/1225 1393/985 64 812. 编写MA TLAB 程序计算,根据程序运算结果填空(要求保留两位小数)(1)58.135.0+=-ea =( 4.29 ); 程序为: a=exp(-0.5)+sqrt(13.58)>> format bank;>> a=exp(-0.5)+sqrt(13.58) a =4.29(2)已知两个圆的半径分别为cm 5.31=r , cm 62=r ,则两个圆的周长分别为=1l ( 21.99 )cm, =2l ( 37.70 )cm;面积分别为=1s ( 38.48 )2cm ,=2s ( 113.10 )2cm程序:>> format bank; >> r=[3.5 6]; >> [2*pi*r;pi*r.*r]提示:圆的周长计算公式为:r L π2=,圆的面积2r S π=,其中r 为圆的半径。

  1. 1、下载文档前请自行甄别文档内容的完整性,平台不提供额外的编辑、内容补充、找答案等附加服务。
  2. 2、"仅部分预览"的文档,不可在线预览部分如存在完整性等问题,可反馈申请退款(可完整预览的文档不适用该条件!)。
  3. 3、如文档侵犯您的权益,请联系客服反馈,我们会尽快为您处理(人工客服工作时间:9:00-18:30)。
相关文档
最新文档