《《自动控制理论》课程设计指导应用matlab.doc》由会员分享，可在线阅读，更多相关《《自动控制理论》课程设计指导应用matlab.doc（53页珍藏版）》请在taowenge.com淘文阁网|工程机械CAD图纸|机械工程制图|CAD装配图下载|SolidWorks_CaTia_CAD_UG_PROE_设计图分享下载上搜索。
1、【精品文档】如有侵权，请联系网站删除，仅供学习与交流自动控制理论课程设计指导应用matlab.精品文档.自动控制理论课程设计指导书霍爱清 薛朝妹自动化教研室目 录1 控制系统的数学描述11.1 微分方程11.2 传递函数61.3 状态空间描述91.4 模型转换112 控制系统的校正152.1 单变量系统的两种主要校正方式152.2 PI、PD、PID校正152.3 串联校正举例183 MATLAB在自动控制系统中的应用243.1 概述243.2 控制系统函数全集303.3 应用实例314 SIMULINK简介414.1 引导414.2 SIMULINK模型的构建455 自动控制系统设计任务50
2、5.1 任务一505.2 任务二525.3 任务三535.4 任务四546 附录566.1 时域分析例题及程序566.2 根轨迹分析例题及程序646.3 频域分析例题及程序66第一章 控制系统的数学描述1.1 微分方程1.1.1 物理系统的微分方程利用机械学 、电学、流体力学和热力学等的物理规律，我们可以得到物理系统的动态方程。它们通常用常系数线性微分方程来描述。1.1.2 数值解通过拉普拉斯变换和反变换，可得到线性时不变方程的解析解，也可用状态转移矩阵 (t)求解。这些分析方法通常只限于常系数的线性微分方程。解析解是精确的，然而通常寻找解析解是困难的，甚至是不可能的。而数值分析方法直接在时域
3、里求解微分方程，不仅适用于线性时不变方程，也适用于非线性以及时变微分方程。MATLAB提供了两个求微分方程数值解的函数，它们采用龙格-库塔(Runge-kutta)法。Ode23和ode45分别表示采用2阶和4阶龙格-库塔公式，后者具有更高的精度。n阶微分方程必须化为n个首1的一阶微分方程组，且放入M-文件中，以便返回方程状态变量的导数，下面的例子介绍这些函数的用法。 例1.1 对图1-1的机械系统，已知三个量拉力、摩擦力、以及弹簧力都影响质量M的加速度。利用牛顿运动定理，建立系统的力平衡方程式令 ，有设质量M=1kg，摩擦系数B=5N/m/sec，弹簧常数K=25N/m。在t=0时刻，施加2
4、5N的拉力。上述方程及已知量在M-文件mechsys.m中定义如下：function xdot=mechsys(t, x);F=25;M=1;B=5;K=25;xdot=x(2);1/M*(F-B*x(2)-K*x(1);下面的M-文件使用ode23对系统在零初始条件下进行仿真：t0=0; tfinal=3; 时间间隔03秒x0=0,0; 零初始条件tol=0.001; 精度trace=0; 如果非零，则打印出每一步的计算值t, x=ode23(mechsys,t0,tfinal,x0,tol,trace)subplot(211),plot(t, x);title (Time response
5、 of mechanical translational system)xlabel (Time-sec)text (2,1.2,displacement)text (2,.2,veloclty)d=x(:,1);v=x(:,2);subplot(212),plot(d,v);title (velocity versus displacement)xlabel (displacement)ylabel (velocity)subplot(111)仿真结果如图1-2。图1-2例1.2 电路图见图1-3。R=1.4，L=2H，C=0.32F。初始状态：电感电流为零，电容电压为0.5V。时刻接入伏的
6、电压。求0tP=1 9 31.25 61.25 67.75 14.75 15;r=roots(P)多项式的根从列向量r中得到r =-4.0000-3.0000-1.0000 + 2.0000i-1.0000 - 2.0000i-0.0000 + 5.0000i-0.0000 - 5.0000i例1.5 多项式的根为-1,-2,-3j4。求多项式方程。为了输入复数，必须首先建立虚数单位。然后在行/列向量中输入根。使用poly得到多项式方程。 i=sqrt(-1); r=-1;-2;-3+4*i;-3-4*i; p=poly(r)多项式的系数从行向量中得到p =1 9 45 87 50因此，多项式
7、方程为s4+9s3+45s2+87s+50=0例1.6 求下列矩阵的特征方程的根用ploy求矩阵的特征方程的根。用roots求方程的根。 A=0 1 -1;-6 -11 6;-6 -11 5; P=poly(A) r=roots(P)结果如下：P= 1.0000 6.0000 11.0000 6.0000r =-3.0000-2.0000-1.00001.2.2 传递函数的零点和极点l 函数tf2zp求传递函数的零点，极点和增益。例1.7 求下列传递函数的零点，极点和增益。 num=1 11 30 0; den=1 9 45 87 50； z, p, k=tf2zp(num, den)z =0
8、-6.0000-5.0000p =-3.0000 + 4.0000i-3.0000 - 4.0000i-2.0000-1.0000因而有l 函数zp2tf根据给定零点，极点和增益求传递函数。例1.8 系统的零点为-6,-5,0,极点为-3j4,-2,-1,增益为1。求其传递函数。 z=-6;-5;0;k=1; i=sqrt(-1); p=-3+4*i;-3-4*i;-2;-1; num, den=zp2tf(z, p, k)上面程序的结果为num =0 1 11 30 0den =1 9 45 87 50因此，传递函数为1.2.3 部分分式展开函数r, p, k=residue(b, a)，对
9、两个多项式的比进行部分分式展开，如 (1.3)向量b, a以s的降幂顺序排列多项式的系数。部分分式展开后余数送入列向量r，极点送入列向量p，常数项送入k。例1.9 对F(s)进行部分分式展开 b=2 0 9 1; a=1 1 4 4; r, p, k=residue(b, a)结果如下：r =0.0000 - 0.2500i0.0000 + 0.2500i-2.0000p =-0.0000 + 2.0000i-0.0000 - 2.0000i-1.0000k =2因而，部分分式展开为函数b, a=residue(r, p, k)将部分分式转化为多项式比P(s)/Q(s)。1.3 状态空间描述集
10、总参数的线性网络可用微分方程表示为 (1.4)该系统的一阶微分方程即为状态方程，X是状态向量。状态空间方法易采用数字或模拟计算机求解。另外，状态空间方法容易拓展到非线性系统。状态方程可从n阶微分方程得到，或者在系统模型中选用合适的状态变量直接写出。1.3.1 将微分方程化成状态方程设一个n阶线性系统由微分方程描述，我们讨论如何选择状态变量，使该系统由状态方程描述。 (1.5)y(t),u(t)分别为系统的输出、输入。状态模型不是唯一的，它依赖于状态变量的选择。一个有用的选择方法是定义那么，y的n阶导数可用各状态变量的线性组合来表示。于是有写成矩阵形式为 (1.7)输出方程为y=1 0 0 0x
11、 (1.8)即 (1.9) D=0例1.10 给定系统的微分方程描述为y(3)+2y(2)+3y(1)+4y=5u(t)利用(1.9)即可写出其相应的一个状态空间描述为：上述状态方程，称为能控标准型。1.3.2 矩阵的对角化为什么要研究A矩阵的对角化问题呢？原因之一是将矩阵对角化后，对角元即为矩阵A的特征值1,2,n。从而，状态转移矩阵亦可化为对角元为e1t,e2t,ent的对角矩阵。（假定A有n个互不相同的特征值）。给定线性系统。当A有n个互不相同的特征值时，我们可以找到非奇异变换矩阵P，令x(t)=py(t) (1.10)将上面状态方程化为对角线规范型 (1.11)其中a=p-1AP, b
12、=p-1B有很多方法可以求得P，如利用A的特征向量可构造P。例1.11 定系统的状态空间描述为求它的对角规范型和变换矩阵P。 A=0 1 -1;-6 -11 6;-6 -11 5; B=0;0;1 P,L=eig(A); L是一个对角元为特征值的对角矩阵P是一个变换矩阵，其列是相应于特征值的特征向量 a=inv(P)*A*P A矩阵对角化 b=inv(P)*B结果为：p = -0.7071 -0.2182 -0.0921 -0.0000 0.4364 -0.5523 -0.7071 0.8729 -0.8285a = -1.0000 -0.0000 -0.0000 0.0000 -2.0000
13、 -0.0000 -0.0000 0.0000 -3.0000b = 2.8284 13.7477 10.86281.4 模型转换1.4.1 传递函数向状态空间描述的转换控制系统工具箱包含一组模型转换的函数。A,B,C,D=tf2ss(num, den)将传递函数转换为状态空间描述。例1.12 求下面传递函数的状态空间描述 num=1 7 2; den=1 9 26 24; A, B, C, D=tf2ss(num, den)状态方程各矩阵如下： D=01.4.2 状态空间描述向传递函数的转换已知状态方程和输出方程 (1.12)y=Cx+Du (1.13)采用拉普拉斯变换有Y(s)=C(sI-
14、A)-1Bu(s)+Du(s)则 (1.14)函数ss2tf(A,B,C,D,i)是将状态空间描述转换为对第一个输入的传递函数。num,den=ss2tf(A,B,C,D,i)是将状态空间描述化为分子、分母多项式形式的传递函数。z,p,k=ss2zp(A,B,C,D,i)将状态空间描述化为零极点形式表示的传递函数。例1.13 一个系统的状态空间描述如下y=1 0 0x求传递函数G(s)=Y(s)/U(s)A=0 1 0; 0 00 1; -1 -2 -3; B=10; 0; 0;C=1 0 0;D=0;num,den=ss2tf(A,B,C,D,1)z,p,k=ss2zp(A,B,C,D,1)
15、其中，ss2tf(A,B,C,D,1)中“1”表示对第一个输入。传递函数的分子、分母多项式系数如下：num=0 10.0000 30.0000 20.0000den=1.0000 3.0000 2.0000 1.0000传递函数的零、极点如下：z=-1-2p=-0.3376+0.5623i-0.3376-0.5623i-2.3247k=10因而传递函数为1.4.3 由方框图求状态空间描述和传递函数控制系统工具箱中提供了函数A,B,C,D=connect(a, b, c, q, iu, iy)。将方框图描述转换成状态空间描述和传递函数。其中q矩阵规定了各框之间的连接关系。其每一行的第一个元素是框
16、号，其余的元素依次是于该框连接的框号，iu,iy分别表示输入，输出施加的框号。例1.14 将图1-7由框图表示的系统转换成状态空间描述和传递函数。n1=1; d1=1; n2=0.5; d2=1; n3=4; d3=1 4;n4=1; d4=1 2; n5=1; d5=1 3; n6=2; d6=1;n7=5; d7=1; n8=1; d8=1;nblocks=8; blkbuildq=1 0 0 0 0 q矩阵表示框图的结构。2 1 -6 -7 -8 如第2个框于第1个框按3 2 0 0 0 1的关系连接，于第6.7.84 3 0 0 0 个框按-1关系连接，依次类推。5 4 0 0 06
17、3 0 0 07 4 0 0 08 5 0 0 0;iu=1; 输入施加于第1个框上iy=8; 由第8个框输出A, B, C, D=connect(a, b, c, d, q, iu, iy)num, den=ss2tf(A,B,C,D,1) 转换成传递函数结果为A=-8.0 -2.5 -0.50.4 -2.0 00 1.0 -3.0B=0.500C=0 0 1D=0num=0 0 0 2den=1.0 13.0 56.0 80.0即第二章 控制系统的校正2.1 单变量系统的两种主要校正方式单变量系统常用的校正主要有两种方式。一种是校正装置与被控对象串联，如图2-1所示。这种校正方式称为串联校
18、正。另一种校正方式是从被控对象中引出反馈信号，与被控对象或其一部分构成反馈回路，并在局部反馈回路设置校正装置。这种校正方式称为局部反馈校正，如图2-2所示。图2-1图2-2串联校正和局部反馈校正应用都相当普遍，究竟选择哪一种，取决于系统中信号的性质，可供采用的元件以及其他条件。两种校正方式结合起来可以收到更好的效果。2.2 PI、PD、PID校正2.2.1 PD校正（超前校正）超前校正（亦称PD校正）的传递函数为 （2.1）其对数频率特性如图2-3所示。超前校正能够产生相位超前角，它的强度可由参数表征。超前校正的相频特性函数是（）=arctgT-arctgT （2.2）最大相移点位于对数频率的
19、中心点，即 (2.3)最大相移量为 （2.4）或者 （2.5）容易求出，在点有L()=10lg (2.6)图2-3基于频率法综合超前校正的步骤是：1. 首先根据静态指标要求，确定开环比例系数K，并按已确定的K画出系统固有部分的Bode图。2. 根据静态指标要求预选c，从Bode图上求出系统固有部分在c点的相角。3. 根据性能指标要求的相角裕量，确定在c点是否需要提供相角超前量。如需要，算出需要提供的相角超前量m。4. 如果所需相角超前量不大于60,按（2.5）求出超前校正强度。5. 令，从而求出超前校正的两个转折频率1/T和1/T。6. 计算系统固有部分在c点的增益Lg(dB)及超前校正装置在
20、c点的增益Lc(dB)。如果Lg+Lc0，则校正后系统的截止角频率c比预选的值要高。如果高出较多，应采用滞后超前校正，如果只是略高一些，则只需核算c点的相角裕量，若满足要求，综合完毕，否则转第3步。如果Lg+Lc0，令Lg(c)=20lg，求出，就是滞后校正的强度，如果Lg(c)0，则无须校正，且可将开环比例系数提高。4. 选择2=1/T=（1/51/10）c，进而确定1=1/(T)。5. 画出校正后系统的Bode图，校核相位裕量。滞后校正的主要作用是降低中频段和高频段的开环增益，但同时使低频段的开环增益不受影响，从而达到坚固静态性能和稳定性。它的副作用是会在c点产生一定的相角滞后。2.2.3
21、 PID校正（滞后超前校正）超前校正的主要作用是增加相角裕量，改善系统的动态响应特性。滞后校正的主要作用是改善系统的静态特性，两种校正结合起来就能同时改善系统的动态和静态特性特性。滞后超前校正（亦称PID校正）综合了前面两种校正的功能。滞后超前校正的传递函数为 (2.8)它相当于一个滞后校正与一个超前校正相串联，其对数频率特性如图2-5所示。图2-5基于频率法的滞后超前校正的综合步骤是：1. 首先根据静态指标要求确定开环比例系数，按照所确定的K画出系统固有部分的Bode图。2. 按指标要求确定c，检查系统固有部分在c的对数幅频的斜率是否为-2，如果是，求出c点的相角。3. 按综合超前校正的步骤
22、36综合超前部分Gc1(s)(注意在确定时要计入滞后校正将造成的50120的相角滞后量)。在第6步时注意，通常Lg(c)+Lc(c)比0高出很多，所以要引进滞后校正。4. 令20lg= Lg(c)+Lc(c),求出。5. 按综合校正的步骤45综合滞后部分Gc2(s)。6. 将滞后校正与超前校正串联在一起，构成滞后超前校正：Gc(s)=Gc1(s)*Gc2(s)。2.3 串联校正举例对于一结构如图2-6所示的系统，给定固有部分的传递函数Gg(s)和性能指标要求，试设计串联校正装置K(s)。设1. 若要求开环比例系数K100 1/秒，相角裕量30，c45 1/秒。2. 若要求开环比例系数K100
23、1/秒，相角裕量40，c=5 1/秒。3. 若要求开环比例系数K100 1/秒，相角裕量40，c20 1/秒，分别设计校正装置，并比较它们的作用效果。图2-6对于1，需进行超前校正，这里给出一个参考的K(s)=(0.05s+1)/(0.005s+1)，于是系统的开环传递函数为闭环传递函数为其单位阶跃响应的实现为：num=1000,000 20,000,000;den=1 310 23,000 1200,000 20,000,000;t=xx: xx: xx;（起始时间，步长，终止时间）c=step(num, den, t);plot(t, c);xlabel (t-sec.), ylabel(
24、c(t),grid, pause对于2需进行滞后校正，这里给出一参考的仿照上述，可求出系统的阶跃响应。对于3，需进行滞后超前校正，这里给出一参考的同理可得到系统的阶跃响应。比较上述三种不同的校正装置K(s)所产生的不同效果，可以看出三种校正方式的作用。这里需要特别指出阶跃响应的起始时间，步长及终止时间的给定问题。起始时间通常给0，终止时间要根据系统的过渡过程时间的长短来给定。我们可按估计过渡过程时间的经验公式ts(49)/c给定此值：如2，c=5 1/秒，则过渡过程时间大约在0.81.8秒之间，于是终止时间可给2秒，我们希望在0至2秒内显示出20个点，则步长为0.1秒。若已知系统的c，建议步长
25、按1/(2c)给定。对c未知，也不便估计过渡过程时间的系统只好用试探法，待曲线显示出以后再调整步长和终止时间。图2-7为系统固有部分的频率特性图（幅频特性和相频特性）可以看出该系统已接近临界稳定状态了。图2-8至图2-11分别表示不加校正，加超前校正，加滞后校正，和加滞后超前校正的系统阶跃响应。图2-7num=100;den=conv(1,0, conv(0.1,1, 0.01,1);bode(num, den)图2-8numo=100;deno=conv(1,0,conv(0.1,1,0.01,1);t=0:0.01:2;numc=numo;denc=zeros(1,length(deno)
26、-length(numo),numo+deno;c=step(numc,denc,t);plot(t,c)xlabel(t-sec),ylabel(c(t),grid,pause图2-9numo=conv(100,0.05,1);deno=conv(0.005,1,conv(1,0,conv(0.1,1,0.01,1);t=0:0.01:2;numc=numo;denc=zeros(1,length(deno)-length(numo),numo+deno;c=step(numc,denc,t);plot(t,c)xlabel(t-sec),ylabel(c(t),grid,pause图2-1
27、0numo=conv(100,0.5,1);deno=conv(10,1,conv(1,0,conv(0.1,1,0.01,1);t=0:0.01:2;numc=numo;denc=zeros(1,length(deno)-length(numo),numo+deno;c=step(numc,denc,t);plot(t,c)xlabel(t-sec),ylabel(c(t),grid,pause图2-11numo=conv(conv(0.25,1,0.1,1),100);deno=conv(conv(1.25,1,0.02,1),conv(1,0,conv(0.1,1,0.01,1);t=0
28、:0.01:2;numc=numo;denc=zeros(1,length(deno)-length(numo),numo+deno;c=step(numc,denc,t);plot(t,c)xlabel(t-sec),ylabel(c(t),grid,pause第三章 MATLAB在控制理论中的应用3.1 概述MATLAB提供了大量的控制工程计算、设计库函数。其中，控制系统软件包包括复数运算、特征值计算、方程求解、矩阵变换以及FFT等重要计算工具及举例。MATLAB的线性代数处理，矩阵运算和数值分析的能力为控制系统工程设计及其它学科研究提供了可靠的基础和强有力的研究工具。控制系统软件包利用M
29、ATLAB矩阵功能提供了适用于控制工程的专用函数，这些函数大部分用M文件表示。控制系统软件包可以方便地用于控制系统设计、分析和建模。 在控制系统软件包中，控制系统通常采用传递函数与状态空间两种形式建模，允许“经典”和“现代”技术并用，既可处理连续时间系统也可处理离散时间系统，并且可以进行不同模型表示形式之间的相互转换，也可以计算和绘制时间响应、频率响应及根轨迹图。此外M文件还能够进行极点配置和最优控制器的参数计算。即使在软件包中没有提供的功能，也可以通过编写新的M文件方式来构造。3.1.1系统模型 控制系统软件包可用于线性时不变（简称LTI）系统模型。时不变系统模型包括：状态空间模型；传递函数模型；零一极点增益模型；部分分式模型；离散时间模型。 下面介绍如何用矩阵表示不同类型的线性时不变系统模型。 1状态空间模型 一个LTI微分方程系统总可以用一组一阶微分方程组来表示。按矩阵或状态空间表示形式，LTI系统模型的一般形式如下 =AxBu y=CxDu其中：u是一个nu维控制输入向量，x是一个ns维状态向量，y是ny维输出向量。 采用MATLAB表示状态空间系统十分容易。A，B，C，D都是矩阵，均作为独立变量处理。 例