数学建模论文-制动器试验台的控制方法分析3.pdf

1 制动器试验台的控制方法分析 摘 要 本文围绕制动器试验台的控制方法进行研究， 基于物理学刚体定轴转动的理想模型，利用 bp 神经网络进行拟合，建立了精确控制电流的计算机控制方法。 针对问题一，根据能量守恒定律以及物理量间的相互转化关系，求出等效的转动惯量为 51.992kg m。 针对问题二，根据机械惯量的定义求出各飞轮的机械惯量，分析题意可以得出八种组合方式，筛选出属于题中电动机补偿惯量范围内的两组数据：11.98462kg m，-18.02372kg m。 针对问题三，由附件提供的可观测的离散化瞬时转速与瞬时扭矩，在物理公式的基础上建立关于驱动电流的差分方程模型。结合问题一、二的结果计算出驱动电流值 i=174.9a。 针对问题四，根据附件数据求出计算机控制方式中每个可观测时间段的能量，得出试验台中能量的波动区间，再与路试时平均能量进行差值比较，以此建立以区间能量为指标的评价体系，对此体系进行劳斯判据检测，发现其稳定性较好。 针对问题五，将基于问题三的差分方程模型通过 z 变换，根据变换结果选用 logistic 模型以及 gompertz 模型进行拟合，利用最小二乘法对二者进行误差分析，最终选择 gompertz 模型作为改进的模型。 针对问题六，将得到的 gompertz 模型通过 bp 神经网络进行多次训练，拟合出精度更高的 gompertz 模型。 关键词：差分方程模型 最小二乘法 logistic 模型 gompertz 模型 bp 神经网络 2 1 问题重述 汽车行车制动器 （以下简称制动器） 的作用是在行驶时使汽车降速或者停止，在制动器的设计过程中路试是检验优劣的重要方式， 但在车辆设计阶段无法进行路试，只能在专门的制动器试验台上对其进行模拟实验。 模拟实验的原则与路试过程基本一致，即：当车辆加速到指定速度时，断开发动机的输出，同时以恒定的力施加制动，使车辆完全停止下来或车速降到某数值以下。在此过程中假设轮胎与地面的摩擦力为无穷大，即仅有静摩擦力，而不考虑滑动摩擦力。 路试车辆制动过程可简化考虑只有轮胎的平动，即不考虑轮胎自身的转动。在制动过程中车轮会承受一定的载荷， 而将此载荷在车辆平动时具有的能量等效地转化为实验台上飞轮和主轴等机构转动时的能量， 与此能量相应的转动惯量在本文中称为等效的转动惯量。 试验台上的主轴等不可拆卸机构的惯量称为基础惯量。飞轮的惯量加上基础惯量称为机械惯量，即试验台上的实验值。但在实际应用中一般不能精确地用机械惯量模拟实验， 等效的转动惯量与机械惯量之间存在一定的差值；为此，在制动过程中，用电动机提供由于机械惯量不足而缺少的能量，称为电动机补偿惯量。 本文中假设试验台采用电动机的驱动电流与其产生的扭矩成正比， 比例系数取为 1.5；且实验台工作时主轴的瞬时转速与瞬时扭矩是可观测的离散量。由于制动器性能的复杂性和可观测量的有限， 造成电动机驱动电流与时间之间产生误差。工程上实际常用计算机控制方法来控制制动过程，即把整个时间离散化为许多小的时间段，再根据前一时间段观测到的瞬时转速/或瞬时扭矩设计本时间段驱动电流的输出值，以此逐次进行完成制动。我们用离散的可观测量建立初步的电动机驱动电流的差分方程模型，并计算最初的驱动电流。 能量误差的大小是评价控制方法优劣的一个重要能量指标。 由附件中提供的数据拟合得出工程上实际常用计算机的控制方法， 并对电动机驱动电流的数学模型进行评价。 本文中能量误差是指所设计的路试时的制动器与实验台上制动器在制动过程中消耗的能量之差， 将等效的转动惯量与机械惯量之差的理论值与用此控制方法试验得到的数据比较，得出误差分析。 由电动机驱动电流的差分方程通过递推得到由前一时间段的观测值来控制本时间段电流值的计算机控制方法。绘制电动机驱动电流的差分方程的图像，与logistic 曲线及 gompertz 曲线的特点相似，由此拟合出 gompertz 函数，对该计算机控制方法进行了改进。 改进后的控制方法仍有不足之处， 因此重新设计出一个尽量完善的计算机控制方法，并作出评价。 2 条件的假设与符号的约定 2.1 条件的假设 1、假设路试时车轮与地面的摩擦力为无穷大，即仅有静摩擦力，而不考虑滑动摩擦力。 2、为了计算的方便，忽略车轮自身转动具有的能量，即将车辆载荷时平动时具有的能量等效地转化为制动器实验台上飞轮和主轴等转动时具有的能量。 3、假设飞轮质量是连续分布的， 且制动装置安装在车轮上对其质量分布不影响。 3 4、假设飞轮是均匀分布的刚体，我们可以把飞轮看成有一定厚度的圆筒。 5、假设试验台采用电动机的驱动电流与其产生的扭矩成正比。 2.2 符号的约定 w： 能量 ke ： 动能 i： 电流 m： 扭矩 m： 载荷质量 v： 汽车行驶时轮胎线速度 ： 汽车行驶时轮胎角速度 r： 汽车滚轮半径 ： 角加速度 ih ： 各飞轮的厚度()1,2,3i = ij ： 各飞轮的转动惯量 ij ： 电动机可能补偿的惯量()1,2,i =? it ： 第i个时间段的时间 i： 第i个时间段的飞轮角速度 in： 第i个时间段的飞轮转速 ： 计算机控制方法的控制效率 3 问题的分析 3.1 问题一 转动惯量是指刚体上所有质点的质量与质点到轴的垂直距离平方乘积的和，它只决定于刚体的形状、 质量分布和转轴的位置， 而同刚体绕轴的转动状态无关。而本题中，将车轮制动时承受的载荷等效为车轮的质量，且假设车轮是规则形状的均质刚体。根据题意，将载荷在车辆平动时具有的能量等效地转化为试验台上飞轮和主轴等机构转动时具有的能量，平动时车辆具有的能量为平动动能，在试验台上时具有的能量为转动动能。 根据能量的守恒原理及相关物理量间的转化关系，即可求得单个前轮等效的转动惯量。 3.2 问题二 所谓机械惯量就是指飞轮组的转动惯量之和再加上基础惯量。 假设飞轮是均匀分布的刚体，我们可以把飞轮看成有一定厚度的圆筒。根据圆筒的分布情况和惯量的定义，可以求出每个飞轮的惯量（详细分析见正文） ，由于 3 个飞轮可以4 都不固定在主轴上，可以任选一个固定在主轴上，可以任选二个在主轴上，还可以全放在主轴上。则可能有的组合为012333338ncccc=+=种，而电动机的补偿是因为机械惯量的不足，即将等效的转动惯量减去机械惯量。 3.3 问题三 实验台工作时可以测得瞬时扭矩及瞬时转速， 一般假设试验台采用的电动机的驱动电流与其产生的扭矩成正比， 那么要建立驱动电流的数学模型就要利用扭矩的物理公式mj=，实际工作中只能观测到离散的转速值，求出各个时间段的角加速度。由这些离散的数据可以建立相应驱动电流的差分方程模型。 利用所建的差分方程求解在问题1、2的条件下，求出初始速度150km h，末速度为0，制动时间为5s的驱动电流。 3.4 问题四 根据题目所给数据，求出测试过程中由电动机补偿的能量w以及计算机控制方法中每一小时间段的能量w。 由于整体的数据不能代表局部的情况 （例如：如果一个时间段的数据过大，而另一时间段的数据过小，但整体的数据不会有什么变化） ，但局部数据可以反映整体的情况，所以本题中评价不能是求和控制方法中的每一个时间段的能量（即4631iw）与测试的补偿能量（即w）比值，而是要比较每一个时间段能量的相差值。由于路试过程是一个实际测试过程，每一个阶段的能量变化比较合理，所以在时间相差很短的范围内，能量应该相差不大。故可以用平均能量w代表每一个阶段的能量。即比较w与w相差多少。因为对于本题来说：w有463个，所以取w中的最小值min和最大值max中与w相差较大的数（即max, minsmaxww=） ，与w比较（即sw=） ，即看的大小。 3.5 问题五 对问题三的差分方程模型进行z变换，将离散的数据连续化，根据变换结果引入logistic曲线以及gompertz曲线，通过拟合得出在此物理背景下的logistic模型以及gompertz模型。为了比较二者的优劣，采用最小二乘法对其进行误差分析，由误差分析结果选出更合理的模型。 3.6 问题六 基于问题五的gompertz模型的拟合初期效果不佳，为解决此不足，采用bp神经网络拟合模型，得出更为精确的gompertz模型。 5 4 模型的建立及求解 4.1 问题一 当汽车在行驶时，车轮的运动包括平动和转动两种。而本题中假设路试时车轮与地面的摩擦力为无穷大，即仅有静摩擦力，而不考虑滑动摩擦力。因此在制动时只考虑车轮的平动，而不考虑其转动。为了计算的方便，忽略车轮自身转动具有的能量， 即将车辆载荷时平动时具有的能量等效地转化为制动器实验台上飞轮和主轴等转动时具有的能量。 即能量w等于动能ke.根据动能的定义： 212kemv=（m为载荷质量，v为行驶时的线速度） (1) 在转动时可得飞轮的转动动能： 212kej=（j为转动惯量，为角速度） (2) 根据能量守恒定律：平动时具有的动能等于转动时具有的动能。即 kkee= (3) 根据线速度与角速度的关系： vr=（r为半径） (4) 由公式（1）（4）可得： 2jmr= (5) 由载荷质量6230635.71429.8gmg=(单位：kg),0.286r =（单位：m）代入上式得等效的转动惯量： 51.9989j =（单位：2kg mi） 4.2 问题二 钢制飞轮在绕固定轴转动时，假设飞轮的质量连续分布，则飞轮对转动轴的转动惯量为： 22i iijmrr dm= (6) 对于外径为1m，内径为0.2m，厚度为(1 2 3)ih i =、 ，密度为78102kg m的均匀分布的飞轮，我们可以把飞轮看成半径为r，宽度为dr，厚度为(1 2 3)ih i =、的圆筒，如图1（120.5,0.1rr=） ，并令这个筒的转动惯量为idj， 6 图1 则2iidjr dm= （7） 而2iiidmdvhrdr= （8） 即32iidjhr dr= （9） 于是，该飞轮对中心轴的转动惯量为 0.530.12iiijdjhr dr= （10） 当厚度10.0392h =m时，0.5310.11923.608130.0083jr dr=（单位：2kg m） ； 厚度20.0784h =m时，0.5320.13847.219460.0166jr dr=（单位：2kg m） ； 厚度30.1568h =m时，0.5330.17831.8396122.1767jr dr=（单位：2kg m） 。 则可以组成8种机械惯量， 分别为10，40.0083，70.0166，132.1767，110.0249，192.1933，100.0249，192.19332kg m. 由上题求出等效的转动惯量51.9989j =2kg m， 而电动机能补偿的能量相应的惯量的范围为2kg m30,30。计算出需要电动机补偿的惯量只有2种可能： 即151.992940.008311.9846j=（单位：2kg m） ； 251.992970.016618.0237j= （单位：2kg m） 。 其它的补偿惯量均超出电动机的补偿范围，在这里不再一一列出。 4.3问题三 由于实际中的计算机控制方法是把整个制动过程离散化为许多小的时间段，所以在试验台上得到主轴的瞬时转速和瞬时扭矩都是可观测的离散量， 而不是连续量。为了研究离散性变量的变化率，我们采用差分方程模型。 记第i个时间段的时间差为1iiittt+ =， （单位：s） ； 第i个时间段的角速度变化量为1iii+=， （单位：1s） 符号说明：电流( )i，扭矩()m，等效的转动惯量( )j，机械惯量()*j，角加速度( )，角速度( )，线速度( )v，旋转半径( )r，a为比例系数（本题中7 1.5a =）. 由*()()iammjjjjtvr = (11) 得到最终的电动机驱动电流的差分方程模型为 *()a jjvir t= (12) 在问题1和问题2的条件下，由初始速度1050vkm h=，末速度 10tvkm h=，5ts =，由本题所建的模型计算驱动电流为 *()1.5 (51.99929.993 10) 50174.90.286 5 3.6a jjviar t= 4.4 问题四 当飞轮在测试时所具有的始末转速in（1 2i= 、）分别为514转/分钟，257转/分钟，根据转速与角速度的关系得：2iin= (13) 将分钟转化为秒 112251453.82596060n=（单位：/rad s） ； 222221726.91306060n=（单位：/rad s） 。 由刚体的转动动能可得：初动能21112kej=； (14) 末动能22212kej=. (15) 由于路试等效的转动惯量为248kg m，机械惯量为235kg m，则由电动机提供的补偿惯量213jkg m=. 代入数据得： 118831.9788ke=，24708.0122ke=. 由定轴转动的动能定理得： 2114123.9666kkwee=（单位：j） (16) 由于飞轮测试的结果是比较理想的，所以每一个小时间段中能量是相差不大的。8 所以可以用平均的能量来近似表示测试时每个阶段的能量，即 14123.966630.5053( )463463wwwj= 则路试时汽车每个阶段的能量为30.5033j. 而在实际的计算机控制方法中是将整个制动时间简化为许多个离散的时间段，由附件中提供的数据可以得出，每个离散的小时间段等于0.01s,并且实际操作中观测到了各个时间点的瞬时转速与瞬时扭矩。 根据第1i+个时间点测得的瞬时转速和第i个时间点瞬时转速以及各物理量间的的关系，编程实现实际控制方法中各时间段的能量w（实现程序见附录1） 。取w的最大值与最小值，构成一个区间min,max，取max, minsmaxww=，则控制效率sw=. 对于得到的462个时间段的能量28.0932,33.3451w，则 33.3451 30.5053*100%9.4727%30.5053= 利用劳斯判据得出结论：“conclusion = 是 2”由此得出稳定性好， 评价体系合理。 （劳斯判据程序见附录2） 4.5问题五 根据问题3建立的关于驱动电流的差分方程模型*()a jjviajr t=. 关于变量的声明，无论从物理角度还是实际轮胎制动的过程中，都是关于时间的单调递减的一次函数，故将其拟合为一次函数（见图2） ，实现程序见附录3。 图2 角速度拟合 9 可以知道：第i个时间段的输出电流应为11iiiiiiiajtt+=，此处ij为随电流变化的不可控量，i由可观测量瞬时转速求得，0.01ts =，所以当第i个 时 间 点 的 输 出 电 流 已 知 ， 则 第1i+个 时 间 点 的 瞬 时 转 速 也 为 定 值11()iiiiiii ttaj+=+。由此可以做到根据前一个时间段观测到的瞬时转速或瞬时扭矩来控制本时间段的瞬时转速。 对差分方程模型进行z变换，并利用原始数据，用matlab绘制出角速度关于时间的图像，如图3 图3 并由该图发现与阻滞增长模型(logistic模型)类似， 在横坐标大于50以后总体呈振荡增长的趋势。在此基础上拟合出此物理背景下的logistic模型：( )1ktli tae=+，其中422.0829l=，9.2217a=，7.8874k=。 再考虑与logistic曲线非常相似的gompertz曲线，同样在本题的物理背景下得gompertz函数：( )ktbei tl e=i， 其中423.0078l=，3.097b=，5.5961k=。 将0( )i t、( )i t、( )i t三条曲线绘制在一个图像中，如图4 10 图4 从样本数据出发，用最小二乘法计算拟合误差，得到gompertz函数的拟合误差小于logistic函数。 由公式212kkkkiw=，取w的最大值与最小值，构成一个区间min,max，取max, minsmaxww=，则控制效率sw=. 对于得到的462个时间段的能量29.728,31.7968w，则 31.796830.5053*100%4.2337%30.5053= 由此可以按照此模型利用计算机对电流进行控制。 4.6 问题六 上文所建的gompertz的模型，由于制动器试验台在制动初期需要效应时间，故制动初期拟合效果并不理想。由此，通过bp神经网络进行拟合，对其多次训练。 （神经网络拟合的程序见附录4） 11 图5 训练误差 （其中实线为1b，长虚线为2b，短虚线为3b） 从训练效果看， 用其他训练方法即使很多步也难以达到levenberg-marquardt算法(b1)可以以较小的步就能得到满意的效果。用图像直观的体现经过多次训练后的拟合效果 图6 12 本程序也对隐层组合进行试验，采用tansig（）函数的训练的效果更好，所以本次拟合使用sigmoid函数； 并且从数据上分析增加隐层节点个数会使训练步数降低，但是神经网络泛化出现了问题；对于如何选择节点个数，只能根据实际情况用试凑方式选择。 最终得出更为精确的gompertz的模型： ( )5.559723.08913422.04781.7424tei tee= 同理由问题四建立的评价体系， 由公式212kkkkiw=，取w的最大值与最小值，构成一个区间min,max，取max, minsmaxww=，则控制效率sw=. 对于得到的462个时间段的能量30.4390,30.5908w，则 30.590830.5053*100%0.2802%30.5053= 5 模型的评价 5.1 模型的优点： 1、本文从简单粗糙的差分方程模型入手，通过拟合得到较为精确的gompertz模型，最后通过bp神经网络模型拟合出误差更小的模型。 2、文中能量误差的评价标准是让设计路试的每个阶段与实验所得各阶段的能量进行比较，基于物理知识，此标准较为合理； 3、 利用劳斯判据对本文所得评价体系进行稳定性分析， 认为其稳定性较好，该体系合理； 4、文中通过多次训练改进了bp神经网络，利用求出的局部极小点部分解决了bp神经网络的缺陷。 5.2 模型的缺点： 1、由于在制动器试验台模拟制动初期，飞轮减速过程不稳定，造成gompertz模型在初期的拟合效果不是很理想； 2、因为本题数据较少，只经过一次实验，使模型可能不适应更多不同条件下的实验数据。 参考文献 1马文蔚，物理学(第5版)m，北京，高等教育出版社，2006年4月。 2姜启源，谢金星，叶俊，数学模型（第三版）m,北京：高等教育出版社,2003年8月。 3薛定宇，陈阳泉，高等应用数学问题的matlab求解m,北京,清华大学出版社：2004年8月第1版。 13 4胡守仁，神经网络应用技术m，长沙，国防科技大学出版社，1995. 5陈建军， 制动器试验台机械惯量电模拟控制方法j， 起重运输机械，12，2007：27-30 6杨勇，基于神经元动态补偿的神经网络自适应控制及应用j，华北电力技术，2，2009：30-32 附录 附录1（问题四的程序） ： %a为导入的数据矩阵 n=a(:,2); %转速 m=a(:,1);%扭矩 for i=1:463 t(i)=(2*pi*n(i)./60;%角速度 end for i=1:462 w(i)=(t(i+1)-t(i)./0.01;%角加速度 end for i=1:462 if(w(i)=0) j(i)=0; else j(i)=m(i)./w(i);%转动惯量 end end t=t; w=w; j=j; for i=1:461 e(i)=0.5*(j(i+1)*t(i+1).2-j(i)*t(i).2); end e=e; w=e 附录2（劳斯判据） function routh_list,conclusion = routh(chara_equ) % = % 自编劳斯判据求解系统稳定性函数 % 输入： % chara_equ = 特征方程向量 % 输出： % routh_list = 劳斯表 14 % conclusion = 给出系统是否稳定或存在多少个不稳定的根的结论 % author:xianfa110 % example: % routh_list,con = routh(1 2 3 4 5); % return: % routh_list = % % 1 3 5 % 2 4 0 % 1 5 0 % -6 0 0 % 5 0 0 % con = % % there is 2 unstable roots! % = n=length(chara_equ); chara_equ=reshape(chara_equ,1,n); if mod(n,2)=0 n1=n/2; else n1=(n+1)/2; chara_equ=chara_equ,0; end routh=reshape(chara_equ,2,n1); routh_list=zeros(n,n1); routh_list(1:2,:)=routh; i=3; while 1; % =特殊情况1(第一列为0，其余列不为0)= if routh_list(i-1,1)=0 & sum(routh_list(i-1,2:n1)=0 chara_equ = conv(chara_equ,1 3); n=length(chara_equ); if mod(n,2)=0 n1=n/2; else n1=(n+1)/2; chara_equ=chara_equ,0; end routh=reshape(chara_equ,2,n1); 15 routh_list=zeros(n,n1); routh_list(1:2,:)=routh; i=3; end % =计算劳斯表= ai=routh_list(i-2,1)/routh_list(i-1,1); for j=1:n1-1 routh_list(i,j)=routh_list(i-2,j+1)-ai*routh_list(i-1,j+1); end % =特殊情况2(全0行)= if sum(routh_list(i,:)=0 k=0; l=1; f=zeros(1,n1); while n-i-k=0 f(l)=n-i+1-k; k=k+2; l=l+1; end routh_list(i,:)=routh_list(i-1,:).*f(1,:); end % =更新= i=i+1; if in break; end end % =outhput= r=find(routh_list(:,1)0); if isempty(r)=1 conclusion=the system is stable!; else n2=length(r); m=n2; for i=1:n2-1 if r(i+1)-r(i)=1 m=m-1; end end str1=there is ; 16 if r(n2)=n str2=num2str(m*2-1); else str2=num2str(m*2); end str3= unstable roots!; conclusion = str1,str2,str3; end 附录3(对的拟合程序) f=inline(a(1)*t+a(2),a,t); xx,res=lsqcurvefit(f,10 54,t,t) 结果为：xx=-6.0138 54.8998 plot(t,t) t1=-6.0138*t+54.8998; hold on plot(t,t1) logistic程序： f=inline(a(1)./(1+a(2).*exp(-a(3).*t),a,t); xx,res=lsqcurvefit(f,500 9 8,t,i) 结果：xx =422.0829 9.2217 7.8874 i2=422.0829./(1+9.2217*exp(-7.8874.*t); gompertz程序： f=inline(a(1).*exp(-a(2).*exp(-a(3).*t),a,t); xx,res=lsqcur