《数学建模论文制动器试验台的控制方法分析 .doc》由会员分享，可在线阅读，更多相关《数学建模论文制动器试验台的控制方法分析 .doc（22页珍藏版）》请在三一办公上搜索。
1、制动器试验台的控制方法分析摘要：本论文研究的是与汽车制造工业紧密相连的有关汽车制动器试验台系统的模拟仿真问题，围绕着如何在制动器试验台上精准的模拟路试试验展开。第一问在基本假设的前提下，运用物理学相关知识，算得的等效惯量为51.99889 ，第二问主要运用转动惯量的计算公式，算出三个飞轮各自的转动惯量后再进行组合，得到8个机械惯量，随后求出问题一中的补偿惯量为12或-18。第三问对试验的机理进行分析，建立起电流I依赖于主轴瞬时转速n和主轴的瞬时扭矩M的数学模型，再根据问题一、问题二的条件和题中所给数据，求出驱动电流为=174.825A（对应补偿惯量为12）, = -262.238A（对应补偿惯
2、量为-18）。第四问中首先跟据路试与试验台模拟时制动器消耗能量一致原则确立了的评价指标。在模型一中运用矩形公式通过积分的方法求出，再运用能量公式求出，解得a=94.298%。但考虑到该模型可能会低估模拟误差，模型二中改进的方法为，对时间轴t进行划分，每段用梯形公式积分，算出每一小段对应的、，进而求的比值b，然后对所有的b值求期望和方差，求得=91.32%，=0.1834。第五问我们首先用灰色理论建立了以主轴瞬时扭矩为影响因素的GM(1,1)模型，随后又尝试了以主轴瞬时扭矩、瞬时转速为影响因素的GM(1，2)模型，虽然经过多次调整，但两种模型的残差均偏大，遂舍弃该模型，应用广义预测控制（GPC）
3、的基本方法原理，创建一套广义预测控制系统，用来控制电流，数值模拟精度高达96.01%，结果令人满意，但是该算法无法避免转速出现短期震荡的现象。针对第五问模型的不足，第六问我们引入机械惯量系统动力方程建立模型研究双闭环直流电机系统的转速变化规律。改进的优点在于：可以通过改变电机控制系统的各种给定信号，使电机输出转速变化能准确地模拟机械惯量系统的转速变化，模拟效果很好。关键词：梯形公式 广义预测控制（GPC） 机械惯量系统动力方程 双闭环直流电机系统一、模型的基本假设：1 假设路试时 轮胎与地面的摩擦力为无穷大，因此轮胎与地面无滑动。2 假设模拟实验中，主轴的角速度始终与车轮角速度一致。3 在计算
4、等效惯量的时候，忽略车轮自身转动具有的能量。4 假设试验台采用的电动机的驱动电流与其产生的扭矩成正比，且试验台工作时主轴的瞬时转速与瞬时扭矩是可观测的离散量。5 不考虑观测误差、随机误差和连续问题离散化所产生的误差。二、符号说明：J转动惯量m质量M扭矩角速度角加速度I电流F载荷v路试车轮的前进速度注：上述符号为基本符号，题中其余符号在建模时均有说明。三、模型的建立和求解：问题一：由题意可知，路试车辆的指定车轮在制动时承受载荷，将这个载荷在车辆平动时具有的能量（忽略车轮自身转动具有的能量）等效地转化为试验台上飞轮和主轴等机构转动时具有的能量，与此能量相应的转动惯量(以下转动惯量简称为惯量)在本题
5、中称为等效的转动惯量。因此我们可以得到以下的方程：经过整理我们得到：带入数值，其中r=0.286m，F=6230N,算得等效惯量为51.99889 问题二：由物理学知识我们得到：飞轮的转动惯量计算公式为：经过整理可得：其中是钢材密度，取值为，R代表外直径，取值为1m，r代表内直径，取值为0.2m，代表飞轮的厚度，取之分别为0.0392m、0.0784m、0.1568m，i=1，2，3。代入数值计算得 ：可以组合成种机械惯量，分别为：10，40，70，100，130，160，190，220 。题1中的等效惯量为51.99889 ，近似为52，考虑到电动机能补偿的能量相应的惯量范围为-30, 30
6、 kgm2 ,当机械惯量取值为40时，补偿惯量为12，当机械惯量取值为70时，补偿惯量为 -18。问题三：图3.1 试验台电机工作示意图由图可知，在试验台模拟仿真，制动阶段开始，当需要电动机补偿惯量时，设电动机产生的扭矩为，制动力产生的扭矩为M，n为飞轮的转速， m为飞轮质量，J为飞轮组实际的转动惯量，为等效惯量，为补偿惯量，有：整理得到：设电动机的驱动电流为I，由题前假设可知，电动机的驱动电流与成正比，即：，,K为比例系数，取1.5，最终可得： （ 3.1.1式） （ 3.2.1式） （ 3.3.1式）可根据不同情况选择具体的模型，以上三个公式相互之间均等价。另外，当为变量，没有办法直接求解
7、时，我们可以用一阶差商来代替倒数，具体做法如下：将时间段t分为i个小段，每段长度为，有，算出其中每一小段的变化量和与其相对应的电流，并对整个区间变化量和电流求期望E，有整理可得： （ 3.1.2式） （ 3.2.2式）在第一问和第二问的条件下，假设制动减速度为常数时，由（ 3.2.1式）：带入具体数值，r=0.286m, =50km/h, =12,= -18解得=174.825A, = -262.238A问题四：由与所设计的路试等效的转动惯量为48 kgm2，机械惯量为35 kgm2，我们可以得到电流补偿惯量为13 kgm2。对于题中给出的控制方法，我们用数量指标a来评价：其中：M为试验台上制
8、动器产生的扭矩，为t时间内飞轮组转过的角度，W为t时间内路试时制动器所做功，、分别为t时间内飞轮组的初角速度、末角速度,J为飞轮组机械惯量，为等效惯量，为电流补偿惯量。模型一：对附表内的数据进行分析，得到飞轮组角速度与时间t的关系图线（见图4.1.1）和制动力扭矩M与时间t的关系图线（见图4.1.2），图4.1.1图4.1.2由图可知，与t近似呈线性关系，但是在t=0.1s处出现剧烈抖动，M的值在的区间内迅速上升，在t=0.1s后在某一均值处上下频繁波动，我们推测出现此现象的原因为：在制动器刚开始制动的时候，所加载荷不可能一下增加到最大值，有一个上升的过程，此阶段为不稳定阶段，而在t=0.1s
9、后我们看到M的值虽有波动但幅度不大，可认为此时制动力已达最大值，处于相对稳定阶段。而的抖动很可能是由于共振等物理原因，难以避免。我们将时间轴t划分为一个个小区间，区间长度为，取出t时刻的一个区间：如图，与t轴围成的面积即为时间内飞轮组转过的角度，当无限小时，可近似认为，用矩形公式可算得=，积分可得：小区间长度取为10ms,用Matlabb编程求得a=94.298%,可见第四问中的控制方法效果还是比较令人满意的。模型二：模型一中采用积分的方法求出了区间内飞轮组转过的角度，但是这种方法不够精确，可能低估了模拟误差。好比第一小段的误差为 2，第二小段的误差为-2，但是相加以后，误差反而为零了(可参见
10、下图)。为了尽量避免这种情况，我们设计了模型二针对模型一的不足加以改进。具体做法为利用梯形公式，算出每个小区间所对应的，然后分析得到的，求出期望和方差，其中：为一阶线性拟合下的角加速度，取为100ms，求得b的期望值为91.32%，方差为0.1834，可见第四问的控制方法效果只能算作良好，并不十分理想。第五问：我们首先用灰色理论先后建立了以主轴瞬时扭矩为影响因素的GM(1,1)模型，随后又尝试了以主轴瞬时扭矩、瞬时转速为影响因素的GM(1，2)模型，虽然经过多次调整，但两种模型的残差均偏大，遂舍弃该模型，应用广义预测控制（GPC）的基本方法原理，创建一套广义预测控制系统，用来控制电流，具体方程
11、如下：1. 模型建立我们在GPC中建立有关电流I（被控对象）的数学模型，用受控自回归积分滑动平均（CARIMA）模型来描述。 （2.1）式中式中，表示顺延一个采样周期的相应量，称为后移算子；为差分算子；是一个不相关的随机序列，表示白噪声的影响。A、B、C都是的多项式，其中多项式的若干首项元素，表示电流相应的纯滞后的步数，可以为零。假设，用以更好的突出方法原理。式（2.1）也可用z传递函数给出电流模型的描述。由输入u对输出y之间的传递函数为 （2.2）为了利用模型式（2.1）导出j步之后输出的预测值，首先引入下述Diophantine方程（恒等式）： （2.3）式中，由和预测长度j唯一确定在式（
12、2.1）中两端同乘，得利用式（2.3），可以得出时刻的输出量关注，的形式，可知与，有关，只与，有关，而与，有关。由于白噪声未知，所以最合适的预测值可由下式表示 （2.4）由此，式（2.4）即为GPC的预测模型，我们根据已知的瞬时扭矩、瞬时转速，及未来的扭矩和转速，就可以预测出电动机电流。在式（2.4）中，记，结合式（2.2），式（2.3）可得因此多项式中前j项的系数正是装置阶跃响应前j项的采样值，记做，。若把展开写作则有。为了由式（2.4）计算预测值，必须首先知道，对于不同的这相当于并行的求解一组Diophantine方程（2.3），计算量是相当惊人的。为此克拉克给出了一个，的递推算法。首先，
13、根据式（2.3）可写出恒等式两式相减得到记：则可得显然，必须有并且由于的首项系数为1，必然有这一系数的递推关系亦可用向量式记为式中此外还可得的递推公式为当时，方程（2.3）为 故可取，为，的初值。这样，便可按照下式递推计算 （2.5）2反复循环优化在GPC中，t时刻的优化性能指标具有以下形式 （2.6）式中，E为数学期望；为装置输出的期望值；和分别为优化时域的始值和终值；NU为控制时域，即在NU步后控制量不再变化为控制加权系数，为简化，一般常可假设其为常识。该性能指标采用了长时间段预测的概念，把所要优化的方差从一个时间点扩展到一段时域，其中应该大于装置的时滞数，而应大到对象的动态特性能充分表现
14、出来。由于以多步预测优化代替了一步预测优化，即使对时滞估计不当或时滞发生变化，仍能从总体优化中得到合理的控制，这是GPC对模型不精确性具有鲁棒性的重要原因。为了更清楚的表达，在下面的式子中将假设。利用预测模型式(2.4)，可得其中使上述性能指标最优的解为：式中，为对象输出的期望值，其形式如下：，即时最优控制量则由以下式给出 （2.7）式中，是矩阵的第一行。至此GPC算法的主要步骤已经表达完毕，但需要指出的是在GPC算法中，对模型参数和的辨识用的是递推最小二乘法，这一步必须在每个循环周期的第一步完成。3仿真模型的离散化题中制动器开始制动后，载荷需要一段时间才能上升到最大值，因此存在一个不稳定阶段
15、，并且由于制动器本身性能复杂，再加上测量装置本身的限制，因此可观测数据读入有一定的滞后性，在仿真时将其简化为一个带有纯滞后的一阶惯性环节，其传递函数为：式中：K为静态增益；T为惯性时间；为纯滞后时间。为了进行计算机控制，我们将连续系统离散化以便获得离散信号。分析离散系统可以用古典的Z变换法求出系统的脉冲传递函数G(z),再从脉冲传递函数G(z)求出离散系统的动态方程差分方程，就能方便地对离散系统进行数字仿真及控制了。用z变换法求连续系统离散化后的脉冲传递函数及差分方程，要在系统中加入采样器和保持器。在这里我们采用零阶保持器。设零阶保持器的传递函数为：式中，为采样周期。则广义对象传递函数为，对其
16、进行Z变换可得到与G(s)等价的离散模型：由上式得上式两边进行z反变换得令,， 则上式可整理为4.仿真试验研究 为突出方法原理，在仿真实验中我们略去了参数辨识这一步，直接选取对象的传递函数为，采样周期为0.1秒，并加入一个不相关的随机序列，表示一类白噪声，离散化后为，仿真运行后的结果如图.再用第四问中模型二的评价指标b来检验模拟效果：用Matlab计算可得，b=96.01%，可见控制的效果令人满意。图中曲线一为原数据曲线，曲线二为模拟的数据曲线 仿真结果表明:本模型采用的广义控制算法对仿真中加入的白噪声具有很强的抗干扰能力，对模型参数的变化也有较强的鲁棒性, 符合题目的具体要求，能够实现基于该
17、模型的广义预测控制。但是缺点在于该算法无法避免转速出现短期震荡的现象，因此在第六问中我们又重新设计了一个新的模型。第六问：针对第五问模型中存在的不足，我们引入机械惯量系统动力方程建立模型。1.本题中电动机用以模拟机械装置传动系转动，当飞轮以角速度稳定旋转时，设在零时刻受到方向与大小恒定的力矩M作用，则根据电机拖动的动力学方程（G为重力，D为惯性直径）得在t时刻的转速变化n(t)为: （6.1.1式）式中为零时刻飞轮初始转速，n(t)为t时刻飞轮转速，由式(6.1)经拉氏变换，可得: （6.1.2式）2. 双闭环直流电机系统动力方程当受到恒定力矩M的作用时，其负载电流在瞬时发生变化，从而导致电机
18、转速发生变化。把电流环作为一个环节考虑，可得到图6.1（转速、电流双闭环直流调速系统的动态结构变形图）7。图6.1 双闭环直流调速系统的动态结构变形图由该图可得零初始条件下，系统输出转数变化与负载电流变化之间的闭环传递函数: ( 6.2.1 式)其中，为电机额定励磁下电动势转速比，为电机额定励磁下转矩电流比，且，为负载力矩作用下电机转速变化量的拉氏变换，为调速系统负载电流变化量的拉氏变换，为转速反馈系数，为电流闭环传递函数。在设计转速调节器ASR时，可把已设计好的电流环看作是转数调速系统中的一个环节，它的等效传递函数为: （6.2.2式）把转速环校正成典型11型系统，ASR应采用PI调节器6，
19、其传递函数为: （6.2.3式）式中为转速调节器的比例系数，为转速调节器的超前时间常数。按照典型H型系统的参数选择方法，采用“振荡指标法”中所用的闭环幅频特性峰值最小准则5可得: （6.2.4式） （6.2.5式）式中为转速环等效时间常数，其为小惯性环节和合并成的时间常数，即，h为调速系统中频宽，为系统机电时间常数。 将（6.2.2式）、（6.2.3式）、（6.2.4式）、（6.2.5式）带入（6.2.1式），并考虑，可得 （6.2.6式）由上式经拉氏反变换，即得到双闭环直流调速系统在恒定负载转矩作用下的转速变化规律:（6.2.7式），式中，在零时刻恒定力矩M开始作用，对应于不同的h值，各系数
20、C、p都有相应的取值。3. 电惯量控制的实现方法8由式6.1.1与式6.2.7的和可以看出，由于电机控制系统的结构形式不同，在同样的恒定转矩作用下，他励直流电机系统与双闭环直流电机系统的转速变化规律大相径庭。在传统的制动器试验台中，电机控制系统不存在闭环反馈，只有飞轮与制动恒定力矩两个环节，因此系统的转速变化规律是按照式6.1.1进行的。当引入转速、电流反馈双闭环控制后，系统数学模型即输入与输出之间的关系发生改变，导致系统转速变化规律改变，即得式6.2.7。设想减小或取消传统制动器试验台的飞轮转动惯量，得到的新系统为符合传统试验台的真实工作情况，可以通过改变电机控制系统的各种给定信号，使电机输
21、出转速变化能准确地模拟传统机械惯量系统的转速变化，由此即构成了实现电惯量的制动器试验台。四、模型的评价与分析1. 问题一、问题二、问题三均用常规的物理学知识入手，列出方程组，建立起数学模型；2. 问题四以 为评价指标，第一个模型用矩形公式积分的方法求出进而得出a值，为了避免由于连续积分造成的误差，更精确的体现实验模拟的细节，我们推出了模型二，将时间轴细分为一个个小段，在每一段上用梯形公式进行积分，分别求出每一段对应的指标值b,最后求出b的期望和方差；3. 问题五我们首先尝试用灰色理论建立模型进行模拟，首先建立了以主轴瞬时扭矩为影响因素的GM(1,1)模型，后又尝试了以主轴瞬时扭矩、瞬时转速为影
22、响因素的GM(1，2)模型，虽然经过多次调整，但两种模型的残差均偏大，遂舍弃该模型。而后我们应用广义预测控制（GPC）的基本方法原理，创建一套广义预测控制系统，用来控制电流。在进行数据仿真后我们及时对模拟效果进行评价，结果证明本模型采用的广义控制算法对仿真中加入的白噪声具有很强的抗干扰能力，对模型参数的变化也有较强的鲁棒性, 符合题目的具体要求，能够实现基于该模型的广义预测控制，但是缺点在于该算法无法避免转速出现短期震荡的现象。4. 针对第五问模型中存在的不足，我们引入机械惯量系统动力方程建立模型研究双闭环直流电机系统的转速变化规律。在传统的制动器试验台中，电机控制系统不存在闭环反馈，只有飞轮
23、与制动恒定力矩两个环节，当引入转速、电流反馈双闭环控制后，系统数学模型即输入与输出之间的关系发生改变，导致系统转速变化规律改变，即得式6.2.7。改进的优点在于：可以通过改变电机控制系统的各种给定信号，使电机输出转速变化能准确地模拟传统机械惯量系统的转速变化，由此即构成了实现电惯量的制动器试验台，模拟结果说明效果很好。参考文献1Clarke D W, Mohtadi C, Tuffs P S, Generalized predictive controlPartand Part ,Automatica, 23(2):137160, 19872阮学斌，广义预测控制在空调监控系统中的应用，控制工程
24、，14(1):18-20，20073钱积新，赵均，徐祖华，预测控制，北京：化学工业出版社，20074师五喜，广义预测控制中Diophantine方程的显式解，自动化学报，32(1):148-153，20065张青、张鸣，转速的双闭环PID控制J，湖北工学院学报，2004，第19卷第3期;6.朱晓锦、张为公等，基于计算机测控的汽车同步器试验系统研制J，测控技术，2000，第19卷第11期:45487.胡寿松，自动控制理论M，国防工业出版社，1994:1751888鹏战争，制动器试验台转速控制研究J，4042附录一、第五问中模型的Matlab程序：M=4;P=6;A=1 -0.496585;A1=
25、1 -1.28585 0.19685;B=0 0.93 1.21;la=length(A);lb=length(B);na=la-1;nb=lb-1;S=zeros(P,na+1);S(1,:)=1.10085 -1.1;for j=1:P-1 for i=1:na S(j+1,i)=S(j,i+1)-A1(i+1)*S(j,1); end S(j+1,na+1)=-S(j,1)*A1(na+2);endG=zeros(P,P+nb);for j=1:nb+1 G(1,j)=B(j);endfor j=1:P-1 for i=1:j+nb if ij bij=0; else bij=B(i-j
26、+1); end G(j+1,i+1)=G(j,i+1)+S(j,1)*bij; endendfor j=2:P G(j,1)=B(1);endG1=zeros(P,P);for i=1:P for j=1:i G1(i,j)=G(i-j+1,i-j+1); endendG11=zeros(P,M);for i=1:P for j=1:M G11(i,j)=G1(i,j); endendG2=zeros(P,nb);for i=1:P for j=1:nb G2(i,j)=G(i,nb+1+i-j); endendu1=zeros(P,1); u0=zeros(nb,1);y0=zeros(n
27、a+1,1);y=;u01=zeros(na+2,1);u=;for j=1:na+1u(j)=0;endu=u;yr=;yr(1)=0;yr1=zeros(P,1);sy=zeros(P,1);d1=zeros(1,M);d1(1,1)=0.9;f=zeros(P,1);w=-27;r=0.9;a=0.9;for t=1:200 for j=1:P yr(t+j)=a*yr(t+j-1)+(1-a)*w;yr1(j)=yr(t+j);end for j=1:P syy=0; for i=1:na+1 syy=syy+S(j,i)*y0(i); end sy(j)=syy;endf=G2*u0
28、+sy;lamda=r*eye(M);lamda(1,1)=1;u1=inv(G11*eye(P)*G11+lamda)*G11*eye(P)*(yr1-f);ut1=d1*u1;ut=ut1+u01(na+2);for i=1:na+1 u01(i)=u01(i+1);endu01(na+2)=ut;u=u;ut;for i=1:nb-1u0(i)=u0(i+1);endu0(nb)=ut1; if t30 yyy1=0; yyy2=0; for i=1:lb yyy1=yyy1+B(i)*u(nb+1+t-i); end for i=2:la yyy2=yyy2+A(i)*y0(i); e
29、nd y1=yyy1-yyy2;else yyy1=0; yyy2=0; for i=1:lb yyy1=yyy1+B(i)*u(nb+1+t-i); end for i=2:la yyy2=yyy2+A(i)*y0(i)-0.1; end y1=yyy1-yyy2;endy=y;y1;for i=la:-1:2 y0(i)=y0(i-1);endy0(1)=y1; end tem=y;nu=length(y);for i=i:nu-1 re(i)=(y(i)+y(i+1)/2;endre(nu)=y(nu);by=1;for i=1:2:2*nu-1 qe(i)=y(by); qe(i+1)=re(by); by=by+1;end t=0.01:0.01:4;figure;plot(t,-qe+40);grid on
