控制系统仿真matlab第六章习题答案
自动控制原理习题全解及MATLAB实验 第6章习题解答
系统开环传递函数为 G0 s
s0.1s
K
10.2s
1
,要求:
(1)系统响应斜坡信号 r(t)=t 时,稳态误差 ess 0.01 ;
(2) 系统相位裕量 ' 40 。
试用分析法设计一个串联滞后-超前校正装置。
解:(1)系统为Ⅰ型系统,在单位斜坡信号下
分稳态误差为
essr
1 k
令 essr
稳态性能与动态性能? 答:PID 兼有 PI、PD 控制的特点,它相当于提供了一个积分环节与两个一阶微分环节。
积分环节改善稳态性能,两个一阶微分环节改善动态性能。 试分别叙述利用比例负反馈和微分负反馈包围振荡环节所起到的作用。
答:二阶振荡环节的频率特性为
1
T 2S 2 2 S 1
用比例负反馈 H(s)=h
0.2s 1 0.0143s 1
(5) Gc (s)
s 1 14s 1
0.2s 1 0.0143s
(6)
G
k
(s)
s(14s
100(s 1) 1)(0.1s 1)(0.0143s
1)
' 180 [90 arctan 7 arctan(14 7) arctan(0.1 7) arctan(0.0143 7)] 41.9 40
10lg( 12)
6dB
,
最后得出 c' m 4.47rad/s>4.4rad/s
(4) 确定校正装置的转折频率
1 m
2.2rad/s ,2 m
8.8rad/s ,T 1 0.45s , 1
G(s)=
s
2.2 s
1 1
0.45s 0.11s
1 1
控制系统的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( 对于难以用文字或数字给出的情况,将提供MATLAB程序或Simulink模型。
第 1 章1.1 系统是被研究的对象,模型是对系统的描述,仿真是通过模型研究系统的一种工具或手段。
1.2 数学仿真的基本工具是数字计算机,因此也称为计算机仿真或数字仿真。
将数学模型通过一定的方式转变成能在计算机上实现和运行的数学模型,称之为仿真模型。
1.3 因为仿真是在模型上做试验,是一种广义的试验。
因此,仿真基本上是一种通过试验来研究系统的综合试验技术,具有一般试验的性质。
而进行试验研究通常是需要进行试验设计。
1.4 解析法又称为分析法,它是应用数学推导、演绎去求解数学模型的方法。
仿真法是通过在模型上进行一系列试验来研究问题的方法。
利用解析法求解模型可以得出对问题的一般性答案,而仿真法的每一次运行则只能给出在特定条件下的数值解。
,解析法常常是围绕着使问题易于求解,而不是使研究方法更适合于问题,常常因为存在诸多困难而不能适用。
从原则上讲,仿真法对系统数学模型的形式及复杂程度没有限制,是广泛适用的,但当模型的复杂程度增大时,试验次数就会迅速增加,从而影响使用效率。
1.5 仿真可以应用于系统分析、系统设计、理论验证和训练仿真器等方面。
1.6,8,20,71,,,,,,,,,x,100x,0u,,,, ,,,,0100,,,,y,,,002x注:本题答案是用MATLAB中tf2ss()函数给出的,是所谓“第二能控标准型”(下同)。
11.7,3,3,11,,,,,,,,,x,100x,0u,,,, ,,,,0100,,,,y,,,013x1.82s,3s,3G(s), 32s,4s,5s,21.91.368,0.36801,,,,,,,,x(k,1),100x(k),0u(k),,,, ,,,,0100,,,,y(k),,,00.3680.264x(k)1.10 仿真模型见praxis1_10_1.mdl;MATLAB程序见praxis1_10_2.m。
MATLAB基础及其应用教程-周开利-邓春晖课后答案 第六章
第6章习题解答1. 在命令窗口中输入,一次只能执行一行命令,脚本文件是命令的集中,只需输入脚本文件名,就可以按顺序一次执行所有命令。
2.脚本文件是命令的集中,只需输入脚本文件名,就可以一次执行所有命令,没有输入参数,没有输出参数。
文件名可任意取,只需符合变量的定义要求即可。
脚本文件中生成的变量是全局变量。
函数文件可以接受参数,也可以返回参数,一般需要由其他语句调用。
首行语句必须以function开始,有函数名,文件名必须与函数名一致。
函数空间中生成的变量是局部变量。
3.function result=project1(year)if mod(year,400)==0 || (mod(year,4)==0 && mod(year,100)~=0)result='是闰年';elseresult='不是闰年';end4.function min=comp(x1,x2)min=x1;if x1>x2min=x2;end5.(1)在命令窗口输入:k=0;for i=-1000:1000;k=k+1;end;k 即可得出:2001步。
(2)方法同上,10步6.(1)循环进行的条件是var不是10的倍数,循环结束时,var必是10的倍数,即10(2)每次循环的结果是:4,16,256,此时,不满足循环条件而退出循环,即var为256(3)循环条件是var>100,而var的值为3,不满足条件,循环不执行,所以var值不变。
7.方法1:clear;A=[1 2 3;4 5 6];fid=fopen('c:\a1.txt','wt');fprintf(fid,'%d ',A(1,:));fprintf(fid,'\n');fprintf(fid,'%d ',A(2,:));fclose(fid);方法2:clear;a=[1 2];A=[1 2 3 4 5 6]B=reshape(A,2,2) dlmwrite('c:\a2.txt',B)。
Matlab第六章作业答案
第六章
6.2 如何进行一下操作:
(1)翻转模块
(2)给模型窗口加标题
(3)指定仿真时间
(4)设置示波器的显示刻度
解:(1)点击模型窗口的Format ,在下拉菜单中中有两个选项:Flip block 和Rotate block 。
Flip block 可使模块旋转180度(快捷键Ctrl+I ),Rotate block 可以模块顺时针旋转90度(快捷键Ctrl+R )。
(2)在模型窗口的上方双击鼠标左键,会出现的一个文本输入框,编辑输入标题。
(3)单击Simulation ,在下拉菜单中选择configuration parameters ,打开仿真环境参数对话框,在Solver 选项的start time 设置仿真起止时间。
(4)双击示波器Scope 调出显示屏幕,选择Axes properties 选项,改变Y-min 和Y-max 可改变显示刻度。
6.3用Simulink 建立如下控制系统的仿真模型,并对系统进行阶跃响应仿真模型:
22()48
G s s s =++ 解:
6.4在Simulink 环境下,设计一个PID 控制器,实现下面被控制对象的控制,并观察选择不同的PID 参数时对控制效果的影响:
220()212G s s s =++ 其中,系统输入信号分别选择阶跃信号和正弦信号。
解:当输入信号为阶跃信号时:
当输入信号为正弦信号时:。
中南大学matlab课后答案-第六章
实验指导二1,(1)>> A=[1/2,1/3,1/4;1/3,1/4,1/5;1/4,1/5,1/6];B=[,,]';x=A\B (2)>> A=[1/2,1/3,1/4;1/3,1/4,1/5;1/4,1/5,1/6];B=[,,]';x=A\B (3)>> cond(A)2,(1)建立函数文件,命令如下;Function fx=funx(x)fx=x^41+x^3+1;$调用fzero函数求根,命令如下;>> z=fzero(@funx,-1)(2) 建立函数文件,命令如下;Function fx=sin(x)fx=x-(sin(x))/x;调用fzero函数求根,命令如下;>> z=fzero(@sin,(3) 建立函数文件,命令如下;function q=myfun(p)x=p(1);、y=p(2);z=p(3);q(1)=sin(x)+y^2+log(z)-7;q(2)=3*x+2^y-z^3+1;q(3)=x+y+z-5;调用fsolve函数求根,命令如下;>> options=optimset('Display','off');x=fsolve(@myfun,[1,1,1]',options) 3,(1) 建立函数文件,命令如下;function yp=funt(t,y)>yp=-y*+sin(10*t));求微分方程,程序如下:>> t0=0;tf=5;y0=1;[t,y]=ode23(@funt,[t0,tf],y0)(2) 建立函数文件,命令如下;function yp=funr(t,y)yp=cos(t)-(y)/(1+t^2);求微分方程,程序如下:>> t0=0;tf=5;y0=1;[t,y]=ode23(@funr,[t0,tf],y0)4,建立函数文件命令如下:{function fx=mymax(x)fx=-1*(1+x^2)/(1+x^4);求最大值,程序如下:>> [x,y]=fminbnd(@mymax,0,2)5,编写目标函数M文件,命令如下:function f=fop(x)f=-1*(x(1)^(1/2)+x(2)^(1/2)+x(3)^(1/2)+x(4)^(1/2));设定约束条件,并调用fmincon函数求解此约束最优化问题,程序如下:>>x0=[200,200,200,200];%A=[1,0,0,0;,1,0,0;,,1,0;,,,1];b=[400,440,484,];Ib=[0,0,0,0];options=optimset('Display','off');[x,y]=fmincon(@fop,x0,A,b,[],[],Ib,[],[],options)思考练习1,(1)矩阵求逆:>> A=[2,3,5;3,7,4;1,-7,1];b=[10,3,5]';x=inv(A)*b矩阵除法:>> A=[2,3,5;3,7,4;1,-7,1];b=[10,3,5]';x=A\b!矩阵分解:>> A=[2,3,5;3,7,4;1,-7,1];b=[10,3,5]';[L,U]=lu(A);x=U\(L\b) (2)矩阵求逆:>>A=[6,5,-2,5;9,-1,4,-1;3,4,2,-2;3,-9,0,2];b=[-4,13,1,11]';x=inv(A)*b矩阵除法:>>A=[6,5,-2,5;9,-1,4,-1;3,4,2,-2;3,-9,0,2];b=[-4,13,1,11]';x=A\b矩阵分解:>> A=[6,5,-2,5;9,-1,4,-1;3,4,2,-2;3,-9,0,2];b=[-4,13,1,11]'; [L,U]=lu(A);x=U\(L\b)2,(1)建立函数文件,命令如下:function fx=fun1(x)fx=3*x+sin(x)-exp(x);调用fzero函数求根,命令如下;>> y=fzero(@fun1,、(2) 建立函数文件,命令如下:function fx=fun2(x)fx=1-(1/x)+5;调用fzero函数求根,命令如下;>> y=fzero(@fun2,1)(3) 建立函数文件,命令如下;function q=fun3(p)x=p(1);y=p(2);q(1)=x^2+y^2-9;'q(2)=x+y-1;调用fsolve函数求根,命令如下;>> options=optimset('Display','off');x=fsolve(@fun3,[1,1]',options) 3,(1)建立函数文件,命令如下:function ydot=fun5(t,y)ydot(1)=(2-3*y(2)-2*t*y(1))/(1+t^2);ydot(2)=y(1);ydot=ydot';求解微分方程,命令如下:'>> t0=0;tf=5;x0=[0,1];[t,y]=ode45(@fun5,[t0,tf],x0);[t,y](2) (1)建立函数文件,命令如下:function ydot=fun6(t,y)ydot(1)=cos(t)+(5*y(1)*cos(2*t))/(t+1)^2-y(2)-y(3)/(3+sin(t)); ydot(2)=y(1);ydot(3)=y(2);ydot=ydot';求解微分方程,命令如下:>> t0=0;tf=5;x0=[0,1];[t,y]=ode45(@fun5,[t0,tf],x0);[t,y]4,】建立函数文件命令如下:function fx=max(x)fx=-1*(sin(x)+cos(x^2));求最大值,程序如下:>> [x,y]=fminbnd(@mymax,0,pi)5,编写目标函数M文件,命令如下:function f=topm(x)f=-1*x*(3-2*x)^2;设定约束条件,并调用fmincon函数求解此约束最优化问题,程序如下:>>x0=[0];A=[1];b=[];Ib=[0];options=optimset('Display','off');[x,y]=fminco n(@top,x0,A,b,[],[],Ib,[],[],options)。
matlab简明教程第六章答案
m a t l a b简明教程第六章答案本页仅作为文档封面,使用时可以删除This document is for reference only-rar21year.March第6章习题参考答案1. 假设x = [-3, 0, 0, 2, 5, 8]且y = [-5, -2,0, 3, 4, 10]。
通过手算得到以下运算的结果,并使用MATLAB检验计算的结果。
(1) z = y<~x(2)z = x&y(3) z = x|y(4) z = xor(x,y)参考答案: (1) z = [1, 1, 1, 0, 0, 0];(2) z = [1, 0, 0, 1, 1, 1];(3) z = [1, 1, 0, 1, 1, 1];(4) z = [0, 1, 0, 0, 0, 0]2. 在MATLAB中使用一个循环确定:如果用户最初在一个银行帐户中存储$10000,并且在每年的年终再存储$10000(银行每年支付6%的利息),那么账户上要积累$1000000要需要多长时间。
参考答案: 33年。
3.某个特定的公司生产和销售高尔夫手推车。
每周周末,公司都将那一周所生产的手推车转移到仓库(库存)之中。
卖出的所有手推车都是从库存中提取。
这个过程的一个简单模型为:I (k + 1) = P(k) + I (k) - S(k)其中:P(k) = 第k周所生产的手推车数量;I (k) = 第k周库存中的手推车数量;S(k) = 第k周所卖出的手推车数量;以下为10周计划中的每周销售额;假设每周的产量都基于前一周的销售额,所以有P(k) = S(k - 1)。
假设第一周的产量为50辆手推车:即,P(1) = 50。
编写一个MATLAB程序计算:10周之内每周库存之中的手推车数量或者计算手推车库存数量减少到0为止的时间,并同时绘制图形。
针对以下两种情况运行该程序:(1)初始库存为50辆手推车,所以I(1)= 50;(2)初始库存为30辆手推车,所以I (1) = 30。
matlab仿真课后习题
第一章习题3.请指出以下的变量名(函数名、M文件名)中,哪些是合法的Abc 2004x lil-1 wu_2004 a&b _xyz 解:合法的变量名有:Abc wu_20044.指令窗操作(1)求[12+2×(7-4)]÷32的运算结果解:>> [12+2*(7-4)]/3^2ans =2(2)输入矩阵A=[1,2,3;4,5,6;7,8,9],观察输出。
解:>> A=[1,2,3;4,5,6;7,8,9]A =1 2 34 5 67 8 9(3)输入以下指令,观察运算结果;clear;x=-8::8;y=x';X=ones(size(y))*x;Y=y*ones(size(x));R=sqrt(X.^2+Y.^2)+eps;Z=sin(R)./R;mesh(X,Y,Z);colormap(hot)xlabel('x'),ylabel('y'),zlabel('z')解:7.指令行编辑(1)依次键入以下字符并运行:y1=2*sin*pi)/(1+sqrt(5))解:>>y1=2*sin*pi)/(1+sqrt(5))y1 =(2)通过反复按键盘的箭头键,实现指令回调和编辑,进行新的计算;y2=2*cos*pi)/(1+sqrt(5))解:>>y2=2*cos*pi)/(1+sqrt(5))y2 =11.编写题4中(3)的M脚本文件,并运行之。
解:第二章习题1.在指令窗中键入x=1::2和y=2::1,观察所生成的数组。
解:>> x=1::2 x =>> y=2::1 y =Empty matrix: 1-by-02.要求在[0,2π]上产生50个等距采样数据的一维数组,试用两种不同的指令实现。
解: y1=0:2*pi/49:2*pi y2=linspace(0,2*pi,50)3.计算e -2t sint ,其中t 为[0,2π]上生成的10个等距采样的数组。
matlab课后习题答案1到6章
欢迎共阅习题二1.如何理解“矩阵是MATLAB 最基本的数据对象”? 答:因为向量可以看成是仅有一行或一列的矩阵,单个数据(标量)可以看成是仅含一个元素的矩阵,故向量和单个数据都可以作为矩阵的特例来处理。
(2) 删除矩阵A 的第7号元素。
答:A(7)=[](3) 将矩阵A 的每个元素值加30。
答:A=A+30;(4) 求矩阵A 的大小和维数。
答:size(A);ndims(A);(5) 将向量 t 的0元素用机器零来代替。
答:t(find(t==0))=eps; (6) 将含有12个元素的向量 x 转换成34⨯矩阵。
答:reshape(x,3,4); (7) 求一个字符串的ASCII码。
答:abs(‘123’); 或E 。
答:B=A(1:3,:); C=A(:,1:2); D=A(2:4,3:4);E=B*C;(2) 分别求E<D 、E&D 、E|D 、~E|~D 和find(A>=10&A<25)。
答:E<D=010001⎡⎤⎢⎥⎢⎥⎢⎥⎣⎦,E&D=110111⎡⎤⎢⎥⎢⎥⎢⎥⎣⎦,11⎡⎤⎢⎥答:student(1).id='0001';student(1).name='Tom'; student(1).major='computer';student(1).grade=[89,78,67,90,86,85]; 8.建立单元矩阵B 并回答有关问题。
B{1,1}=1;B{1,2}='Brenden';B{2,1}=reshape(1:9,3,3); B{2,2}={12,34,2;54,21,3;4 ,23,67};(1)size(B)和ndims(B)的值分别是多少?答:size(B) 的值为2, 2。
ndims(B) 的值为2。
(2)B(2)和B(4)的值分别是(2)建立5×6随机矩阵A,其元素为[100,200]范围内的随机整数。
控制系统的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程序设计教程(第二版)》第6章__MATLAB解方程与最优化问题求解
例6-7 分别用Jacobi迭代和Gauss-Serdel迭代法求解 下列线性方程组,看是否收敛。 命令如下: a=[1,2,-2;1,1,1;2,2,1]; b=[9;7;6]; [x,n]=jacobi(a,b,[0;0;0]) [x,n]=gauseidel(a,b,[0;0;0])
6.2 非线性方程数值求解 6.2.1 单变量非线性方程求解 在MATLAB中提供了一个fzero函数,可以用来求单变量非 线性方程的根。该函数的调用格式为: z=fzero('fname',x0,tol,trace) 其中fname是待求根的函数文件名,x0为搜索的起点。一个 函数可能有多个根,但fzero函数只给出离x0最近的那个根。 tol控制结果的相对精度,缺省时取tol=eps,trace指定迭代 信息是否在运算中显示,为1时显示,为0时不显示,缺省 时取trace=0。
例6-9 求下列非线性方程组在(0.5,0.5) 附近的数值解。 (1) 建立函数文件myfun.m。 function q=myfun(p) x=p(1); y=p(2); q(1)=x-0.6*sin(x)-0.3*cos(y); q(2)=y-0.6*cos(x)+0.3*sin(y); (2) 在给定的初值x0=0.5,y0=0.5下,调用fsolve函数求方程的 根。 x=fsolve('myfun',[0.5,0.5]',optimset('Display','off')) x= 0.6354 0.3734
例6-10 设有初值问题,试求其数值解,并与精确解 相比较。 (1) 建立函数文件funt.m。
function yp=funt(t,y) yp=(y^2-t-2)/4/(t+1); (2) 求解微分方程。 t0=0;tf=10; y0=2; [t,y]=ode23('funt',[t0,tf],y0); %求数值解 y1=sqrt(t+1)+1; t' y' y1' %求精确解
matlab答案第六章
第六章1. 利用MATLAB提供的randn函数生成符合正态分布的10×5随机矩阵A,进行如下操作:(1) A各列元素的均值和标准方差。
(2) A的最大元素和最小元素。
(3) 求A每行元素的和以及全部元素之和。
(4) 分别对A的每列元素按升序、每行元素按降序排列。
答:clear all; close all; clc;A=randn(10, 5);meanA=mean(A); %(1)A各列元素的均值stdA=std(A); %(1)A各列元素的标准方差maxA=max(max(A)); %(2)A的最大元素minA=min(min(A)); %(2)A的最小元素rowsumA=sum(A, 2); %(3)A每行元素的和sumA=sum(rowsumA); %(3)A全部元素之和sort1=sort(A); %(4)A的每列元素按升序排列sort2=sort(A, 2, 'descend'); %(4)A的每行元素按降序排列2. 按要求对指定函数进行插值和拟合。
(1) 按表6.1用3次样条方法插值计算0~90D范围内整数点的正弦值和0~75D范围内整数点的正切值,然后用5次多项式拟合方法计算相同的函数值,并将两种计算结果进行比较。
表6.1 特殊角的正弦和正切值表α度0 15 30 4560 75 90sinα0 0.2588 0.50000.7071 0.8660 0.9659 1.0000 tanα0 0.2679 0.57741.0000 1.7320 3.7320(2) 按表6.2用3次多项式方法插值计算1~100之间整数的平方根。
表6.2 1~100内特殊值的平方根表N 1 4 9 16 25 36 49 64 81 100N的平方根 1 2 3 4 5 6 7 8 9 10答:(1) 程序设计:clear all; close all; clc;alpha1=0:15:90;sin_alpha1=sin(alpha1*pi/180); %精确正弦值plot(alpha1, sin_alpha1, 'k:p'); hold on; %绘精确正弦曲线alpha2=0:90;sin_Y1=interp1(alpha1, sin_alpha1, alpha2, 'spline'); %3次样条正弦插值plot(alpha2, sin_Y1, 'r-*'); hold on; %绘3次样条插值正弦曲线P1=polyfit(alpha1, sin_alpha1, 5); %5次多项式拟合sin_Y2= polyval(P1, alpha2); %5次多项式求值plot(alpha2, sin_Y2, 'b-o'); %绘5次多项式插值正弦曲线legend('精确正弦值', '3次样条正弦插值', '5次多项式正弦插值'); title('正弦值比较'); alpha3=0:15:75;tan_alpha3=tan(alpha3*pi/180); %精确正切值figure, plot(alpha3, tan_alpha3, 'k:p'); hold on; %绘精确正切曲线alpha4=0:75;tan_Y1=interp1(alpha3, tan_alpha3, alpha4, 'spline'); %3次样条正切插值plot(alpha4, tan_Y1,'r-*'); hold on; %绘3次样条正切曲线P2=polyfit(alpha3, tan_alpha3, 5); %5次多项式拟合tan_Y2= polyval(P2, alpha4); %5次多项式求值plot(alpha4, tan_Y2, 'b-o'); %绘5次多项式插值正弦曲线legend('精确正切值', '3次样条正切插值', '5次多项式正切插值'); title('正切值比较');(2)程序设计:clear all; close all; clc;X=[1, 4, 9, 16, 25, 36, 49, 64, 81, 100]; Y=1:10;X1=1:100; Y1=interp1(X, Y, X1,'cubic');plot(X, Y, 'r:o'); hold on; %绘精确平方根曲线plot(X1, Y1, 'k-x'); %绘3次多项式插值平方根曲线legend('精确平方根', '3次多项式插值');3. 已知一组实验数据如表6.3所示。
控制系统的Matlab仿真与设计课后答案
end
end
4.2[s,p]=fcircle(10)
4.3y=0;n=100;
fori=1:n
y=y+1/i/i;
end
y
4.4s=0;
fori=1:5
s=s+factor(i);
end
s
4.5sum=0;i=1;
whilesum<2000
sum=sum+i;
i=i+1;
end;
n=i-2
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.4
theta=0:0.01:2*pi;
rho=sin(2*theta).*cos(2*theta);
polar(theta,rho,'k');
3.5
[x,y,z]=sphere(20);
fg=compose(f,g)
6.3symsx
g=(exp(x)+x*sin(x))^(1/2);
Y=sin(X);
P=polyfit(X,Y,3)
MATLAB教程2012a第6章习题解答-张志涌
MATLAB教程2012a第6章习题解答-张志涌第6章 M 文件和句柄函数习题6及解答1 请分别写出用for 和while 循环语句计算1000000210000002.02.02.012.0+++==∑= i i K 的程序。
此外,还请写出避免循环的数值、符号计算程序。
〖解答〗(1)for 环tics1=0;for k=0:1e6s1=s1+0.2^(k);end;s1tocs1 =1.2500Elapsed time is 1.453482 seconds.(2)while 环tics2=1;k=1;while k<1e6+1 %注意:上限与for 环不同s2=s2+0.2^k;k=k+1;ends2tocs2 =1.2500Elapsed time is 2.716870 seconds.(3)数值求和指令tics3=sum(0.2.^(0:1e6))tocs3 =1.2500Elapsed time is 0.626723 seconds.(4)符号求和指令ticsyms k;s4=vpa(symsum(0.2^k,0,1e6))tocs4 =1.2500000000000000000000000000000Elapsed time is 4.029501 seconds.〖说明〗●使用数值求和指令和“数组运算”的计算速度最快。
●符号计算的精度最高,但速度慢。
●以上程序运行时间仅供参考。
具体时间与所用机器、那程序是否初次运行、在MATLAB指令窗还是在M-book中运行等因素有关。
2 编写一个函数M文件,它的功能:没有输入量时,画出单位圆(见图p6.2-1);输入量是大于2的自然数N时,绘制正N边形,图名应反映显示多边形的真实边数(见图p6.2-2);输入量是“非自然数”时,给出“出错提示”。
此外,函数M文件应有H1行、帮助说明和程序编写人姓名。
〖解答〗(1)函数M文件function prob_solve602(n)% prob_solve602(n) plot a circle or a polygon with n edges % prob_solve602 plot a circle% n 应为大于2的自然数% By ZZY, 2006-2-15if nargin==0t=0:pi/100:2*pi;x=exp(i*t);str='Circle';elseif (nargin~=0)&(n<=2)error('输入量应是大于2的自然数')end;if n-round(n)~=0 %检查非自然数error('输入量应是大于2的自然数')end;t=(0:n)/n*2*pi;x=exp(i*t);str=['Polygon with ', int2str(n),' edges']; % 合成字符串endplot(real(x),imag(x),'r','LineWidth',4)title(str)axis square image offshg(2)各典型运行情况prob_solve602prob_solve602(2)Error using ==> prob_solve602输入量应是大于2的自然数prob_solve602(7.3)Error using ==> prob_solve602输入量应是大于2的自然数3 用泛函指令fminbnd求|]y x--=在x=0附近的极小x|e)sin[cos(x值。
MATLAB程序设计与应用 实验答案 第六章 刘卫国
1.(1)A=randn(10,5)A =-0.4326 -0.1867 0.2944 -0.3999 -1.6041 -1.6656 0.7258 -1.3362 0.6900 0.2573 0.1253 -0.5883 0.7143 0.8156 -1.0565 0.2877 2.1832 1.6236 0.7119 1.4151 -1.1465 -0.1364 -0.6918 1.2902 -0.8051 1.1909 0.1139 0.8580 0.6686 0.5287 1.1892 1.0668 1.2540 1.1908 0.2193 -0.0376 0.0593 -1.5937 -1.2025 -0.9219 0.3273 -0.0956 -1.4410 -0.0198 -2.1707 0.1746 -0.8323 0.5711 -0.1567 -0.0592 B=mean(A)B =0.0013 0.2310 0.0253 0.3588 -0.4197 C=std(A)C =0.9034 0.8829 1.1898 0.7832 1.0821 (2)D=max(max(A))D =2.1832E=min(min(A))E =-2.1707(3)F=sum(A,2)F =-2.3288-1.32870.01056.2215-1.48953.36024.9201-3.6964-3.3998-0.3025G=sum(sum(A))G =1.9666(4)H=sort(A)H =-1.6656 -0.8323 -1.5937 -1.2025 -2.1707-1.1465 -0.5883 -1.4410 -0.3999 -1.6041-0.4326 -0.1867 -1.3362 -0.1567 -1.0565-0.0376 -0.1364 -0.6918 -0.0198 -0.92190.1253 -0.0956 0.2944 0.6686 -0.80510.1746 0.0593 0.5711 0.6900 -0.05920.2877 0.1139 0.7143 0.7119 0.21930.3273 0.7258 0.8580 0.8156 0.25731.1892 1.0668 1.2540 1.1908 0.52871.19092.1832 1.6236 1.2902 1.4151L=-sort(-A,2,'descend')L =-1.6041 -0.4326 -0.3999 -0.1867 0.2944-1.6656 -1.3362 0.2573 0.6900 0.7258-1.0565 -0.5883 0.1253 0.7143 0.81560.2877 0.7119 1.4151 1.6236 2.1832-1.1465 -0.8051 -0.6918 -0.1364 1.29020.1139 0.5287 0.6686 0.8580 1.19090.2193 1.0668 1.1892 1.1908 1.2540-1.5937 -1.2025 -0.9219 -0.0376 0.0593-2.1707 -1.4410 -0.0956 -0.0198 0.3273-0.8323 -0.1567 -0.0592 0.1746 0.57112.(1)x=[0 15 30 45 60 75 90];y=[0 0.2588 0.5000 0.7071 0.8660 0.9659 1.0000];x1=0:90;y1=interp1(x,y,x1,'spline')y1 =Columns 1 through 110 0.0175 0.0349 0.0524 0.0698 0.0872 0.1045 0.1219 0.1392 0.1564 0.1737Columns 12 through 220.1908 0.2079 0.2249 0.2419 0.2588 0.2756 0.2923 0.3090 0.3255 0.3420 0.3583Columns 23 through 330.3746 0.3907 0.4067 0.4226 0.4384 0.4540 0.4695 0.4848 0.5000 0.5150 0.5299Columns 34 through 440.5446 0.5592 0.5736 0.5878 0.6018 0.6157 0.6293 0.6428 0.6561 0.6691 0.6820Columns 45 through 550.6947 0.7071 0.7193 0.7313 0.7431 0.7547 0.7660 0.7771 0.7880 0.7986 0.8090Columns 56 through 660.8191 0.8290 0.8387 0.8480 0.8571 0.8660 0.8746 0.8829 0.8910 0.8987 0.9062Columns 67 through 770.9135 0.9204 0.9271 0.9335 0.9396 0.9454 0.9510 0.9563 0.9612 0.9659 0.9703Columns 78 through 880.9744 0.9782 0.9817 0.9849 0.9878 0.9904 0.9927 0.9946 0.9963 0.9977 0.9987Columns 89 through 910.9995 0.9999 1.0000x=[0 15 30 45 60 75];y=[0 0.2679 0.5774 1.0000 1.7320 3.7320];x1=0:75;y1=interp1(x,y,x1,'spline')y1 =Columns 1 through 110 0.0184 0.0365 0.0545 0.0724 0.0902 0.1079 0.1255 0.1431 0.1607 0.1784Columns 12 through 220.1961 0.2138 0.2317 0.2497 0.2679 0.2863 0.3048 0.3236 0.3427 0.3620 0.3817Columns 23 through 330.4017 0.4221 0.4429 0.4641 0.4858 0.5079 0.5305 0.5537 0.5774 0.6017 0.6266Columns 34 through 440.6520 0.6780 0.7046 0.7317 0.7593 0.7876 0.8163 0.8456 0.8754 0.9058 0.9367Columns 45 through 550.9681 1.0000 1.0325 1.0658 1.1003 1.1364 1.17431.2145 1.2572 1.3028 1.3516Columns 56 through 661.4041 1.4604 1.5211 1.5863 1.6565 1.7320 1.8131 1.9002 1.99362.0937 2.2008Columns 67 through 762.3152 2.4374 2.5675 2.7060 2.85323.0095 3.17523.3506 3.5361 3.7320x=[0 15 30 45 60 75 90];y=[0 0.2588 0.5000 0.7071 0.8660 0.9659 1.0000];y1=polyfit(x,y,5)y1 =0.0000 0.0000 -0.0000 0.0000 0.0174 0.0000x=[0 15 30 45 60 75];y=[0 0.2679 0.5774 1.0000 1.7320 3.7320];y1=polyfit(x,y,5)y1 =0.0000 -0.0000 0.0000 -0.0010 0.0245 0.0000(2)x=[1 4 9 16 25 36 49 64 81 100];y=1:10;x1=1:100;y1=interp1(x,y,x1,'cubic')y1 =Columns 1 through 111.0000 1.3729 1.71252.0000 2.2405 2.4551 2.64942.82923.0000 3.1636 3.3186Columns 12 through 223.4661 3.6069 3.7422 3.87294.0000 4.1237 4.24354.3599 4.4730 4.5832 4.6907Columns 23 through 334.7958 4.89885.0000 5.0993 5.1966 5.2921 5.38575.4777 5.5681 5.6570 5.7446Columns 34 through 445.8309 5.91606.0000 6.0829 6.1647 6.2454 6.32496.4035 6.4810 6.5577 6.6334Columns 45 through 556.7082 6.7823 6.8556 6.92817.0000 7.0712 7.14167.2113 7.2804 7.3487 7.4164Columns 56 through 667.4835 7.5500 7.6159 7.6812 7.7459 7.8102 7.8739 7.9372 8.0000 8.0623 8.1242Columns 67 through 778.1855 8.2464 8.3068 8.3668 8.4263 8.4854 8.5441 8.6024 8.6603 8.7178 8.7749Columns 78 through 888.8317 8.8881 8.9442 9.0000 9.0555 9.1107 9.16559.2201 9.2744 9.3284 9.3821Columns 89 through 999.4354 9.4884 9.5412 9.5935 9.6456 9.6973 9.7486 9.7996 9.8502 9.9005 9.9505Column 10010.00003.xi=[165 123 150 123 141];yi=[187 126 172 125 148];P=polyfit(xi,yi,3)P =1.0e+003 *-0.0000 0.0013 -0.1779 8.4330线性拟合曲线为:p(x)=1.3x^2—177.9x+84334.(1)P1=[0,3,2];P2=[5,-1,2];P3=[1,0,-0.5];P=conv(P1,conv(P2,P3))P =0 15.0000 7.0000 -3.5000 0.5000 -2.0000 -2.0000 (2)P1=[0,3,2];P2=[5,-1,2];P3=[1,0,-0.5];P=conv(P1,conv(P2,P3));Y=roots(P)Y =0.70710.1000 + 0.6245i0.1000 - 0.6245i-0.7071-0.6667(3)P1=[0,3,2];P2=[5,-1,2];P3=[1,0,-0.5];P=conv(P1,conv(P2,P3));k=0:10;xi=0.2*k;Y=polyval(P,xi)Y =-2.0000 -2.3920 -2.6112 -1.7024 2.7104 15.0000 42.1120 94.1408 184.9056 332.5264 560.00005.(1)[U,fmin]=fminsearch('xiti651',[1,1])U =1.0e-004 *-0.0675 0.1715fmin =1.9920e-010(2)f=inline('-sin(x)-cos(x.^2)');fminbnd(f,0,pi)ans =0.73106.(1)x=[pi/6 pi/4 pi/3 pi/2];f=inline('sin(x).^2+cos(x).^2');dx=diff(f([x,5*pi/12]))/(pi/12)dx =0 0 0 0(2)x=1:3;f=inline('sqrt(x.^2+1)');dx=diff(f([x,4]))dx =0.8219 0.9262 0.96087.(1)g=inline('sin(x).^5.*sin(5*x)');I=quadl(g,0,pi)I =0.0982(2)g=inline('(1+x.^2)./(1+x.^4)');I=quad(g,-1,1)I =2.2214(3)g=inline('(x.*sin(x))./(1+cos(x).^2)'); I=quadl(g,0,pi)I =2.4674(4)f=inline('abs(cos(x+y))');dblquad(f,0,pi,0,pi)ans =6.28329.(1)矩阵求逆法:A=[2,3,5;3,7,4;1,-7,1];b=[10;3;5];x=inv(A)*bx =-1.8060-0.53733.0448矩阵除法:A=[2,3,5;3,7,4;1,-7,1];b=[10;3;5];x=A\bx =-1.8060-0.53733.0448矩阵分解法:A=[2,3,5;3,7,4;1,-7,1];b=[10;3;5];[Q,R]=qr(A);x=R\(Q\b)x =-1.8060-0.53733.0448(2)矩阵求逆法:A=[5,1,-1,0;1,0,3,-1;-1,-1,0,5;0,0,2,4]; b=[1;2;3;-1];x=inv(A)*bx =1.4000-5.90000.1000-0.3000矩阵除法:A=[5,1,-1,0;1,0,3,-1;-1,-1,0,5;0,0,2,4]; b=[1;2;3;-1];x=A\bx =1.4000-5.90000.1000-0.3000矩阵分解法:A=[5,1,-1,0;1,0,3,-1;-1,-1,0,5;0,0,2,4];b=[1;2;3;-1];[Q,R]=qr(A);x=R\(Q\b)x =1.4000-5.90000.1000-0.300010.A=[2 1 -1 1;4 2 -2 1;2 1 -1 -1];b=[1;2;1];[x,y]=line_solution(A,b)原方程组有无穷个解,特解为x,齐次方程组的基础解系为yWarning: Rank deficient, rank = 2, tol = 4.3512e-015. > In line_solution at 11方程组无解x =[]y =-0.5000 0.50001.0000 00 1.00000 011.(1)f=inline('x-sin(x)./x');x=fzero(f,0.5)x =0.8767(2)f=inline('(sin(x).^2).*exp(-0.1.*x)-0.5.*abs(x)');x=fzero(f,1.5)x =1.673812.x=fsolve('xiti612',[0.5,0.5],optimset('Display','off')) x =0.5000y =0.5000f =-0.0509f =-0.0509 0.1173 x =0.5000y =0.5000f =-0.0509f =-0.0509 0.1173 x =0.5000y =0.5000f =-0.0509f =-0.0509 0.1173 x =0.6459y =0.3739f =0.0055f =0.0055 0.0044 x =0.6459y =0.3739f =0.0055f =0.0055 0.0044 x =0.6459y =0.3739f =0.0055f =0.0055 0.0044 x =0.6355y =0.3734f =1.9417e-005f =1.0e-004 *0.1942 0.2589 x =0.6355y =0.3734f =1.9425e-005f =1.0e-004 *0.1942 0.2589x =0.6355y =0.3734f =1.9419e-005f =1.0e-004 *0.1942 0.2591x =0.6354y =0.3734f =2.3754e-010f =1.0e-009 *0.2375 0.2957x =0.6354y =0.3734f =7.9432e-009f =1.0e-008 *0.7943 0.5602x =0.6354y =0.3734f =1.8684e-009f =1.0e-007 *0.0187 0.1936x =0.6354 0.373413.x0=0;xf=20;y1=0;[x,y]=ode45('xiti6131',[x0,xf],y1) x =0.50001.00001.50002.00002.50003.00003.50004.00004.50005.00005.50006.00006.50007.00007.50008.00008.50009.00009.500010.000010.500011.000011.500012.000012.500013.000013.500014.000014.500015.000015.500016.000016.500017.000017.500018.000018.500019.000019.500020.0000 y =NaNNaNNaNNaNNaNNaNNaNNaNNaNNaNNaNNaNNaNNaNNaNNaNNaNNaNNaNNaNNaNNaNNaNNaNNaNNaNNaNNaNNaNNaNNaNNaNNaNNaNNaNNaNNaNNaNNaNNaN(2)x0=0;xf=20;[x,y]=ode45('xiti6132',[x0,xf],[-3,2])[x,y]此题结果过长,疑似错误或不会做。
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度。
- 1、下载文档前请自行甄别文档内容的完整性,平台不提供额外的编辑、内容补充、找答案等附加服务。
- 2、"仅部分预览"的文档,不可在线预览部分如存在完整性等问题,可反馈申请退款(可完整预览的文档不适用该条件!)。
- 3、如文档侵犯您的权益,请联系客服反馈,我们会尽快为您处理(人工客服工作时间:9:00-18:30)。
控制系统仿真第六章课后题作业6.1在图6.1中,已知单位负反馈系统被控对象的传递函数为)1001.0)(11.0()(++=s s s K s G 试编写matlab 程序,设计系统的超前矫正器Gc(s),要求:1)在斜坡信号r (t)=2t 作用下,系统的稳态误差ess<=0.002;2)校正后系统的相位裕度Pm 范围为:45~55;3)绘制系统校正后的bode 图和阶跃响应曲线。
程序:>> s=tf('s');>> G=1000/(s*(0.1*s+1)*(0.001*s+1));>> margin(G) % 绘制校正前的bode 图>> figure(2)>> sys=feedback(G,1);>> step(sys) %绘制校正前的单位阶跃响应曲线>> [Gm,Pm]=margin(G); %该句值计算bode图的增益裕量Gm和相位裕量Pm >> [mag,phase,w]=bode(G); %该句只计算bode图上多个频率点w对应的幅值和相位>> QWPm=50; %取矫正后的相位为50>> FIm=QWPm-Pm+5;>> FIm=FIm*pi/180;>> alfa=(1-sin(FIm))/(1+sin(FIm));>> adb=20*log10(mag);>> am=10*log10(alfa);>> wc=spline(adb,w,am);>> T=1/(wc*sqrt(alfa));>> alfat=alfa*T;>> Gc=tf([T 1],[alfat 1]) %校正器的传递函数Transfer function:0.01794 s + 1-------------0.00179 s + 1>> figure(3)>> margin(Gc*G) % 系统矫正后的bode图>> figure(4)>> step(feedback(Gc*G,1)) % 校正后的单位阶跃响应曲线作业6.2在图6.1中,已知单位负反馈系统被控对象的传递函数为)102.0)(11.0()(++=s s s K s G 试编写matlab 程序,设计系统的滞后校正器Gc (s ),要求:1) 在斜坡信号r(t)=t 作用下,系统的稳态误差ess 01.0≤;2) 校正后系统的相位裕度Pm 范围为:40~50;3) 绘制系统矫正前后的bode 图和阶跃响应曲线。
计算K ess=1/K<=0.01 K>=100 取K=100程序如下:>> s=tf('s');>> G=100/(s*(0.1*s+1)*(0.02*s+1));>> margin(G) %绘制校正前的bode 图如下图示,这是一个不稳定的系统。
>> figure(2)>> step(feedback(G ,1))可以看出系统的动态响应不稳定,处于震荡发散的状态。
>> P0=45;>> fic=-180+P0+5;>> [mu,pu,w]=bode(G);>> wc2=spline(pu,w,fic);>> d1=conv(conv([1 0],[0.1 1]),[0.02 1]);>> na=polyval(100,j*wc2);>> da=polyval(d1,j*wc2); %该句是用j*wc2代替多项式中的s >> G1=na/da;>> g1=abs(G1); %abs()函数是取绝对值,这里表示求复数G1的模>> L=20*log10(g1);>> beta=10^(L/20); T=1/(0.1*wc2);>> betat=beta*T;>> Gc=tf([T 1],[betat 1])Transfer function:1.558 s + 1-----------20.27 s + 1>> figure(3)>> margin(Gc*G)>> figure(4)>> step(feedback(Gc*G,1))作业6.3在图6.1中,已知单位负反馈系统被控对象的传递函数为)2)(1()(++=s s s K s G 试编写matlab 程序,设计系统的滞后-超前校正器Gc(s),要求:1) 在单位斜坡信号r(t)=t 作用下,系统的速度误差系数Kv=20s-1;2) 校正后系统的实际相位裕量Pm 范围为:42~58;3) 校正后系统的实际剪切频率wc2>=1.3rad/s;4) 绘制系统校正前后的bode 图和阶跃响应曲线。
由题知 202)2)(1(lim )(lim 00==++==→→K s s s K s s sG Kv s s K=40 程序如下:先绘制校正前的bode 图和单位阶跃响应曲线,再绘制校正后的 第一种>> s=tf('s');>> G0=40/(s*(s+1)*(s+2));>> figure(1)>> margin(G0)>> figure(2)>> step(feedback(G0,1))校正前系统的增益裕量Gm=-16.5,相位裕量Pm=-40.4都是负值,系统不稳定 剪切频率wc1=3.19rad/s,阶跃响应曲线是发散的。
>> wc2=4;>> [Gm,Pm,wc1]=margin(G0);>> beta=9;T1=1/(0.1*wc1);>> betat=beta*T1;>> Gc1=tf([T1 1],[betat 1]) %计算并显示滞后校正器传递函数Transfer function:7.071 s + 1-----------63.64 s + 1>> sope=G0*Gc1; %计算原系统与滞后校正器串联后的传递函数>> num=sope.num{1};den=sope.den{1};>> na=polyval(num,j*wc2);>> da=polyval(den,j*wc2);>> G=na/da;>> g1=abs(G);>> L=20*log10(g1);>> alfa=10^(L/20);>> T=1/(wc2*(alfa)^(1/2));>> alfat=alfa*T;>> Gc2=tf([T 1],[alfat 1])Transfer function:1.018 s + 1-------------0.06139 s + 1>> G=G0*Gc1*Gc2; >> sys=feedback(G,1); >> figure(3)>> margin(G)>> figure(4)>> step(sys)第二种:改进后>> s=tf('s');>> G0=40/(s*(s+1)*(s+2));>> figure(1)>> margin(G0)>> figure(2)>> step(feedback(G0,1))>> wc2=10;>> [Gm,Pm,wc1]=margin(G0);>> beta=9;T1=1/(0.1*wc1);>> betat=beta*T1;>> Gc1=tf([T1 1],[betat 1]); %¼ÆËã²¢ÏÔʾÖͺóУÕýÆ÷´«µÝº¯Êý>> n1=conv([T1 1],[0 40]);>> d1=conv(conv(conv([1 0],[1 1]),[1 2]),[betat 1]);>> na=polyval(n1,1i*wc2);>> da=polyval(d1,1i*wc2);>> G=na/da;>> g1=abs(G);>> L=20*log10(g1);>> alfa=10^(L/20);>> T=1/(wc2*(alfa)^(1/2)); >> alfat=alfa*T;>> Gc2=tf([T 1],[alfat 1]); >> G=G0*Gc1*Gc2;>> sys=feedback(G,1);>> figure(3)>> margin(G)>> figure(4)>> step(sys)作业6-4在图6.1中,已知受控对象为一个带延迟的惯性环节,其传递函数为 e s s s G 71176)(-+= 试编写matlab 程序,用Ziegler-Nichols 经验整定公式中的两种方法,分别计算P ,PI ,PID 控制器的参数,并进行阶跃响应仿真。
由传递函数知:比例系数K=6,惯性时间常数T=17.纯延迟时间常数t=7.可采用Ziegler-Nichols 经验整定公式中阶跃响应整定法。
ττττττ5.0,2,2.1,3,9.0,======PIDTd PIDTi K T PIDKp PITi K T PIKp K T PKp >> K=6;T=17;tau=7;>> s=tf('s');>> Gz=K/(T*s+1);>> [np,dp]=pade(tau,2);>> G1=tf(np,dp);>> G=Gz*G1;>> Pkp=T/(K*tau) % 计算比例控制器的KpPkp =0.4048>> step(feedback(Pkp*G,1)),hold on>> PIKp=0.9*T/(K*tau);>> PITi=3*tau; %计算比例积分控制器的Kp和Ti >> PIGc=PIKp*(1+1/(PITi*s))Transfer function:7.65 s + 0.3643---------------21 s>> step(feedback(PIGc*G,1)),hold on>> PIDKp=1.2*T/(K*tau);>> PIDTI=2*tau;>> PIDTd=0.5*tau;>> PIDGc=PIDKp*(1+1/(PIDTI*s)+PIDTd*s/((PIDTd/10)*s+1))Transfer function:26.18 s^2 + 6.97 s + 0.4857---------------------------4.9 s^2 + 14 s>> step(feedback(PIDGc*G,1)),hold on>> [PIDKp,PIDTI,PIDTd]ans =0.4857 14.0000 3.5000>> gtext('p'); %在曲线上打印文字>> gtext('PI');>> gtext('PID');。