《线性系统理论》课程——倒立摆实验报告
基本情况
实验完成了基本要求,通过pid、极点配置、根轨迹、和ldr方法调试运行一级倒立摆,设计新的pid参数,调试运行状态,逐渐使一级倒立摆稳定,完成了实验的基本要求。
在对一级倒立摆完成实验的基础上,进一步对二级倒立摆进行了分析研究。这其中的工作主要包括针对LDR方法运行demo,观察系统稳定性,快速性,调整系统参数,查看有什么问题,并且针对问题提出修改意见。在多次试验后,对系统有了进一步的了解,便开始着手二级倒立摆极点配置方法的实现问题。
这部分继续学习了极点配置的方法,通过编写m文件,计算K,仿真运行系统,查看系统图像,查看调节时间,超调量等。逐渐调试参数,使系统指标顺利达到。最后是进行试验,进一步调整系统参数。在这一个过程中,经验很重要,同时偶然因素也起到了重要的作用。所以调试一个系统真的不容易。
这一部分的内容在第六节中进行了较为详细的介绍
收获
对倒立摆的系统原理有了更深层次的了解
掌握了pid、极点配置、根轨迹、ldr方法设计系统 学会了一些调试运行系统的经验
加强了和同学之间的交流,锻炼了软件实现编程能力
改进意见
这里我有一个小小的建议,这是我在做实验的时候遇到了问题总
结。
系统参数含义还不是很清楚。在这个方面尤其是参数对应着系统的具体实际含义不明确,只能在尝试凑参数,有时出现了一个问题,不知道是哪个参数引起的,所以影响了效率,结果也不是很明显。 改进意见:共有四次实验,第一次实验安排不变但是试验后,负责人要收集问题,主要是要老师来解决的,在第二次实验前针对上一次的问题进行集体讲解一下,尤其是与物理的联系,不要仅仅是自己做实验吧,第三次和第一次相同,第四次与第二次相同。在这个完成后,如果课堂有时间,可以进行了一个小小的试验心得介绍,和大家交流心得体会。或者是老师统一解决一下这个总体过程中的问题,我觉得这样结果会更好一点。 下面是具体的详细报告
一、倒立摆系统介绍
倒立摆是机器人技术、控制理论、计算机控制等多个领域、多种技术的有机结合,其被控系统本身又是一个绝对不稳定、高阶次、多变量、强耦合的非线性系统,可以作为一个典型的控制对象对其进行研究。倒立摆系统作为控制理论研究中的一种比较理想的实验手段,为自动控制理论的教学、实验和科研构建一个良好的实验平台,以用来检验某种控制理论或方法的典型方案,促进了控制系统新理论、新思想的发展。由于控制理论的广泛应用,由此系统研究产生的方法和技术将在半导体及精密仪器加工、机器人控制技术、人工智能、导弹
拦截控制系统、航空对接控制技术、火箭发射中的垂直度控制、卫星飞行中的姿态控制和一般工业应用等方面具有广阔的利用开发前景。
倒立摆已经由原来的直线一级倒立摆扩展出很多种类,典型的有直线倒立摆环形倒立摆,平面倒立摆和复合倒立摆等,本次实验采用的是直线一级倒立摆。
倒立摆的形式和结构各异,但所有的倒立摆都具有以下的特性: 1)非线性:倒立摆是一个典型的非线性复杂系统,实际中可以通过线性化得到系统的近似模型,线性化处理后再进行控制。也可以利用非线性控制理论对其进行控制。倒立摆的非线性控制正成为一个研究的热点。2)不确定性:主要是模型误差以及机械传动间隙,各种阻力等,实际控制中一般通过减少各种误差来降低不确定性,如通过施加预紧力减少皮带或齿轮的传动误差,利用滚珠轴承减少摩擦阻力等不确定因素。3)耦合性:倒立摆的各级摆杆之间,以及和运动模块之间都有很强的耦合关系,在倒立摆的控制中一般都在平衡点附近进行解耦计算,忽略一些次要的耦合量。4)开环不稳定性:倒立摆的平衡状态只有两个,即在垂直向上的状态和垂直向下的状态,其中垂直向上为绝对不稳定的平衡点,垂直向下为稳定的平衡点。5)约束限制: 由于机构的限制,如运动模块行程限制,电机力矩限制等。为了制造方便和降低成本,倒立摆的结构尺寸和电机功率都尽量要求最小,行程限制对倒立摆的摆起影响尤为突出,容易出现小车的撞边现象。
倒立摆控制器的设计是倒立摆系统的核心内容,因为倒立摆是一个绝对不稳定的系统,为使其保持稳定并且可以承受一定的干扰,需
要给系统设计控制器,本次实验采用LQR控制。
二.直线一级倒立摆的物理模型:
系统建模可以分为两种:机理建模和实验建模。实验建模就是通 过在研究对象上加上一系列的研究者事先确定的输入信号,激励研究对象并通过传感器检测其可观测的输出,应用数学手段建立起系统的输入-输出关系。。机理建模就是在了解研究对象的运动规律基础上通过物理、化学的知识和数学手段建立起系统内部的输入-状态关系。,由于倒立摆本身是自不稳定的系统,实验建模存在一定的困难。但是忽略掉一些次要的因素后,倒立摆系统就是一个典型的运动的刚体系统,可以在惯性坐标系内应用经典力学理论建立系统的动力学方程。
下面我们采用牛顿-欧拉方法建立直线型一级倒立摆系统的数学模型:
在忽略了空气阻力和各种摩擦之后,可将直线一级倒立摆系统抽象成小车和匀质杆组成的系统,如图所示:
我们不妨做以下假设: M 小车质量 m 摆杆质量 b 小车摩擦系数
l 摆杆转动轴心到杆质心的长度 I 摆杆惯量 F 加在小车上的力 x 小车位置
φ 摆杆与垂直向上方向的夹角
θ 摆杆与垂直向下方向的夹角(考虑到摆杆初始位置为竖直向下) 图是系统中小车和摆杆的受力分析图。其中,N 和P 为小车与摆杆相互作用力的水平和垂直方向的分量。
注意:在实际倒立摆系统中检测和执行装置的正负方向已经完全确
定,因而矢量方向定义如图所示,图示方向为矢量正方向。
分析小车水平方向所受的合力,可以得到以下方程:
由摆杆水平方向的受力进行分析可以得到下面等式:
即:
把这个等式代入式(4-1)中,就得到系统的第一个运动方程:
为了推出系统的第二个运动方程,我们对摆杆垂直方向上的合力进行分析,
可以得到下面方程:
力矩平衡方程如下:
注意:此方程中力矩的方向,由于
,故等式前面有负号。
合并这两个方程,约去P 和N,得到第二个运动方程:
设θ=φ+π( φ是摆杆与垂直向上方向之间的夹角),假设φ与1(单位是弧度)相比很小,即φ<<1,则可以进行近似处理:
用u 来代表被控对象的输入力F,线性化后两个运动方程如下:
对式(3-9)进行拉普拉斯变换,得到
注意:推导传递函数时假设初始条件为0。
由于输出为角度φ,求解方程组的第一个方程,可以得到:
或
如果令
则有:
把上式代入方程组的第二个方程,得到:
整理后得到传递函数:
其中
设系统状态空间方程为:
方程组 对,
解代数方程,得到解如下:
整理后得到系统状态空间方程:
由(4-9)的第一个方程为:
对于质量均匀分布的摆杆有:
于是可以得到:
化简得到:
设
则有:
另外,也可以利用MATLAB 中tf2ss 命令对(4-20)式进行转化,求得上述状态方程。 实际系统的模型参数如下: M 小车质量 1.096 Kg m 摆杆质量 0.109 Kg b 小车摩擦系数 0 .1N/m/sec
l 摆杆转动轴心到杆质心的长度 0.2 5m I 摆杆惯量 0.0034 kg*m*m
把上述参数代入,可以得到系统的实际模型。 摆杆角度和小车位移的传递函数:
摆杆角度和小车加速度之间的传递函数为:
摆杆角度和小车所受外界作用力的传递函数:
以外界作用力作为输入的系统状态方程:
以小车加速度作为输入的系统状态方程:
三.直线一级倒立摆的线性二次型最优控制
1.线性二次最优控制LQR 基本原理及分析
线性二次最优控制LQR 基本原理为,由系统方程:
确定下列最佳控制向量的矩阵K:
u(t)=-K*x(t) 使得性能指标达到最小值:
式中 Q——正定(或正半定)厄米特或实对称阵 R——为正定厄米特或实对称阵
图 3-1 最优控制LQR 控制原理图
方程右端第二项是考虑到控制能量的损耗而引进的,矩阵Q和R确定了误差和能量损耗的相对重要性。并且假设控制向量u(t)是无约束的。
对线性系统:
根据期望性能指标选取Q 和R,利用MATLAB 命令lqr 就可以得到反馈矩阵 K 的值。
K=lqr(A,B,Q,R)
改变矩阵Q 的值,可以得到不同的响应效果,Q 的值越大(在一定的范围之内),系统抵抗干扰的能力越强,调整时间越短。但是Q
不能过大
2. LQR 控制参数调节及仿真
前面我们已经得到了直线一级倒立摆系统的比较精确的动力学模型,下面我们针对直线型一级倒立摆系统应用 LQR 法设计与调节控制器,控制摆杆保持竖直向上平衡的同时,跟踪小车的位置。 前面我们已经得到了直线一级倒立摆系统的系统状态方程:
应用线性反馈控制器,控制系统结构如下图。图中 R 是施加在小车上的阶跃输入,四个状态量x,x,φ,φ分别代表小车位移、小车速度、摆杆角度和摆杆角速度,输出y = [x,φ]’ 包括小车位置和摆杆角度。设计控制器使得当给系统施加一个阶跃输入时,摆杆会摆动,然后仍然回到垂直位置,小车可以到达新的指定位置。 假设全状态反馈可以实现(四个状态量都可测),找出确定反馈控制规律的向量K 。在 Matlab 中得到最优控制器对应的K 。Lqr 函数允许你选择两个参数——R 和Q,这两个参数用来平衡输入量和状态量的权重。最简单的情况是假设
R = 1,Q =C’ *C 。当然,也可以通过改变Q 矩阵中的非零元素来调
节控制器以得到期望的响应
其中, Q11 代表小车位置的权重,而Q33 是摆杆角度的权重,输入的权重R是1。
下面来求矩阵K,Matlab 语句为K = lqr(A,B,Q,R) 。下面在MATLAB 中编程计算:
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=1500;Q33=300; 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,Bc,Cc,Dc,U,T); plot(T,X(:,1),'-');hold on; plot(T,X(:,2),'-');hold on; plot(T,X(:,3),'.');hold on; plot(T,X(:,4),'-');
legend('cartpls','cartspd','pendang','pendspd') 令Q11= 1,Q33 =1求得
K=[-1 -1.7855 25.422 4.6849]
在 Simulink 中建立直线一级倒立摆的模型如下图所示:
“LQR Controller”为一封装好的模块,在其上单击鼠标右键,选择
“Look under
mask”打开LQR Controller 结构如下:
双击“Matrix gain K”即可输入控制参数:
点击 执行仿真,得到如下仿真结果:
LQR 控制的阶跃响应如上图所示,从图中可以看出,闭环控制系统响应的超调量很小,但稳定时间和上升时间偏大,我们可以通过增大控制量来缩短稳定时间和上升时间。
可以发现,Q矩阵中,增加Q11 使稳定时间和上升时间变短,并且使摆杆的角度变化减小。经过多次尝试,这里取Q11=1500, Q33 =300,
则K = [ -32.7298 -23.8255 81.6182 14.7098] 输入参数,运行得到响应曲线如下:
从图中可以看出,系统响应时间有明显的改善,增大Q11和Q33系统的响应还会更快,但是对于实际离散控制系统,过大的控制量会引起系统振荡。
3.直线一级倒立摆LQR控制实验 打开直线一级倒立摆LQR 实时控制模块
其中“LQR Controller”为LQR 控制器模块,“Real Control”为实时控制模块,双击“LQR Controller”模块打开LQR 控制器参数设置窗口如下:
在“LQR Controller”模块上点击鼠标右键选择“Look under mask”打开模 型如下:
双击“Real Control”模块打开实时控制模块如下图:
其中“Pendulum”模块为倒立摆系统输入输出模块,输入为小车的速度“Vel ”和“Acc ”,输出为小车的位置“Pos”和摆杆的角度
“Angle ”。
双击“Pendulum”模块打开其内部结构:
其中“Set Cart’s Acc and Vel”模块的作用是设置小车运动的速度和加速度, Get Cart’s Position”模块的作用是读取小车当前的实际位置,“Get Pend’s Angle” 的作用是读取摆杆当前的实际角度。
2) 运行程序实验运行结果如下图所示:
其中图片上半部分为小车的位置曲线,下半部分为摆杆角度的变化曲线,从图中可以看出,小车位置和摆杆角度比较稳定。控制效果很好。
在此实验中,R值固定,R=1,则只调节Q值,Q11 代表小车位置的权重,而Q33是摆杆角度的权重,若Q33增加,使得θ的变化幅度减小,而位移r的响应速度变慢;若Q11增加,使得r的跟踪速度变快,而θ的变化幅度增大。当给系统施加一个阶跃输入后,得到系统的响应结果。从响应曲线可明显看出是否满足系统所要达到的性能指标要求。通过这样反复不断的试凑,选取能够满足系统动态性能要求的Q和R。
四、直线二级倒立摆的物理模型
直线两级倒立摆由直线运动模块和两级倒立摆组件组成。 为简化系统,我们在建模时忽略了空气阻力和各种摩擦,并认为摆杆为刚体。
二级倒立摆的组成如图 6-1 所示:
图 6-1 直线两级倒立摆物理模型
倒立摆参数定义如下: M 小车质量 m1 摆杆 1 的质量 m2 摆杆2 的质量 m3 质量块的质量
l1 摆杆 1 中心到转动中心的距离 l2 摆杆2 中心到转动中心的距离 θ1 摆杆 1 与竖直方向的夹角
θ2 摆杆2 与竖直方向的夹角 F 作用在系统上的外力
利用拉格朗日方程推导运动学方程:
拉格朗日方程为:
L(q,q)=T(q,q)-V(q,q)
其中L为拉格朗日算子,q为系统的广义坐标,T为系统的动能,V为系统的势能。
其中i1,2,,n,fi为系统在第i个广义坐标上的外力,在二级倒立摆系统中,系统的广义坐标有三个广义坐标,分别为x,1,2。
首先计算系统的动能:
其中Tm,Tm1,Tm2,Tm3分别为小车的动能,摆杆1的动能,摆杆2的动能和质量块的动能。
小车的动能:
Tm1 = Tm1' +Tm2 '' 其中Tm1' ,Tm2 ' 分别为摆杆 1 的平动动能和转动动能。
Tm2 = Tm2 ' +Tm2 '' 其中Tm2 ' ,Tm2 ' 分别为摆杆2 的平动动能和转动动能。
对于系统,设以下变量: xpend1 摆杆 1 质心横坐标; yangle1 摆杆 1 质心纵坐标; xpend2 摆杆2 质心横坐标; yangle2 摆杆2 质心纵坐标; xmass 质量块质心横坐标; ymass 质量块质心纵坐标; 又有:
由于系统在θ1,θ2 广义坐标下没有外力作用,所以有:
在Mathematics中计算以上各式。
因其余各项为 0,所以这里仅列举了k12、k13、k17、k22、k23、k27 等7 项,得到结果如下:
运行的到:
五、直线二级倒立摆的线性二次型最优控制
5.1 在MATLAB Simulink 中建立直线两级倒立摆的模型:
其中“State-Space”模块为直线两级倒立摆的状态方程,双击模块打开模型:
“Controller”模块为控制器模块,在“Controller”模块上单击鼠标右键,选择 “ Look under mask”打开模型内部结构:
其中“Matrix Gain K”为反馈矩阵。
双击“Controller”模块打开其参数设置窗口: 先设置参数为“1”。
“Disturbance”模块为外界干扰模块,其作用是给系统施加一个阶跃信号,点击
“
”运行模型进行开环系统仿真。
得到运行结果如下:
从仿真结果可以看出,系统发散,为使系统稳定,需要对其添加控制器。
5.2 LQR 控制器设计及仿真
给系统添加LQR 控制器,添加控制器后的系统闭环图如下图所示
:
下面利用线性二次最优控制 LQR 方法对系统进行控制器的设计 clear;clc;
k12=86.69;k13=-21.62;k17=6.64; k22=-40.31;k23=39.45;k27=-0.088;
a=[0 0 0 1 0 0;0 0 0 0 1 0;0 0 0 0 0 1;0 0 0 0 0 0;0 k12 k13 0 0 0;0 k22 k23 0 0 0]; b=[ 0 0 0 1 k17 k27]';
c=[1 0 0 0 0 0;0 1 0 0 0 0;0 0 1 0 0 0]; d=[0; 0; 0]; q11=1;q22=1;q33=1;
q=[q11 0 0 0 0 0;0 q22 0 0 0 0;0 0 q33 0 0 0;0 0 0 0 0 0;0 0 0 0 0 0;0 0 0 0 0 0]; r=1;
k=lqr(a,b,q,r) aa=a-b*k; b=b*k(1); sys=ss(aa,b,c,d); t=0:0.01:5;
[y,t,x]=step(sys,t);
plot(t,y(:,1),'g',t,y(:,2),'r',t,y(:,3)); grid on
运行得到以下结果: LQR 控制参数为:
K=[ 1 73.818 -83.941 2.0162 4.2791 -13.036] 得到仿真结果如下:
可以看出,系统稳定时间过长,因此增加权重Q的值。 1、设Q11=300;Q22=500;Q33=500; 运行得到仿真结果: LQR 控制参数为:
K=[ 17.321 110.87 -197.57 18.468 2.7061 -32.142]
从图中可以看出,系统可以很好的稳定,在给定倒立摆干扰后,系统在 2.5 秒内可以恢复到平衡点附近。 把以上仿真参数输入 Simulink 模型中 得到运行结果
从图中可知,系统稳定性还不错。
但这未必是最好的参数。所以,下面改变LQR参数,比较结果变化。 确定最合适参数。
2. 设Q11=1000;Q22=500;Q33=500; 运行得到仿真结果: LQR 控制参数为:
k=[31.6228 116.7093 -238.1742 29.1041 1.2221 -39.3596]
可看出位置在2秒左右就可恢复到平衡点位置。而角度依然是在2.5秒内恢复到平衡位置。
把以上仿真参数输入 Simulink 模型中 得到运行结果:
3. 设Q11=1500;Q22=500;Q33=500; 运行得到仿真结果: LQR 控制参数为:
k =[44.7214 121.1834 -272.5934 38.3562 -0.0849 -45.4751]
可看出位置依然在1.5秒就可恢复到平衡点位置。而角度依然是在2.5秒内恢复到平衡位置。
把以上仿真参数输入 Simulink 模型中 得到运行结果:
4.取Q11=500;Q22=500;Q33=500; 运行得到仿真结果: LQR 控制参数为:
k =[17.321 110.87 -197.57 18.468 2.7061 -32.142]
把以上仿真参数输入 Simulink 模型中 得到运行结果:
5.设Q11=500;Q22=7000;Q33=5000; 运行得到仿真结果: LQR 控制参数为:
k =[17.32 115.24 -206.16 18.87 2.74 -33.58]
把以上仿真参数输入 Simulink 模型中 得到运行结果:
可看出位置在1.5内就可恢复到平衡点位置。而角度是在2秒内恢复到平衡位置。
通过对比,第6个参数最合适。
六、二级倒立摆的极点配置控制
系统的状态空间方程与对应的矩阵为
按极点配置步骤进行计算。
1) 检验系统可控性,由3.1.1.4系统可控性分析可以得到,系统的状态完全可控性矩阵的秩等于系统的状态维数(4),系统的输出完全可控性矩阵的秩等于系统输出向量y 的维数(2),所以系统可控。
2) 计算特征值
根据要求,并留有一定的裕量(设调整时间为2 秒),我们选取期望的闭环极点s = (i = 1,2,3,4,5,6)μi ,其中:
u110,u210,u323,u423,u510,u610
系统的主导极点是后两个,因此期望的特征方程为:
(su1)(su2)(su3)(su4)(su5)(su6)1s644s5767s46680s330200s268000s17000s0
144,2767,36680,430200,56800,67000系统的特征方程为
10,229.4,30,40,50,60
K[6a6M5a5M4a4M3a3M2a2M1a1]
极点配置仿真结果
可以看出,在给定系统干扰后,倒立摆可以在2 秒内很好的回到平衡位置,满足系统要求
系统系数
代码如下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 ]';
J=[ -10 0 0 0; 0 -10 0 0; 0 0 -2-2*sqrt(3)*i ; 0 0 0 -2+2*sqrt(3)*I;0 0 -10 0; 0 0 0 -10]; pa=poly(A);pj=poly(J);
M=[B A*B A^2*B A^3*B A^4*B A^5*B];
W=[ pa(6) pa(5) pa(4) pa(3) pa(2) 1; pa(5) pa(4) pa(3) pa(2) 1 0; pa(4) pa(3) pa(2) 1 0 0; pa(3) pa(2) 1 0 0 0;pa(2) 1 0 0 0 0 0]; T=M*W;
K=[ pj(7)-pa(7) pj(6)-pa(6) pj(5)-pa(5) pj(4)-pa(4) pj(3)-pa(3) pj(2)-pa(2)]*inv(T) 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; plot(T,X(:,2),'-.'); hold on; plot(T,X(:,3),'.'); hold on; plot(T,X(:,4),'-')
legend('CartPos','CartSpd','PendAng','PendSpd')
因篇幅问题不能全部显示,请点此查看更多更全内容