控制系统CAD难点.doc

研究生课程考试成绩单（试卷封面）院 系自动化学院专业控制工程学生姓名侯萌学号157175课程名称控制系统CAD授课时间 2016年 2 月至 2016年 6 月周学时36学分2简要评语考核论题总评成绩（含平时成绩）备注任课教师签名： 日期： 注：1. 以论文或大作业为考核方式的课程必须填此表，综合考试可不填。“简要评语”栏缺填无效。2. 任课教师填写后与试卷一起送院系研究生秘书处。3. 学位课总评成绩以百分制计分。被控对象的单位阶跃响应数据在文件matlab_work.mat中(要拷在MATLAB的work目录下才能看到)，t是采样时间点,y是对应于采样时间点处的单位阶跃响应数据。试辨识被控对象模型并用不同的控制方法（包括经典控制方法、现代控制方法和计算机控制方法）设计控制器，使对象的开环增益K=10，输出响应的超调量%=20%，调节时间ts=90%；4 程序及仿真结果。5 讨论控制器的鲁棒性和抗干扰性。1 原始数据分析被控对象的开环单位阶跃响应数据在文件matlab_work.mat中。在MATLAB指令窗中键入： load matlab_work.mat plot(t,y) xlabel(t) ylabel(y) title(原始响应曲线)在matlab中画出在输入为单位阶跃响应数据时输出随时间的变化曲线如下图所示：图1 被控对象的单位阶跃响应2 系统模型的辨识由于被控对象模型结构未知，题中仅给出了被控对象的单位阶跃响应。因此需要根据原始的数据分析得出原系统的基本结构。至于结构以及参数的最终确定，将使用MATLAB的系统辨识工具箱来实现。2.1系统模型结构的估计系统辨识工具箱提供的模型结构选择函数有struc、arxstruc和selstruc。函数struc 生成arx结构参数，调用格式为：nnstruc (Na, Nb, Nk)其中，Na、Nb 分别为arx模型多项式A(q)、B(q) 的阶次范围；Nk为arx模型纯滞后的大小范围；nn为模型结构参数集构成的矩阵。函数arxstruc用来计算arx模型结构的损失函数，即归一化的输出预测误差平方和，调用格式为：varxstruc (ze, zv, nn)其中，zey u为模型辨识的I/O数据向量或矩阵。zvyr ur为模型验证的I/O数据向量或矩阵。nn为多个模型结构参数构成的矩阵，nn的每行都具有格式：nnna nb nk。v的第一行为各个模型结构损失函数值，后面的各行为模型结构参数。函数selstruc 用来在损失函数的基础上进行模型结构选择，调用格式为：nn, vmodselstruc (v, c)其中，v 由函数arxstruc获得的输出矩阵，为各个模型结构的损失函数。c为可选参数，用于指定模型结构选择的方式。根据图2所示曲线的形状初步估测被控对象的模型应该为二阶系统或者更高阶系统。并且可以看出，纯滞后为0，故Nk恒为零。在MATLAB指令窗中键入： u=ones(size(y)Z=y,u v2=arxstruc(Z,Z,struc(2,0:2,0) nn2=selstruc(v2,0) v3=arxstruc(Z,Z,struc(3,0:3,0) nn3=selstruc(v3,0) v4=arxstruc(Z,Z,struc(4,0:4,0) nn4=selstruc(v4,0) v5=arxstruc(Z,Z,struc(5,0:5,0) nn5=selstruc(v5,0) v6=arxstruc(Z,Z,struc(6,0:6,0) nn6=selstruc(v6,0) v7=arxstruc(Z,Z,struc(7,0:7,0) nn7=selstruc(v7,0) v8=arxstruc(Z,Z,struc(8,0:8,0) nn8=selstruc(v8,0)得到如下结果：nn2 = 2 2 0nn3 = 3 3 0nn4 = 4 1 0nn5 = 5 5 0nn6 = 6 1 0nn7 = 7 5 0nn8 = 8 3 0于是，去除积分环节后的模型阶数为：二阶系统2 2 0、三阶模型3 3 0和四阶模型4 1 0，五阶系统模型为5 5 0，六阶系统模型为6 1 0，七阶系统模型为7 5 0，八阶系统模型为8 3 0。2.2系统模型结构的确定为了确定模型结构以及参数，使用MATLAB的系统辨识工具箱中已有的辨识函数arx()。辨识函数arx()的使用方法是：如果输入信号的列向量为u，输出信号的列向量为y，并选定了系统的分子多项式阶次m-1，分母多项式阶次n及系统的纯滞后d，则可以通过下面的指令辨识出系统的数学模型：T=arx(y,u, n,m,d)该函数将直接显示辨识的结果，且所得的T为一个结构体，其中T.A和T.B分别表示辨识得到的分子和分母多项式。由给定的观测数据建立系统数学模型后，还需要进行检验，看模型是否适用，如果不适用，则要修改模型结构，重新进行参数估计等。MATLAB的系统辨识工具箱中用于模型验证和仿真的函数主要有compare、resid、pe、predict 和idsim。此次实验主要用的是函数compare对模型进行验证。函数compare可将模型的预测输出与对象实际输出进行比较。验证过程与结果如下所示。 对二阶系统的验证在MATLAB指令窗中键入： Z=iddata(y,u,0.05) M=arx(Z,2,2,0) compare(M,Z)得到结果如图2所示：图2 二阶模型的匹配结果 对三阶系统的验证在MATLAB指令窗中键入： Z=iddata(y,u,0.05) M=arx(Z,3,3,0) compare(M,Z)得到结果如图3所示：图3 三阶模型的匹配结果 对四阶系统的验证在MATLAB指令窗中键入： Z=iddata(y,u,0.05) M=arx(Z,4,1,0) compare(M,Z)得到结果如图4所示：图4 四阶模型的匹配结果 对五阶系统的验证在MATLAB指令窗中键入： Z=iddata(y,u,0.05) M=arx(Z,5,5,0) compare(M,Z)得到结果如图5所示：图5 五阶模型的匹配结果 对六阶系统的验证在MATLAB指令窗中键入： Z=iddata(y,u,0.05) M=arx(Z,6,1,0) compare(M,Z)得到结果如图6所示：图6 六阶模型的匹配结果 对七阶系统的验证在MATLAB指令窗中键入： Z=iddata(y,u,0.05) M=arx(Z,7,5,0) compare(M,Z)得到结果如图7所示：图7 七阶模型的匹配结果 对八阶系统的验证在MATLAB指令窗中键入： Z=iddata(y,u,0.05) M=arx(Z,8,3,0) compare(M,Z)得到结果如图8所示：图8 八阶模型的匹配结果 模型阶次的选取从图28可以看出，八阶模型的拟合率已高达93%，符合要求。因此，为简便起见，不再选用更高阶系统，本次试验选用八阶系统模型。在MATLAB指令窗中键入： H=tf(M)得到如下结果：H = From input u1 to output y1: 0.4025 z-2 - 1 - 1.082 z-1 - 0.1049 z-2 + 0.3345 z-3 + 0.1519 z-4 + 0.2915 z-5 + 0.02519 z-6 - 0.6381 z-7 + 0.4242 z-8Sample time: 0.05 secondsDiscrete-time transfer function.上述模型为离散时间系统模型，为了得到连续时间系统模型，在MATLAB指令窗中键入： G=d2c(H(1)得到如下结果： G = From input u1 to output y1: -1.958 s7 - 71.66 s6 - 3434 s5 + 4.562e04 s4 + 6.38e06 s3 + 2.581e08 s2 + 5.035e09 s + 4.182e10 - s8 + 17.15 s7 + 4570 s6 + 5.992e04 s5 + 4.939e06 s4 + 4.392e07 s3 + 9.317e08 s2 + 3.74e09 s + 4.177e10Continuous-time transfer function.2.3传递函数模型系数的修正得到了系统传递函数，对其使用阶跃响应step命令，并与原数据线对比，看是否吻合。命令如下：plot(t,y,red);hold on;step(G)得到的结果：图9 估算后的曲线匹配从上图可以看出，由arx模型得到的传递函数的阶跃响应曲线并不能与原始数据线很好的拟合。故考虑在该传递函数的基础上，对其系数修改。1) 考虑忽略传递函数分子的高次项G=tf(5.035e09 4.182e10, 1 17.15 4570 5.992e004 4.939e006 4.392e007 9.317e008 3.74e009 4.177e010);plot(t,y,red);hold on step(G)图10 忽略高次项后的曲线匹配由上图看到，忽略高次项对系统的阶跃响应并无多大影响，但是使系统结构简化，在此基础上继续对传递函数的系数做修改。2) 调整系统的零点。由于系统零点的作用为减小峰值时间，使系统响应速度加快，并且零点越接近虚轴，这种作用越显著。但对分子常数项的修改，则会影响到系统的收敛值，较麻烦，故考虑只对分子一次项的系数进行修改。经过多次尝试，最后得到的系统传递函数如下：G=tf(6.38e006 4.182e010, 1 17.15 4570 5.992e004 4.939e006 4.392e007 9.317e008 3.74e009 4.177e010);plot(t,y,red);hold on step(G)图11 调整零点后的曲线匹配即得到最终的系统传递函数为：Transfer function: 6.38e006 s + 4.182e010-s8 + 17.15 s7 + 4570 s6 + 59920 s5 + 4.939e006 s4 + 4.392e007 s3 + 9.317e008 s2 + 3.74e009 s + 4.177e0103) 计算数据匹配率：计算公式：FIT=100(1-norm(Y-YHAT)/norm(Y-mean(Y) (in %)代码如下：yhat=step(G,t);num=norm(y-yhat,2);ym=mean(y);ymean=y;ymean(:)=ym;den=norm(y-ymean,2);fit=100*(1-num/den)得出数据匹配率为：fit =91.6637满足要求。3设计控制器3.1经典方法（超前-滞后校正）首先，我们必须得到系统的幅值裕度、频率、相位裕度、剪切频率。MATLAB控制系统工具箱中提供了margin函数，可以直接用于系统的幅值与相位裕度的求取。程序与结果如下：gm,pm,wg,wp = margin(G)gm=0.1507 pm=-113.4251 wg=9.8935 wp=14.0522由于系统的相位裕度为负，系统闭环不稳定，且要补偿的角度很大，故考虑引入滞后-超前校正器。将所给时域性能指标要求转化成频域指标，并据此设计滞后-超前校正器。MATLAB的控制系统工具箱中提供了一个控制器设计界面sisotool（）。可以得出几个控制系统设计界面，该界面允许和修改控制器的结构，允许添加零极点，调整增益，从而设计出控制器模型。以下是使用该界面设计超前校正器的过程：调用语句sisotool(G)，选择控制器结构如下图：其中H=1。选择Analysis-Other Loop Response菜单来设置系统响应的形式为阶跃响应。在当前的系统的根轨迹图上按滞后-超前校正的形式选择两个实极点与实零点，而后调整零点、极点值，并使控制器增益尽量大，以满足系统开环增益10，最终得到了较理想的响应曲线。图12 控制器设计界面图13 通过控制器得到根轨迹图14 通过控制器得到的曲线可以得到滞后-超前控制器为：图15 滞后-超前控制器传递函数3.2极点配置校正器设计如果已知对象的状态空间模型，可以采用极点配置技术，将闭环系统的全部极点均移到制定的位置，使得其动态性能达到指定值。常用的极点配置算法有Bass-Gura算法、Ackermann算法和鲁棒极点配置算法。本例采用Ackermann算法设计控制器。控制系统工具箱给出了一个acker()函数来实现该算法，该函数的调用格式为K=acker(A,B,P),A、B为系统矩阵,P为期望闭环极点的位置。由所给性能指标%=20%，ts=2sec,设计闭环主导极点为-3.6+4.8i、-3.6-4.8i,另两个极点取-40、-50。程序如下：num=31000 356800;den=1 18 4130 8300 356800;a,b,c,d=tf2ss(num,den);qc=ctrb(a,b);rank(qc); p=-3.6+4.8i -3.6-4.8i -40 -50;k=acker(a,b,p)aa=a-b*k;sy=ss(aa,b,c,d);step(sy)得到的状态反馈矩阵K=80 -1450 9340 -2848