倒立摆实验报告

倒立摆实验报告 实验人: 2011 年 12 月 29 日 目录 一实验概述 倒立摆简介 3 倒 立 摆 控 制 的 目 标 与 控 制 器 设 计 方 法 3 实验装置简介 3 二实验内容 实验对象建模 4 控制器的设计、仿真与实验结果 频率响应 控制实验 5 直线一级倒立摆 PID 控制实验 13 直线一级倒立摆 LQR 控制实验 17 三小结 21 一实验概述 1.倒立摆简介 倒 立 摆 控 制 系 统 是 一 个 复 杂 的 、 不 稳 定 的 、 非 线 性 系 统 ， 是 进 行 控 制 理 论 教 学 及 开 展 各 种 控 制 实 验 的 理 想 实 验 平 台 。 对 倒 立 摆 系 统 的 研 究 能 有 效 的 反 映 控 制 中 的 许 多 典 型 问 题 ： 如 非 线 性 问 题 、 鲁 棒 性 问 题 、 镇 定 问 题 、 随 动 问 题 以 及 跟 踪 问 题 等 。 通 过 对 倒 立 摆 的 控 制 ， 用 来 检 验 新 的 控 制 方 法 是 否 有 较 强 的 处 理 非 线 性 和 不 稳 定 性 问 题 的 能 力 。 同 时 ， 其 控 制 方 法 在 军 工 、 航 天 、 机 器 人 和 一 般 工 业 过 程 领 域 中 都 有 着 广 泛 的 用 途 ， 如 机 器 人 行 走 过 程 中 的 平 衡 控 制 、 火 箭 发 射 中 的 垂 直 度 控 制 和 卫 星 飞 行 中 的 姿 态 控 制 等 。 2.倒 立 摆 控 制 的 目 标 与 控 制 器 设 计 方 法 。 倒 立 摆 的 控 制 问 题 就 是 使 摆 杆 尽 快 地 达 到 一 个 平 衡 位 置 ， 并 且 使 之 没 有 大 的 振 荡 和 过 大 的 角 度 和 速 度 。 当 摆 杆 到 达 期 望 的 位 置 后 ， 系 统 能 克 服 随 机 扰 动 而 保 持 稳 定 的 位 置 。 本 实 验 的 控 制 对 象 是 一 级 倒 立 摆 ,控 制 目 标 是 实 现 起 摆 后 摆 杆 稳 定 于 倒 立 状 态 ,有 一 定 的 抗 干 扰 能 力 。 倒 立 摆 系 统 的 输 入 为 小 车 的 位 移 （ 即 位 置 ） 和 摆 杆 的 倾 斜 角 度 期 望 值 ， 计 算 机 在 每 一 个 采 样 周 期 中 采 集 来 自 传 感 器 的 小 车 与 摆 杆 的 实 际 位 置 信 号 ， 与 期 望 值 进 行 比 较 后 ， 通 过 控 制 器 处 理 得 到 控 制 量 ， 再 经 数 模 转 换 驱 动 直 流 电 机 实 现 倒 立 摆 的 实 时 控 制 。 直 流 电 机 通 过 皮 带 带 动 小 车 在 固 定 的 轨 道 上 运 动 ， 摆 杆 的 一 端 安 装 在 小 车 上 ， 能 以 此 点 为 轴 心 使 摆 杆 能 在 垂 直 的 平 面 上 自 由 地 摆 动 。 作 用 力 u 平 行 于 铁 轨 的 方 向 作 用 于 小 车 ， 使 杆 绕 小 车 上 的 轴 在 竖 直 平 面 内 旋 转 ， 小 车 沿 着 水 平 铁 轨 运 动 。 当 没 有 作 用 力 时 ， 摆 杆 处 于 垂 直 的 稳 定 的 平 衡 位 置 （ 竖 直 向 下 ） 。 为 了 使 杆 子 摆 动 或 者 达 到 竖 直 向 上 的 稳 定 ， 需 要 给 小 车 一 个 控 制 力 ， 使 其 在 轨 道 上 被 往 前 或 朝 后 拉 动 。 其 中 控 制 器 的 设 计 方 法 有 PID 控制、根轨迹以及频率响应法、状态空间法、最优控制理论、模糊控 制理论、神经网络控制、拟人智能控制、鲁棒控制方法、自适应控制，以及这些控制理论的相互结合组成 更加强大的控制算法。本实验采用频率响应法、PID 控制和 LQR 控制。 3.实验装置简介。 固高科技公司倒立摆实验装置(旋转编码器作为角位移传感器对输出采样，交流伺服电机驱动导轨上的 小车作为执行器) PC机; 基于Matlab simulink平台的固高实验控制软件。 二实验内容 1.实验对象建模。 PNFbxMmg 忽略掉一些次要因素后，倒立摆系统就是一个典型的运动刚体系统，可以在惯性坐标系内应用经典力 学理论建立系统的动力学方程，对倒立摆系统进行机理建模。 在忽略了空气阻力, 各种摩擦之后, 可将直线一级倒立摆系统抽象成小车和匀质杆组成的系统F的作用下， 小车及摆均产生加速运动，根据牛顿第二定律，在水平直线运动方向的惯性力与F平衡和绕摆轴旋转运动 的惯性力矩应与重力力矩平衡，我们可以得到下面运动方程： （1）cossin)( in2 2xmlgllIxbM 设 ， ，进行近似处理： ， ， 。用 u 来代表被控对11cossin0)(2dt 象的输入力 F，对式（5-1）进行线性化后得到下面两个运动方程： （2）xmlglI ubx)(2 整理后得系统的状态方程如下： （3）CXYBA 其中： ， ， 22222010()()100()()ImlbglMIMmAllIIl 220()()ImlMBlI10C1234xX 带入固高系统固件的相关参数，可得到角度作为输出与加速度作为输入的传递函数： 2.控制器的设计、仿真与实验结果。 （ 1）直线一级倒立摆 频率响应 控制实验 系统对正弦输入信号的响应，称为频率响应。在频率响应方法中，我们在一定范围内改变输入信号的频率， 研究其产生的响应。频率响应可以采用以下三种比较方便的方法进行分析，一种为伯德图或对数坐标图，伯德 图采用两幅分离的图来表示，一幅表示幅值和频率的关系，一幅表示相角和频率的关系；一种是极坐标图, 也 常称为奈奎斯特图，奈奎斯特稳定判据使我们有可能根据系统的开环频率响应特性信息，研究线性闭环系统的 绝对稳定性和相对稳定性。 频率响应分析 由直线一级倒立摆的物理模型，实际系统的开环传递函数： 其中输入为小车的加速度V (s) ，输出为摆杆的角度 (s) 。 在 MATLAB 下绘制系统的 Bode 图和奈奎斯特图，程序如下： m=0.109; L=0.25; I=0.0034; g=9.8; num=m*L; den=I+m*L*L 0 -m*g*L; z=roots(num); p=roots(den); subplot(2,1,1); bode(num,den); subplot(2,1,2); nyquist(num,den)； 程序运行结果如下： -100 -50 0 Ma gni tud e (d B) 100 101 102 -181 -180 -179 Pha se (de g) Bode Diagram Frequency (rad/sec) -1 -0.9 -0.8 -0.7 -0.6 -0.5 -0.4 -0.3 -0.2 -0.1 0-1 -0.5 0 0.5 1 Nyquist Diagram Real Axis Ima gin ary Ax is 可以得到，系统没有零点，但存在两个极点，其中一个极点位于右半 s 平面，根据奈奎斯特稳定判据， 闭环系统稳定的充分必要条件是：当 从 到 + 变化时，开环传递函数 G( j) 沿逆时针方向包围- 1 点 p 圈，其中 p 为开环传递函数在右半 S 平面内的极点数。对于直线一级倒立摆，由图 3-21 我们可 以看出，开环传递函数在 S 右半平面有一个极点，因此 G( j) 需要沿逆时针方向包围-1 点一圈。可 以看出，系统的奈奎斯特图并没有逆时针绕-1 点一圈，因此系统不稳定，需要设计控制器来镇定系统。 选取系统静态误差为10，相角裕量50，增益裕量不小于10分贝。Matlab仿真程序如下： m=0.109; L=0.25; I=0.0034; g=9.8; %输入实际物理参数 Kv=10; Ko=Kv*g; %确定静差和开环增益 num=Ko*m*L; den=I+m*L*L 0 -m*g*L; z=roots(num); p=roots(den); %根据已知传递函数计算零极点 subplot(2,1,1) bode(num,den) subplot(2,1,2) nyquist(num,den) %画伯德图和奈奎斯特曲线 程序运行结果如下： -50 0 50 Ma gni tud e (d B) 100 101 102 -181 -180 -179 Pha se (de g) Bode Diagram Frequency (rad/sec) -10 -9 -8 -7 -6 -5 -4 -3 -2 -1 0-1 -0.5 0 0.5 1 Nyquist Diagram Real Axis Ima gin ary Ax is 可以看出，系统相角裕量为0，由设计要求，取最大超前相角为55，用如下Matlab程序得出超前控 制器参数并绘制相应伯德图及奈奎斯特曲线： m=0.109; L=0.25; I=0.0034; g=9.8; %输入实际物理参数 Kv=10; Ko=Kv*g; %确定静差和开环增益 a=(1+sin(0.96)/(1-sin(0.96); %求校正参数a wc1=(Ko*m*L*a0.5)-m*g*L)/(I+m*L2)0.5; %求校正后幅值穿越频率 p1=-(a0.5)*wc1; z1=p1/a; num=Ko*m*L; den=I+m*L*L 0 -m*g*L; z=roots(num); p=roots(den); za=z;z1; pa=p;p1; %计算校正后系统零极点 Kc1=Ko*a; %校正后开环增益 sys=zpk(za,pa,Kc1); sysc=sys/(1+sys); subplot(2,1,1) bode(sys) subplot(2,1,2) nyquist(sys) %输出伯德图和奈氏曲线 运行结果如下： -100 0 100 Ma gni tud e ( dB ) 10-1 100 101 102 103 104 -180 -135 -90 Ph ase (d eg) Bode Diagram Frequency (rad/sec) -4 -3.5 -3 -2.5 -2 -1.5 -1 -0.5 0-1 -0.5 0 0.5 1 Nyquist Diagram Real Axis Ima gin ary Ax is 加上命令行： figure sysc=sys/(1+sys); t=0:0.005:5; impulse(sysc,t) 可得其冲击响应如下 0 0.5 1 1.5 2 2.5 3 3.5 4 4.5 5-2 0 2 4 6 8 10 12 Impulse Response Time (sec) Am plitu de 可以看出，系统在遇到干扰后，在 1 秒内可以达到新的平衡，但是超调量比较大。 为提高稳态精度，缩短响应时间，在已完成的超前校正基础上进行滞后校正，设计15为新的穿越频率. Matlab仿真程序如下： m=0.109; L=0.25; I=0.0034; g=9.8; %输入实际物理参数 Kv=10; Ko=Kv*g; %确定静差和超前校正开环增益 a=(1+sin(0.96)/(1-sin(0.96); wc1=(Ko*m*L*a0.5)-m*g*L)/(I+m*L2)0.5; wc2=10; %新穿越频率 p1=-(a0.5)*wc1; z1=p1/a; b=(wc22+p12)0.5/(wc22+z12)0.5*(I+m*L2)*wc22+m*g*L)/(m*L*Ko); z2=-wc2/10; p2=b*z2; num=Ko*m*L; den=I+m*L*L 0 -m*g*L; z=roots(num); p=roots(den); za=z;z1;z2; pa=p;p1;p2; %计算校正后系统零极点 Kc2=Ko*a*b; %确定新的开环增益 sys=zpk(za,pa,Kc2); sysc=sys/(1+sys); subplot(2,1,1) bode(sys) subplot(2,1,2) nyquist(sys) %输出伯德图和奈氏曲线 运行结果如下： -100 -50 0 50 Ma gni tud e (d B) 10-2 10-1 100 101 102 103 104 -180 -150 -120 Pha se (de g) Bode Diagram Frequency (rad/sec) -4 -3.5 -3 -2.5 -2 -1.5 -1 -0.5 0-1 -0.5 0 0.5 1 Nyquist Diagram Real Axis Ima gin ary Ax is 显然系统稳定. 实验室调试: 按理论设计得到的控制器零点为:Za=-8.9353;-0.02 ; 极点为Pa=-89.902;-0.020149 开环增益Kc2=993.35 实验现象:倒立摆无法保持平衡,小车会向一个方向运动直到碰到挡板. 考虑存在静差. 将Kv改为100，并适当降低穿越频率wc2=5，重新由Matlab算得Za=-28.667;-0.5 ; 极点为Pa=-288.43;-0.096954 开环增益Kc2=1912 仿真结果如下: -100 0 100 Ma gni tud e (d B) 10-3 10-2 10-1 100 101 102 103 104 -270 -180 -90 Pha se (de g) Bode Diagram Frequency (rad/sec) -40 -35 -30 -25 -20 -15 -10 -5 0-20 -10 0 10 20 Nyquist Diagram Real Axis Ima gin ary Ax is 单位阶跃响应： 0 0.5 1 1.5 2 2.5 3 3.5 4 4.5 50 0.2 0.4 0.6 0.8 1 1.2 1.4 1.6 Step Response Time (sec) Amp litud e 连接倒立摆装置输入调整后参数： 立摆成功,达到控制目的.监控软件输出如下: （ 2）直线一级倒立摆 PID 控制实验 PID 理论控制分析 经典控制理论的研究对象主要是单输入单输出的系统，控制器设计时一般需要有关被控对象的较精确 模型。PID 控制器因其结构简单，容易调节，且不需要对系统建立精确的模型，在控制上应用较广。 首先，对于倒立摆系统输出量为摆杆的角度，它的平衡位置为垂直向上的情 况。系统控制结构框图如下： 图中 KD(s) 是控制器传递函数， G( s) 是被控对象传递函数。 考虑到输入 r(s) = 0 ，结构图可以很容易的变换成： 该系统的输出为： 其中： num 被控对象传递函数的分子项 den 被控对象传递函数的分母项 numPID PID 控制器传递函数的分子项 denPID PID 控制器传递函数的分母项 通过分析上式就可以得到系统的各项性能。 由(3-13) 可以得到摆杆角度和小车加速度的传递函数： PID 控制器的传递函数为： 需仔细调节 PID 控制器的参数，以得到满意的控制效果。 小车位置输出为： 通过对控制量 v 双重积分即可以得到小车位置。 PID 实验参数设定及仿真 由实际系统的物理模型： 在 Simulink 中建立如图所示的直线一级倒立摆模型： 参数设定过程： a.先设置 PID 控制器为 P 控制器，令 Kp = 1,Kii= 0,Kd= 0 ，得到以下仿真结果： Kp = 1,Kii= 0,Kd= 0 由上图可见，控制曲线不收敛，必须增大 Kp。 b. 令 Kp = 50,Kii= 0,Kd= 0 Kp = 50,Kii= 0,Kd= 0 从图中可以看出，闭环控制系统持续振荡，周期约为 0.6s。为消除系统的振 荡，增加微分控制参数 Kd。 c. 令 Kp = 50,Kii= 0,Kd= 8 Kp = 50,Kii= 0,Kd= 8 从上图可以看出，系统在 1.2秒后达到平衡，但是存在一定的稳态误差。 为消除稳态误差，我们增加积分参数 Ki d. 令 Kp = 50,Ki i= 15,Kd= 8 Kp = 50,Kii= 15,Kd= 8 从上图可以看出，系统最终达到稳态，且误差近似为 0，但稳定时间较长，约 10 秒，因此增大积分参 数 Ki e. 令 Kp = 50,Ki i= 27,Kd= 8 Kp = 50,Kii= 27,Kd= 8 由上图可见，稳定时间有所缩短，在 6.2s 处达到稳定，且误差约为 0，各方面参数满足要求。 小车位置输出曲线如下： 可以看出，由于 PID 控制器为单输入单输出系统，所以只能控制摆杆的角度，并不能控制小车的位置， 所以小车会往一个方向运动。 PID 控制参数验证 将仿真所得的参数:Kp=50,Ki=27,Kd=8,输入到PID控制器中，编译，连接，运行程序。一切正常时，缓 慢提起倒立摆的摆杆到竖直向上的位置，在程序进入自动控制后松开，当小车运动到正负限位的位置时， 用工具挡一下摆杆，使小车反向运动。 试验结果如下图： 从图中可以看出，倒立摆可以实现较好的稳定性，摆杆的角度在 3.142（弧度左右。同仿真结果，PID 控制器并不能对小车的位置进行控制，小车会沿滑杆有稍微的移动。 （ 3）直线一级倒立摆 LQR 控制实验 线性二次最优控制 LQR 基本原理及分析 线性二次最优控制 LQR 基本原理为，由系统方程: 确定下列最佳控制向量的矩阵 K： 使得性能指标达到最小值: 式中 Q正定(或正半定)厄米特或实对称阵 R为正定厄米特或实对称阵 方程右端第二项是考虑到控制能量的损耗而引进的，矩阵 Q 和 R 确定了误差和能量损耗的相对重要性。 并且假设控制向量 u(t)是无约束的。 对线性系统： 根据期望性能指标选取 Q 和 R，利用 MATLAB 命令 lqr 就可以得到反馈阵K 的值： 改变矩阵 Q 的值，可以得到不同的响应效果， Q 的值越大 (在一定的范围之内)，系统抵抗干扰的能力越 强，调整时间越短。 LQR 控制参数调节及仿真 前面我们已经得到了直线一级倒立摆系统的比较精确的动力学模型，下面我们针对直线型以及倒立摆系统 应用LQR法设计与调节控制器，控制摆杆保持竖直向上平衡的同时，跟踪小车的位置。 前面我们已经得到了直线一级倒立摆系统的系统状态方程： 假设全状态反馈可以实现，四个状态量都可以测，找出确定反馈控制规律的向量K。在Matlab中得到最优 控制器对应的K。LQR函数允许选择两个参数：R，Q，这两个参数用来平衡输入量和状态量的权重。最简单的情 况是假设R=1， 。当然也可以通过改变Q矩阵中非零元素来调节控制器以得到期望的响应。 其中，Q11代表小车位置的权重，Q33代表摆杆角度的权重，输入的权重R是1。 下面是求矩阵K的Matlab程序 clear; A= 0 1 0 0; 0 0 0 0; 0 0 0 1; 0 0 29.4 0; B= 0 1 0 3; C= 1 0 0 0; 0 0 1 0; D= 0 0 ; Q11=1; Q33=1; Q=Q11 0 0 0; 0 0 0 0; 0 0 Q33 0; 0 0 0 0; R = 1; K = lqr(A,B,Q,R) Ac = (A-B*K); Bc = B; Cc = C; Dc = D; T=0:0.005:5; U=0.2*ones(size(T); Cn=1 0 0 0; Nbar=rscale(A,B,Cn,0,K); Bcn=Nbar*B; Y,X=lsim(Ac,Bcn,Cc,Dc,U,T); plot(T,X(:,1),-);hold on