基于MATLAB（矩阵实验室)的倒立摆控制系统仿真 摘要 自动控制原理（包括经典部分和现代部分）是电气信息工程学院学生的一门必修专业基础课，课程中的一些概念相对比较抽象，如系统的稳定性、可控性、收敛速度和抗干扰能力等。倒立摆系统是一个典型的非线性、强耦合、多变量和不稳定系统，作为控制系统的被控对象，它是一个理想的教学实验设备，许多抽象的控制概念都可以通过倒立摆直观地表现出来。本文以一级倒立摆为被控对象，用经典控制理论设计控制器（PID控制器）的设计方法和用现代控制理论设计控制器（极点配置）的设计方法，通过MATLAB仿真软件的方法来实现。 关键词：一级倒立摆PID控制器极点配置
Inverted pendulum controlling system simulation based on the MATLAB ABSTRACT Automatic control theory (including classical parts and modern parts) is a compulsory specialized fundamental course of the students majored in electrical engineering. Some of the curriculum concept is relatively abstract, such as the stability, controllability, convergence rate and the anti-interference ability of system. Inverted pendulum system is a typical nonlinear, strong coupling, multivariable and unstable system. It is an ideal teaching experimental equipment as a controlled object, by which many abstract control concepts can be came out directly. This paper chose first-order inverted pendulum as the controlled object. First, the PID controller was designed with classical control theory. Then pole-assignment method was discussed with modern control theory. At last, the effectness of the two methods was verified by MATLAB simulation software. KEY WORDS: First-order inverted pendulum PID controller pole-assignment
MATLAB中如何直接曲线拟合，而不使用cftool的GUI 界面 我们知道在MATLAB中有个很方便的曲线拟合工具：cftool 最基本的使用方法如下，假设我们需要拟合的点集存放在两个向量X和Y中，分别储存着各离散点的横坐标和纵坐标，则在MATLAB中直接键入命令 cftool(X,Y) 就会弹出Curve Fitting Tool的GUI界面，点击界面上的fitting即可开始曲线拟合。 MATLAB提供了各种曲线拟合方法，例如：Exponential, Fourier, Gaussing, Interpolant, Polynomial, Power, Rational, Smoothing Spline, Sum of Functions, Weibull等，当然，也可以使用 Custom Equations. cftool不仅可以绘制拟合后的曲线、给出拟合参数，还能给出拟合好坏的评价 参数(Goodness of fit)如SSE, R-square, RMSE等数据，非常好用。但是如果我们已经确定了拟合的方法，只需要对数据进行计算，那么这种GUI的操作方式就不太适合了，比如在m文件中就不方便直接调用cftool。 MATLAB已经给出了解决办法，可以在cftool中根据情况生成特定的m文件，让我们直接进行特定的曲线拟合并给出参数。具体方法在帮助文件的如下文档中" \ Curve Fitting Toolbox \ Generating M-files From Curve Fitting Tool " ，以下简单举例说明： 以双色球从第125期到第145期蓝球为Y值： Y=[12 15 4 1 7 11 5 7 1 6 16 1 1 14 2 12 9 13 10 12 11]; X=1:1:21; cftool(X,Y); 点击Fitting选择最常用的多项式拟合(Polynomial)，选择3次多项式拟合(cubic)，然后就会出现如下拟合图形： 然后在Curve Fitting Tool窗口中点击 " \ File \ Generate M-file " 即可生成能直接曲线拟合的m函数文件，其中使用的拟合方法就是刚才使用的三次多项式拟合，文件中这条语句证明了这一点： ft_ = fittype('poly3'); 保存该m文件（默认叫做createFit.m），调用方法和通常的m文件一样，使用不同的X和Y值就能拟合出不同的曲线。但是，这种调用方法只能看到一个拟合出的图形窗口，拟合参数以及Goodness of fit参数都看不到了，因此需要在刚才的m文件中稍作修改。 找到这句话： cf_ = fit(X(ok_),Y(ok_),ft_); 修改为： [cf_,gof] = fit(X(ok_),Y(ok_),ft_); 然后将函数声明 function createFit(X,Y) 修改为 function [cf_,gof] = createFit(X,Y) ，这样我们再调用试试看： Y=[12 15 4 1 7 11 5 7 1 6 16 1 1 14 2 12 9 13 10 12 11]; X=1:1:21;
第一章绪论 1.1倒立摆系统的简介 1.1.1倒立摆系统的研究背景及意义 倒立摆系统的最初分析研究开始于二十世纪五十年代，是一个比较复杂的不稳定、多变量、带有非线性和强耦合特性的高阶机械系统，它的稳定控制是控制理论应用的一个典型范例[1]。倒立摆系统存在严重的不确定性，一方面是系统的参数的不确定性，一方面是系统的受到不确定因素的干扰。通过对它的研究不仅可以解决控制中的理论问题，还将控制理论涉及的相关主要学科：机械、力学、数学、电学和计算机等综合应用。在多种控制理论与方法的研究和应用中，特别是在工程中，存在一种可行性的实验问题，将其理论和方法得到有效的验证，倒立摆系统可以此提供一个从控制理论通过实践的桥梁。近些年来，国内外不少专家、学者一直将它视为典型的研究对象，提出了很多控制方案，对倒立摆系统的稳定性和镇定问题进行了大量研究，都在试图寻找不同的控制方法实现对倒立摆的控制，以便检查或说明该方法的严重非线性和绝对不稳定系统的控制能力，其控制方法在军工、航天、机械人领域和一般工业过程中都有着广泛的用途，如精密仪器的加工、机器人行走过程中的平衡控制、火箭发射中的垂直度控制、导弹拦截控制、航空对接控制、卫星飞行中的姿态控制等方面均涉及到倒置问题。因此，从控制这个角度上讲，对倒立摆的研究在理论和方法论上均有着深远意义。倒立摆系统是一个典型的自不稳定系统，其中摆作为一个典型的振动和运动问题，可以抽象为许多问题来研究。随着非线性科学的发展，以前的采用线性化方法来描述非线性的性质，固然无可非议，但这种方法是很有局限性，非线性的一些本质特征往往不是用线性的方法所能体
现的。非线性是造成混乱、无序或混沌的核心因素，造成混乱、无序或混沌并不意味着需要复杂的原因，简单的非线性就会产生非常的混乱、无序或混沌。在倒立摆系统中含有极其丰富和复杂的动力学行为，如分叉、分形和混沌动力学，这方面的问题也值得去探讨和研究。 无论哪种类型的倒立摆系统都具有如下特性[2]： (1)非线性倒立摆是一个典型的非线性复杂系统。实际中可以通过线性化得到系统的近似模型，线性化处理后再进行控制，也可以利用非线性控制理论对其进行控制，倒立摆的非线性控制正成为一个研究的热点。 (2)不确定性主要是指建立系统数学模型时的参数误差、量测噪声以及机械传动过程中的减速齿轮间隙等非线性因素所导致的难以量化的部分。 (3)欠冗余性一般的，倒立摆控制系统采用单电机驱动，因而它与冗余机构，比如说冗余机器人有较大的不同。之所以采用欠冗余的设计是要在不失系统可靠性的前提下节约经济成本或者节约有效的空间。研究者常常是希望通过对倒立摆控制系统的研究获得性能较为突出的新型控制器设计方法，并验证其有效性及控制性能。 (4)耦合特性倒立摆摆杆和小车之间，以及多级倒立摆系统的上下摆杆之间都是强耦合的。这既是可以采用单电机驱动倒立摆控制系统的原因，也是使得控制系统的设计、控制器参数调节变得复杂的原因。 (5)开环不稳定性倒立摆系统有两个平衡状态：垂直向下和垂直向上。垂直向下的状态是系统稳定的平衡点(考虑摩擦力的影响)，而垂直向上的状态是系统不稳定的平衡点，开环时微小的扰动都会使系统离开垂直向上的状态而进入到垂直向下的状态中。 (6)约束限制由于实际机构的限制，如运动模块行程限制，电机力矩限制等。为制造方便和降低成本，倒立摆的结构尺寸和电机功率都尽量要求最小，行程限制对于倒立
第1 页共11 页 倒立摆系统的建模及Matlab仿真 1.系统的物理模型 考虑如图(1)所示的倒立摆系统。图中，倒立摆安装在一个小车上。这里仅考虑倒立摆在图面内运动的二维问题。 图(1)倒立摆系统 假定倒立摆系统的参数如下。 摆杆的质量：m=0.1g l=1m小车的质量：摆杆的长度：2重力加速度：g=9.8m/M=1kg s摆杆的质量在摆杆的中心。 设计一个控制系统，使得当给定任意初始条件(由干扰引起)时，最大超调量?≤10%，调节时间ts ≤4s ，通过小车的水平运动使倒立摆保持在垂直位置。 2.系统的数学模型 2.1建立倒置摆的运动方程并将其线性化。 为简化问题，在数学模型中首先假设：1)摆杆为刚体；2）忽略摆杆与支点之间的摩擦；3）忽略小车与接触面间的摩擦。 ?）,在u设小车瞬时位置为z,摆心瞬时位置为（作用下，小车及摆均产生加速远 动，sin?lz根据牛顿第二定律，在水平直线远动方向的惯性力应与u平衡，于是有 22dzd?)?sinu?M?m(zl22dtdt???2????z(M?mml?)cos?mlusin? 即：??①
绕摆轴转动的惯性力矩与重力矩平衡，因而有． 第2 页共11 页 2??d??? sin??lcosm(z?lsinmgl)??2dt?????22???????即： nis?l?ocgcosincoszs?ls??② 以上两个方程都是非线性方程，为求得解析解，需作线性化处理。由于控制的目的是保持倒立摆直?2?????且可忽略则，立，在试驾合适的外力条件下，假定θ很小，接近于零时合理的，1sincos??,项。于是有 ???M?zm?u?ml??)(③ ????g?z?l??④联立求解可得1mg?u?z????MM 1)?m(M????u??MlMl 列写系统的状态空间表达式。2.2??T xx,x,x,，选取系统变量则 xx,x,xx?,42134123xx??211mgux???x?32MM x?x?431)(M?mu?x?x? 34MlMl 即00100????z??1mg??????000?z?????d MM??Bu?Ax?xux????????00001???dt????1gm?(M)????000??????? MlMl??????Cx?0?y?xx1001代入数据计算得到：0100????000?1??????T0D,?0??1BA?,?001,C100??1000??00011?? 11 页3 页共第 3.设计控制器3.1判断系统的能控性和稳定性 1100????0011????23BBAABAB?Q?故被控对象完全可控， rank()=4,Q kk??11?0?10??011?10???22???11?。出现大于零的特征值，故被，，0 解得特征值为 0由特征方程0??11I?A?)(控对象不稳定3.2确定希望的极点， 另一对为远极点，认为系统性能主要由主导，选其中一对为主导极点和希望的极点n=4ss21极点决定，远极点只有微小影响。根据二阶系统的关系式，先确定主导极点???42??1????10.?e??t1.67?有，闭环可得；取误差带，于是取，则6.?059?0.02.?0? pns??n2????1?js??=-10.8j,远极点选择使它和原点的距离大于主导极点与原点 距离主导极点为?n,21s??15倍，取的54,33.3采用状态反馈方法使系统稳定并配置极点 ??kkkk?k；状态反馈系统的状态方程，馈状态反的控制规律为为kxu??3102?，其
MATLAB机械工程 最小二乘法曲线拟合的应用实例 班级: 姓名: 学号: 指导教师:
一，实验目的 通过Matlab上机编程，掌握利用Matlab软件进行数据拟合分析及数据可视化方法 二，实验内容 1.有一组风机叶片的耐磨实验数据，如下表所示，其中X为使用时间，单位为小时h，Y为磨失质量，单位为克g。要求： 对该数据进行合理的最小二乘法数据拟合得下列数据。 x=[10000 11000 12000 13000 14000 15000 16000 17000 18000 19000 2 0000 21000 22000 23000]; y=[24.0 26.5 29.8 32.4 34.7 37.7 41.1 42.8 44.6 47.3 65.8 87.5 137.8 174. 2] 三，程序如下 X=10000:1000:23000; Y=[24.0,26.5,29.8,32.4,34.7,37.7,41.1,42.8,44.6,47.3,65.8,87.5,137.8,17 4.2] dy=1.5; %拟合数据y的步长for n=1:6 [a,S]=polyfit(x,y,n); A{n}=a;
da=dy*sqrt(diag(inv(S.R′*S.R))); Da{n}=da′; freedom(n)=S.df; [ye,delta]=polyval(a,x,S); YE{n}=ye; D{n}=delta; chi2(n)=sum((y-ye).^2)/dy/dy; end Q=1-chi2cdf(chi2,freedom); %判断拟合良好度 clf,shg subplot(1,2,1),plot(1：6,abs(chi2-freedom),‘b’) xlabel(‘阶次’)，title(‘chi2与自由度’) subplot(1,2,2),plot(1：6,Q,‘r’,1：6,ones(1,6)*0.5) xlabel(‘阶次’),title(‘Q与0.5线’) nod=input(‘根据图形选择适当的阶次（请输入数值）’)； elf,shg, plot(x,y,‘kx’)；xlabel(‘x’),ylabel(‘y’)； axis([8000,23000,20.0,174.2])；hold on errorbar(x,YE{nod},D{nod},‘r’)；hold off title(‘较适当阶次的拟合’) text(10000,150.0,[‘chi2=’num2str(chi2(nod))‘~’int2str(freedom(nod))])
1、一元线性拟合 求HNO 3的正常沸点温度T b 及摩尔汽化热。 程序如下： >> t=[0 20 40 50 70 80 90 100]; >> t=t+273.15; >> p=[1919.52 6385.07 17728.9 27726.4 62251.1 89311 124902.1 170890.6] p = 1.0e+005 * 0.0192 0.0639 0.1773 0.2773 0.6225 0.8931 1.2490 1.7089 >> subplot 121 >> plot(t,p,'o',t,p) >> t1=1./t;p2=log(p); >> pp=polyfit(t1,p2,1) pp = 1.0e+003 * -4.5691 0.0243 >> subplot 122 >> plot(t1,p2,'o',t1,p2) >> gtext('p/pa'),gtext('T/K'),GTEXT('lnP/Pa'),gtext('T^-^1/K') 由克拉贝龙-克劳修斯方程式，~ ln v H P C RT ?=-+ 作1 ln ~P T -得一直线：3 1 ln 4.5691024.30P T -=-?+ 斜率为：~ 3 4.56910v H R ?-?=-
所以摩尔汽化热为：~ 314.569108.31437.99()v H kJ mol -?=??=? 并根据拟合方程，求得一大气压时 1 32.8010T --=? 则正常沸点为：357b T K = 2、多元线性拟合： 某气体混合物由四种气体组成，在常压或低压下其粘度η与各组分摩尔分数x 1,x 2,x 3,x 4之间有如下线性关系：011223344b b x b x b x b x η=++++ 试根据下表所列实验数据用最小二乘法确定上式中的各个系数，并计算其复相关系数。 Matlab 程序如下： >> a=[1.0 0.402 0.153 0.058 0.387;1.0 0.503 0.301 0.183 0.013; 1.0 0.306 0.109 0.224 0.361; 1.0 0.296 0.365 0.009 0.330; 1.0 0.309 0.405 0.109 0.177; 1.0 0.055 0.153 0.506 0.289] a = 1.0000 0.4020 0.1530 0.0580 0.3870 1.0000 0.5030 0.3010 0.1830 0.0130 1.0000 0.3060 0.1090 0.2240 0.3610 1.0000 0.2960 0.3650 0.0090 0.3300 1.0000 0.3090 0.4050 0.1090 0.1770 1.0000 0.0550 0.1530 0.5060 0.2890 >> y=[0.00625 0.00826 0.01182 0.01944 0.02372 0.03243]' y = 0.0063 0.0083 0.0118 0.0194 0.0237 0.0324 >> b=a.'*a
MATLAB软件提供了基本的曲线拟合函数的命令。 曲线拟合就是计算出两组数据之间的一种函数关系，由此可描绘其变化曲线及估计非采集数据对应的变量信息。 1.线性拟合函数：regress() 调用格式： b = regress(y,X) [b,bint,r,rint,stats] = regress(y,X) [b,bint,r,rint,stats] = regress(y,X,alpha) 说明：b=[ε; β]，regress(y,X)返回X与y的最小二乘拟合的参数值β、ε，y=ε+βX。β是p′1的参数向量；ε是服从标准正态分布的随机干扰的n′1的向量；y为n′1的向量；X为n′p矩阵。 bint返回β的95%的置信区间。 r中为形状残差，rint中返回每一个残差的95%置信区间。Stats向量包含R2统计量、回归的F值和p值。 例： x=[ones(10,1) (1:10)']; y=x*[10;1]+normrnd(0,0.1,10,1); [b,bint]=regress(y,x,0.05) 结果得回归方程为：y=9.9213+1.0143x 2.多项式曲线拟合函数：polyfit() 调用格式： p = polyfit(x,y,n) [p,s] = polyfit(x,y,n) 说明：n：多项式的最高阶数； x，y：将要拟合的数据，用数组的方式输入； p：为输出参数，即拟合多项式的系数； 多项式在x处的值y可用下面程序计算： y=polyval(p,x) 例： x=1:20; y=x+3*sin(x); p=polyfit(x,y,6) xi=linspace(1,20,100); z=polyval(p,xi); % 多项式求值函数
倒立摆系统的建模及Matlab 仿真 1.系统的物理模型 考虑如图(1)面内运动的二维问题。 图(1)倒立摆系统 假定倒立摆系统的参数如下。 摆杆的质量：m=0.1g 摆杆的长度：l =1m 小车的质量： M=1kg 重力加速度：g=9.8m/2s 摆杆的质量在摆杆的中心。 设计一个控制系统，使得当给定任意初始条件(由干扰引起)时，最大超调量δ ≤10%，调节时 间ts ≤4s ，通过小车的水平运动使倒立摆保持在垂直位置。 2.系统的数学模型 2.1建立倒置摆的运动方程并将其线性化。 为简化问题，在数学模型中首先假设：1)摆杆为刚体；2）忽略摆杆与支点之间的摩擦；3）忽略小车与接触面间的摩擦。 设小车瞬时位置为z,摆心瞬时位置为（θsin l z +）,在u 作用下，小车及摆均产生加速远动，根据牛顿第二定律，在水平直线远动方向的惯性力应与u 平衡，于是有 u l z dt d m dt z d M =++)sin (22 22θ 即： u ml ml z m M =-++θθθθsin cos )(2&&&&& ① 绕摆轴转动的惯性力矩与重力矩平衡，因而有
θθθsin cos )sin (22mgl l l z dt d m =??? ????+ 即： θθθθθθθsin cos sin cos cos 22g l l z =-+&&&&& ② 以上两个方程都是非线性方程，为求得解析解，需作线性化处理。由于控制的目的是保持倒立摆直 立，在试驾合适的外力条件下，假定θ很小，接近于零时合理的，则1cos ,sin ≈≈θθθ，且可忽略θ θ2&项。于是有 u ml z m M =++θ&&&& )( ③ θθg l z =+&&&& ④ 联立求解可得 u Ml Ml m M u M M mg z 1)(1 -+=+- =θθθ&&&& 2.2列写系统的状态空间表达式。 选取系统变量4321,,,x x x x ， []T x x x x x 4321,,,=则 u Ml x Ml m M x x x u M x M mg x x x 1 )(134433221-+= =+-==&&&& 即 []Cx x x y Bu Ax u Ml M x Ml g m M M mg z z dt d x ===+=?????? ? ???????-+?????????? ??? ? +- =???? ????????=000110100)(0 010 0000000 1 1θθ&&& 代入数据计算得到： [][]0,0001,1010,01100 1000010000 1 0==-=? ? ??? ? ??? ???-=D C B A T
m a t l a b实现插值法和 曲线拟合
插值法和曲线拟合 电子科技大学 摘要：理解拉格朗日多项式插值、分段线性插值、牛顿前插，曲线拟合，用matlab编程求解函数，用插值法和分段线性插值求解同一函数，比较插值余项；用牛顿前插公式计算函数，计算函数值；对于曲线拟 合，用不同曲线拟合数据。 关键字：拉格朗日插值多项式；分段线性插值；牛顿前插；曲线拟合 引言： 在数学物理方程中，当给定数据是不同散点时，无法确定函数表达式，求解函数就需要很大的计算量，我们有多种方法对给定的表格函数进行求解，我们这里，利用插值法和曲线拟合对函数进行求解，进一步了解函数性质，两种方法各有利弊，适合我们进行不同的散点函数求解。 正文： 一、插值法和分段线性插值 1拉格朗日多项式原理 对某个多项式函数，已知有给定的k + 1个取值点： 其中对应着自变量的位置，而对应着函数在这个位置的取值。 假设任意两个不同的x j都互不相同，那么应用拉格朗日插值公式所得到的拉格朗日插值多项式为： 其中每个为拉格朗日基本多项式（或称插值基函数），其表达式为： [3] 拉格朗日基本多项式的特点是在上取值为1，在其它的点 上取值为0。 2分段线性插值原理 给定区间[a,b], 将其分割成a=x 0 3.1 曲线拟合的线性最小二乘法及其MATLAB 程序 例3.1．1 给出一组数据点),(i i y x 列入表3-1中，试用线性最小二乘法求拟合曲线，并估计其误差，作出拟合曲线. 表3-1 例3.1.1的一组数据),(y x 解 （1）在MATLAB 工作窗口输入程序 >> x=[-2.5 -1.7 -1.1 -0.8 0 0.1 1.5 2.7 3.6]; y=[-192.9 -85.50 -36.15 -26.52 -9.10 -8.43 -13.12 6.50 68.04]; plot(x,y,'r*'), legend('实验数据(xi,yi)') xlabel('x'), ylabel('y'), title('例3.1.1的数据点(xi,yi)的散点图') 运行后屏幕显示数据的散点图（略）. （3）编写下列MA TLAB 程序计算)(x f 在),(i i y x 处的函数值，即输入程序 >> syms a1 a2 a3 a4 x=[-2.5 -1.7 -1.1 -0.8 0 0.1 1.5 2.7 3.6]; fi=a1.*x.^3+ a2.*x.^2+ a3.*x+ a4 运行后屏幕显示关于a 1,a 2, a 3和a 4的线性方程组 fi =[ -125/8*a1+25/4*a2-5/2*a3+a4, -4913/1000*a1+289/100*a2-17/10*a3+a4, -1331/1000*a1+121/100*a2-11/10*a3+a4, -64/125*a1+16/25*a2-4/5*a3+a4, a4, 1/1000*a1+1/100*a2+1/10*a3+a4, 27/8*a1+9/4*a2+3/2*a3+a4, 19683/1000*a1+729/100*a2+27/10*a3+a4, 5832/125*a1+324/25*a2+18/5*a3+a4] 编写构造误差平方和的MATLAB 程序 >> y=[-192.9 -85.50 -36.15 -26.52 -9.10 -8.43 -13.12 6.50 68.04]; fi=[-125/8*a1+25/4*a2-5/2*a3+a4, -4913/1000*a1+289/100*a2-17/10*a3+a4, -1331/1000*a1+121/100*a2-11/10*a3+a4, -64/125*a1+16/25*a2-4/5*a3+a4, a4, 1/1000*a1+1/100*a2+1/10*a3+a4, 27/8*a1+9/4*a2+3/2*a3+a4, 19683/1000*a1+729/100*a2+27/10*a3+a4, 5832/125*a1+324/25*a2+18/5*a3+a4]; fy=fi-y; fy2=fy.^2; J=sum(fy.^2) 运行后屏幕显示误差平方和如下 J= (-125/8*a1+25/4*a2-5/2*a3+a4+1929/10)^2+(-4913/1000*a1+2 89/100*a2-17/10*a3+a4+171/2)^2+(-1331/1000*a1+121/100*a2-11/10*a3+a4+723/20)^2+(-64/125*a1+16/25*a2-4/5*a3+a4+663/25)^2+(a4+91/10)^2+(1/1000*a1+1/100*a2+1/10*a3+a4+843/100)^2+(27/8*a1+9/4*a 2+3/2*a3+a4+328/25)^2+(19683/1000*a1+729/100*a2+27/10*a3+a4-13/ 2)^2+(5832/125*a1+324/25*a2+18/5*a3+a4-1701/25)^2 为求4321,,,a a a a 使J 达到最小，只需利用极值的必要条件0=??k a J )4,3,2,1(=k ， 2,3,4次曲线拟合matlab程序 【程序代码】 clf reset H=axes('unit','normalized','position',[0,0,1.5,1],'visible','off'); set(gcf,'currentaxes',H); str='\fontname{微软雅黑}2,3,4次曲线拟合程序'; text(0.17,0.9,str,'fontsize',15);%这是设置字体位置的 h_fig=get(H,'parent'); set(h_fig,'unit','normalized','position',[0.1,0.2,0.8,0.5]);%这是设置出现窗口的大小的 h_axes=axes('parent',h_fig,'unit','normalized','position',[0.1,0.15,0.55,0.7],'xlim',[015],'ylim',[0 1.8],'fontsize',8); h_text=uicontrol(h_fig,'style','text','unit','normalized','position',[0.69,0.90,0.24,0.03],'horizontal','left','s tring',{'左区间'}); h_text1=uicontrol(h_fig,'style','text','unit','normalized','position',[0.69,0.75,0.24,0.03],'horizontal','left',' string',{'右区间'}); h_text2=uicontrol(h_fig,'style','text','unit','normalized','position',[0.69,0.62,0.24,0.03],'horizontal','left',' string',{'步长'}); h_text3=uicontrol(h_fig,'style','text','unit','normalized','position',[0.69,0.48,0.24,0.03],'horizontal','left',' string',{'拟合矩阵'}); h_edit=uicontrol(h_fig,'style','edit','unit','normalized','position',[0.69,0.82,0.24,0.08], 'horizontal','left','callback',['a=str2num(get(gcbo,''string''));','t=a:n:b;','x=x;','p2=polyfit(t,x,2);','f2=poly val(p2,t);','p3=polyfit(t,x,3);','f3=polyval(p3,t);','p4=polyfit(t,x,4);','f4=polyval(p4,t);','plot(t,x,t,f2,t,f3,t, f4)']); h_edit1=uicontrol(h_fig,'style','edit','unit','normalized','position',[0.69,0.67,0.24,0.08], 'horizontal','left','callback',['b=str2num(get(gcbo,''string''));','t=a:n:b;','x=x;','p2=polyfit(t,x,2);','f2=poly val(p2,t);','p3=polyfit(t,x,3);','f3=polyval(p3,t);','p4=polyfit(t,x,4);','f4=polyval(p4,t);','plot(t,x,t,f2,t,f3,t, f4)']); h_edit2=uicontrol(h_fig,'style','edit','unit','normalized','position',[0.69,0.54,0.24,0.08], 'horizontal','left','callback',['n=str2num(get(gcbo,''string''));','t=a:n:b;','x=x;','p2=polyfit(t,x,2);','f2=poly val(p2,t);','p3=polyfit(t,x,3);','f3=polyval(p3,t);','p4=polyfit(t,x,4);','f4=polyval(p4,t);','plot(t,x,t,f2,t,f3,t, f4)']); h_edit3=uicontrol(h_fig,'style','edit','unit','normalized','position',[0.69,0.38,0.24,0.1], 'horizontal','left','callback',['x=str2num(get(gcbo,''string''));','t=a:n:b;','x=x;','p2=polyfit(t,x,2);','f2=poly val(p2,t);','p3=polyfit(t,x,3);','f3=polyval(p3,t);','p4=polyfit(t,x,4);','f4=polyval(p4,t);','plot(t,x,t,f2,t,f3,t, f4)']); h_push1=uicontrol(h_fig,'style','pushbutton','unit','normalized','position',[0.69,0.24,0.12,0.08],'string',' grid on','callback','grid on'); h_push2=uicontrol(h_fig,'style','pushbutton','unit','normalized','position',[0.69,0.15,0.12,0.08],'string',' grid off','callback','grid off'); h_push3=uicontrol(h_fig,'style','pushbutton','unit','normalized','position',[0.81,0.15,0.12,0.08],'string','退出','callback','exit'); h_push4=uicontrol(h_fig,'style','pushbutton','unit','normalized','position',[0.81,0.24,0.12,0.08],'string','关闭','callback','close(gcbf)'); 【操作界面】 ECE451 Controll Engineering Inverted pendulum 09/29/2013 Introduction: Inverted pendulum is a typical fast, multi-varaibles, nonlinear, unstable system, it has significant meaning. We choose the PID controller to fot the inverted pendulum. Assume the input is a step signal , the gravitational acceleration g=9.8m/s^2 and linearize the nonlinear model around the operating point. 1.Mathematic Modling M mass of the car 0.5 kg m mass of the pendulum 0.2 kg b coefficient of friction for cart 0.1 N/m/sec l length to pendulum center of mass 0.3 m I mass moment of inertia of the pendulum 0.006 kg.m^2 F force applied to the cart x coordinate of cart position θpendulum angle from vertical (down) N and F are the force from horizontal and vertical direction. ) Force analysis Consider the horizontal direction cart force, we get the equation: Consider the horizontal direction pendulum force, we get the equation: To get rid of P and N, we get this equation: Merge these two equations, about to P And N, to obtain a second motion equation: u to represent the controlled object with the input force F, linearized two motion equations Apply Laplace transform to the equation above The transfer function of angle and position Matlab 线性回归（拟合） 对于多元线性回归模型： e x x y p p ++++=βββ 110 设变量12,,,p x x x y 的n 组观测值为 12(,,,) 1,2,,i i ip i x x x y i n = ． 记 ??????? ? ?=np n n p p x x x x x x x x x x 2 1 222211121111 1，?? ?? ??? ??=n y y y y 2 1 ，则???? ?? ? ??=p ββββ 10 的估计值为 y x x x b ')'(?1-==β (11.2) 在Matlab 中，用regress 函数进行多元线性回归分析，应用方法如下： 语法：b = regress(y, x) [b, bint, r, rint, stats] = regress(y , x) [b, bint, r, rint, stats] = regress(y , x, alpha) b = regress(y, x)，得到的1+p 维列向量b 即为(11.2)式给出的回归系数β的估计值． [b, bint, r, rint, stats]=regress(y , x) 给出回归系数β的估计值b ，β的95％置信区间（(1)2p +?向量）bint ，残差r 以及每个残差的95％置信区间（2?n 向量）rint ；向量stats 给出回归的R 2 统计量和F 以及临界概率p 的值． 如果i β的置信区间（bint 的第1i +行）不包含0，则在显著水平为α时拒绝0i β=的假设，认为变量i x 是显著的． [b, bint, r, rint, stats]=regress(y , x, alpha) 给出了bint 和rint 的100(1-alpha)%的置信区间． 三次样条插值函数的MATLAB 程序 matlab 的spline x = 0:10; y = sin(x); %插值点 xx = 0:.25:10; %绘图点 yy = spline(x,y ,xx); plot(x,y,'o',xx,yy) 基于matlab的倒立摆的仿真与设计姓名:贾永伟专业:测控技术与仪器学号:1123105950 年级:2011级 摘要：倒立摆系统是一个典型的快速、多变量、非线性、不稳定系统，对倒 立摆的控制研究无论在理论上和方法上都有深远的意义。 本论文以实验室原有的直线一级倒立摆实验装置为平台，重点研究其PID控制方法，设计出相应的PID控制器，并将控制过程在MATLAB上加以仿真。 关键词：一级倒立摆，PID，MATLAB仿真 一、倒立摆模型的研究意义 倒立摆系统的研究能有效的反映控制中的许多典型问题：如非线性问题、鲁棒性问题、镇定问题、随动问题以及跟踪问题等。通过对倒立摆的控制，用来检验新的控制方法是否有较强的处理非线性和不稳定性问题的能力。同时，其控制方法在军工、航天、机器人和一般工业过程领域中都有着广泛的用途，如机器人行走过程中的平衡控制、火箭发射中的垂直度控制都有重要意义 倒立摆控制系统是一个复杂的、不稳定的、非线性系统，是进行控制理论教学及开展各种控制实验的理想实卫星飞行中的姿态控制等。故其研究意义广泛。 二、倒立摆模型的数学建模 质量为m的小球固结于长度为L的细杆（可忽略杆的质量）上，细杆又和质量为M的小车铰接相连。由经验知：通过控制施加在小车上的力F（包括大小和方向）能够使细杆处于θ＝0的稳定倒立状态。在忽略其他零件的质量以及各种摩擦和阻尼的条件下，推导小车倒立摆系统的数学模型 分析过程如下： 如图所示，设细杆摆沿顺时针方向转动为正方向，水平向右方向为水平方向上的正方向。当细杆摆顺时针往右运动时水平方向施加的力应该为水平向右。 现对小车和细杆摆分别进行隔离受力分析： 1 曲线拟合的线性最小二乘法及其MATLAB 程序 例7.2.1 给出一组数据点),(i i y x 列入表7–2中，试用线性最小二乘法求拟合曲线，并用（7.2），（7.3）和（7.4）式估计其误差，作出拟合曲线. 表7–2 例7.2.1的一组数据),(y x 解 （1）在MATLAB 工作窗口输入程序 >> x=[-2.5 -1.7 -1.1 -0.8 0 0.1 1.5 2.7 3.6]; y=[-192.9 -85.50 -36.15 -26.52 -9.10 -8.43 -13.12 6.50 68.04]; plot(x,y,'r*'), legend('实验数据(xi,yi)') xlabel('x'), ylabel('y'), title('例7.2.1的数据点(xi,yi)的散点图') 运行后屏幕显示数据的散点图（略）. （3）编写下列MA TLAB 程序计算)(x f 在),(i i y x 处的函数值，即输入程序 >> syms a1 a2 a3 a4 x=[-2.5 -1.7 -1.1 -0.8 0 0.1 1.5 2.7 3.6]; fi=a1.*x.^3+ a2.*x.^2+ a3.*x+ a4 运行后屏幕显示关于a 1,a 2, a 3和a 4的线性方程组 fi =[ -125/8*a1+25/4*a2-5/2*a3+a4, -4913/1000*a1+289/100*a2-17/10*a3+a4, -1331/1000*a1+121/100*a2-11/10*a3+a4, -64/125*a1+16/25*a2-4/5*a3+a4, a4, 1/1000*a1+1/100*a2+1/10*a3+a4, 27/8*a1+9/4*a2+3/2*a3+a4, 19683/1000*a1+729/100*a2+27/10*a3+a4, 5832/125*a1+324/25*a2+18/5*a3+a4] 编写构造误差平方和的MATLAB 程序 >> y=[-192.9 -85.50 -36.15 -26.52 -9.10 -8.43 -13.12 6.50 68.04]; fi=[-125/8*a1+25/4*a2-5/2*a3+a4, -4913/1000*a1+289/100*a2-17/10*a3+a4, -1331/1000*a1+121/100*a2-11/10*a3+a4, -64/125*a1+16/25*a2-4/5*a3+a4, a4, 1/1000*a1+1/100*a2+1/10*a3+a4, 27/8*a1+9/4*a2+3/2*a3+a4, 19683/1000*a1+729/100*a2+27/10*a3+a4, 5832/125*a1+324/25*a2+18/5*a3+a4]; fy=fi-y; fy2=fy.^2; J=sum(fy.^2) 运行后屏幕显示误差平方和如下 J= (-125/8*a1+25/4*a2-5/2*a3+a4+1929/10)^2+(-4913/1000*a1+2 89/100*a2-17/10*a3+a4+171/2)^2+(-1331/1000*a1+121/100*a2-11/10*a3+a4+723/20)^2+(-64/125*a1+16/25*a2-4/5*a3+a4+663/25)^2+(a4+91/10)^2+(1/1000*a1+1/100*a2+1/10*a3+a4+843/100)^2+(27/8*a1+9/4*a 2+3/2*a3+a4+328/25)^2+(19683/1000*a1+729/100*a2+27/10*a3+a4-13/ 2)^2+(5832/125*a1+324/25*a2+18/5*a3+a4-1701/25)^2 为求4321,,,a a a a 使J 达到最小，只需利用极值的必要条件0=??k a J )4,3,2,1(=k ，曲线拟合的线性最小二乘法及其MATLAB程序
2,3,4次曲线拟合matlab程序
Invertedpendulum倒立摆的matlab建模
Matlab线性回归(拟合)
基于matlab的倒立摆仿真设计
曲线拟合_线性最小二乘法及其MATLAB程序
相关主题
文本预览