列车制动系统辨识和自适应控制

1、1、系统辨识部分1.1 列车制动系统介绍制动系统是列车操纵系统的重要组成部分，它是用来调节车速和进行停车，是列车安全可靠运行的基本保障。现代列车通常用电空混合制动系统，它通过制动控制器协调制动和空气制动的介入比例，基本能够保证在车速变化情况下一致的制动性能20当司机操纵控制手柄时，制动控制器收到司机请求的目标加速度，它根据车重和车速进行相应补偿，然后计算出所需制动力和分配方案，最后驱动执行机构完成制动力的输出。这是一个反馈控制的过程，但是由于它并不补偿由于坡道和弯道所带来的附加影响。1.2 收集先验知识ATO(自动驾驶系统)的功能是代替司机驾驶，它根据目标速度、行车许可及线路情况等自动产生牵引

2、或制动指令。ATO控制器和司机对车辆的控制都是通过车辆本身的牵引控制系统和制动控制系统完成的。对于ATO来说制动模型是包括列车制动控制系统在内的车辆运行模型，它描述从输入到制动控制器的制动指令到车辆运行状态：加速度、速度间的动态关系。图1是制动模型框图12o图1制动模型框图1.3 数学模型分析3制动控制器通过反馈调节实现对目标加速度的跟踪，这是一个动态的过程。这个过程可以用一阶动态系统近似，另外考虑到系统传输延时，则可以用下面方程进行描述U11.a(t)=-a(t)+At。)(1)TTA式中，a(t)为控制加速度，它是由制动系统控制器的作用而是列车产生的加速度；A为目标加速度；丁为系统响应时间

3、常数；仃为传输延迟。车辆实际加速度a(t)由控制加速度a、t)和环境(如弯道、坡道)造成的附加加速度|_a：t)构成，车辆速度v由实际速度决定。(2)a(t)=a(t)La(t)Lv(t)=a(t)(3)目标加速度A(t)为列车制动控制器的输入指令，它通过驾驶员或ATO控制(4)u(t)(制动指令)产生，它们之间的关心可以用式(4)静态函数关系描述。At)=F(u(t)以上方程描述了从ATO控车角度看到的列车制动模型，它可以图2框图表示。a(t)图2制动模型框图其中输入为驾驶员或ATO控制指令u(t),输出为车辆实时速度v(t),式(1)中e和仃是需要辨识的参数。辨识式(1)中参数需要知道控制

4、加速度3t)，但是由于|_a(t)的干扰，无法从测量的实时速度中获得。我们假设列车是在平直轨道上运行，环境影响可忽略不计的理想状态，即：Ua(t)-0(5)AL那么a(t)定V(t)(6)1.4稳态响应参数估计6暂态响应消失后a(t)=A(to)=F(u(t-cr)ta仃(7)可以看出，如果u(t)为一恒值，则A(t)也为恒值。根据式(5),可以利用测量的速度来计算在各制动级别上的稳态加速度，得到如下关系：Lv(t)=A(t-0)=F(u(t-二)=F(u)t>仃i=1,2,3,，N(8)实验测量数据为:U(t)12345A(t)-20.7-39.5-58.4-77.396.167-11

5、5.0-133.8用matlab对其进行线性拟合，代码为：clearu=1,2,3,4,5,6,7;A=-2O7,-39.5,-58.4,-77.3,-96.1,-115Q-133.8;fori=1:7;Z(i)=A(i)；endZL=Z'HL=u(1),1;u(2),1;u(3),1;u(4),1;u(5),1;u(6),1;u(7),1c1=HL'*HL;c2=inv(c1);c3=HL'*ZL;c4=c2*c3a=c4(1)b=c4(2)程序运行结果：ZL=-20.7000-39.5000-58.4000-77.3000-96.1000-115.0000-133.

6、8000HL=c4 =-18.8571-1.8286-18.8571-1.8286即可得1111A t ) = 1 8 . 8i6t ( )1(9)1.5动态响应参数辨识与仿真A为了方便计算，我们将A(t)作为输入，a(t)作为输出。式(1)描述了从目标加速度到控制加速度以7和仃为参数的传递函数，表示为：Gs；i,仃)，s为拉普拉斯变换变量。最优的参数应该是使输入u得到响应与实际响应误差e(t)最小。其中:(10)(11)e-s0s,Le(t)=g(t；,二)LF(u(t)-v(t)式中g(s；工,*是0；工,。)的反变换。误差大小用式(11)进行评估13。J(.,二)=e2(t )dt(12

7、)取邙=1.2s,0.4s,贝ij-1.2sGs,二)二(13)e-0.4s1对式(13)进行Z变换，取采样时间为T=1,得到离散系统的传递函数:Gz)2 0.8647 z 0.05325z -z 一 0.08208(14)则有y(k)0.08208y(k1)=0.8647u(k-2)+0.05325Uk-3)+v(k)(15)式中，v(k)是方差为0.1的白噪声干扰。制动控制系统为时滞系统，数据不能突变，将输入信号定为阶跃信号，用递推最小二乘法进行参数估计。算法步骤：1 .设置3(0)和P(0),输入初始数据；2 .采样当前输出y(k)和输入u(k);3 .利用最小二乘法估计$递推公式，计算

8、K(k),Q(k)和P(k);4 ,kTk+1返回第二步，继续循环4。仿真程序如下：clearall;closeall;a=1-0.08208;b=0.86470.05325'd=2;%寸象参数na=length(a)-1;nb=length(b)-1;%na、nb为A、B的阶次L=500;%仿真长度y1=1;y2=1;y3=1;y4=1;%生成M序列fori=1:L;%生成M序列x1=xor(y3,y4);x2=y1;x3=y2;x4=y3;y(i)=y4;ify(i)>0,u(i)=y(i);elseu(i)=0.5;endy1=x1;y2=x2;y3=x3;y4=x4;en

9、duk=zeros(d+nb,1);%输入初值：uk（i）表示u（k-i）yk=zeros(na,1);%输出初值xi=sqrt(0.1)*randn(L,1);%白噪声干扰theta=a(2:na+1);b;%对象参数真值%thetae 初值thetae_1=zeros(na+nb+1,1);P=10A6*eye(na+nb+1);fork=1:Lphi=-yk;uk(d:d+nb);%此处phi 为列向量y(k)=phi'*theta+xi(k);%采集输出数据%递推最小二乘法K=P*phi/(1+phi'*P*phi);thetae(:,k)=thetae_1+K*(y(

10、k)-phi'*thetae_1);P=(eye(na+nb+1)-K*phi')*P;%更新数据thetae_1=thetae(:,k);fori=d+nb:-1:2uk(i)=uk(i-1);enduk(1)=u(k);fori=na:-1:2yk(i)=yk(i-1);endyk(1)=y(k);endplot(1:L,thetae);xlabel('k');ylabel('2?ey1a?a?0b');); axis(0 L -2 4);legend('a_1','b_0','b_1'对系统进

11、行RIS辨识仿真,参数辨识如下:3210q盘较期O50O400025050图3白色噪声方差为0.01时参数辨识结果图4白色噪声方差为0.1时参数辨识结果q百工±=若卷2J,，1J.1.1.1.LL.050100150200250300350400450500k图5白色噪声方差为1时参数辨识结果以上参数曲线说明随着扰动强度的增强，辨识的效果会有所下降。2、自适应控制部分2.1基于最小方差间接自校正控制基于前节讨论的列车制动模型，采用最小方差间接自校正控制。当被控对象的参数未知时，可首先利用递推增广最小二乘法在线实时估计对象参数，然后再设计最小方差控制律，即将对象参数估计器和控制器的设计

12、分开进行，就形成了最小方差自校正控制间接算法网4。2.3控制器设计根据前文，现对最小方差自适应控制算法进行推导。1)、单步输出预测45：结合本系统，令输出为y(k),输入为u(k)。设系统数学模型为：-4,-d_1一，一、A(z)y(k)=zB(z)u(k)+C(z卜(k)(16)式中，C(z')为Hurwitz多项式，即其零点完全位于z平面的单位圆内；u(k)和y(k)表示系统的输入输出，”k)为方差仃2的白噪声，d>1是系统纯延时，且有Ia(z,)=1,a1z,-a2z2,，anzzJaJb(z,)=b0+b1z，+b2z/+.+bnbz"b,b0#0(17)C(z

13、,)=1CiZJ-az，c%z"c对象(2.1)基于k时刻和之前时刻的输入输出数据记作Yk,Uk=y(k),y(k1),.,u(k),u(k1),.(18)基于Wk,Uk对k+d时刻的输入记作?(k+d|k)(19)(输出)预测误差记作?(k+d|k)=y(k+d)-?(k+d|k)(20)那么，关于提前d步最小方差预测输出可由如下定理给出。2)、定理(最优d步预测输出)45使如下性能指标(预测误差的方差)Ey2(k+d|k)为最小的d步最优预测输出y*(k+d|k)必满足方程*1C(z)y(kd|k)=G(z)y(k)F(z)u(k)式中1_d1C(z)=A(z)E(z)z&quo

14、t;G(z)_1F(z)=B(z)E(z)且Rz,)=1+ez，+.+,zf,(ne=d1)G(z')=g0gz'.gnz,(ng=na-1)gF(z)=f°f1z.z,(%=踊d-1)此时，最有预测误差的方差为d1Ey*(kd|k)2=(1八e2)-23)、最小方差控制律45(21)(22)(23)(24)i1设B(z,)也是Hurwitz多项式，则其零点完全位于z平面单位圆内，即对象是最小相位或者逆稳，可得如下定理：设控制目标是使实际输出y(k+d)跟踪期望输出yr(k+d),使性能指标J=Ey(k+d)-y.(k+d)2(25)为最小。由上面定理可知：.J-*-

15、y(k+d)=E«k+d)+y(k+d|k)(26)将上式代入式(25)中，得2.J二或E(kd)y(kd|k)yr(kd)=EE(k+d)2+Ey*(k+d|k)-yr(k+d)2(27)上式右边第一项不可控，若使J最小，需使上式右边第二项为0,即*一_-,_.y(k+d|k)=yr(k+d)(28)将式(28)代入最优预测输出方程式(21)得Qz)y*(k+d|k)=Gz,)y(k)+F(z)u(k)(29)得Ctz)yr(k+d)=G(z,)y(k)+F(z)Uk)(30)上式经过简单的变形即得最小方差控制律为F(z')Uk)=az)y.(k+d)G(z,)y(k)(3

16、1)2.4仿真验证用matlab对上面自适应控制进行仿真，其步骤如下5：已知：模型阶次na，几，几及纯延时do第一步：设置初值日(0)和R0),输入初始数据；第二步：采样当前实际输出y(k)和期望输出yr(k+d)；第三步：利用递推增广最小二乘法在线估计被控对象参数即A,B和C;第四步：求解Sindiophantine方程，得到多项式E、F、G的系数；第五步：利用式(22)计算并实施u(k)；第六步：返回第二步(ktk+1),继续循环。被控对象为y(k)一0.08208y(k-1)=0.8647u(k-2)0.05325Uk-3)(k)式中qk)为白噪声。6取初值R0)=10I、40)=0.0

17、005,期望输出yr(k)为幅值为5的万波信号，采用最小方差自校正控制间接算法。由于参数估计初期，所获得的估计参数不准确，很可能偏离实际参数真值较远，使多项式B(z,)不稳定，可以在实施自校正控制时，将B(z)的参数估计值限制在稳定范围内6。仿真程序如下：clearall;closeall;clc;%模型参数%输出阶数a=1-0.08208;b=10.5;c=1;d=2;na=length(a)-1;nb=length(b)-1;nc=length(c)-1;nf=nb+d-1;%nf为多项式的阶次L=400;%控制步数uk=0.91*ones(d+nb,1);脸入初值；u(k)表示u(k-1

18、)yk=4.8*ones(na,1);%出初值yrk=1*ones(nc,1);%11望输出初值xik=zeros(nc,1);%白噪声初值xiek=zeros(nc,1);yr=5*0.8*ones(L/2,1);-0.6*ones(L/2+d,1);%期望输出xi=sqrt(0.1)*randn(L,1);硒噪声序列%RELS初值thetae_1=0.0005*ones(na+nb+1+nc,1);%lt匕处不能为0P=10A6*eye(na+nb+1+nc)；fork=1:Ltime(k)=k;%采集输出数据y(k)=-a(2:na+1)*yk+b*uk(d:d+nb)+c*xi(k);

19、%递推增广最小二乘法phie=-yk;uk(d:d+nb);K=P*phie/(1+phie'*P*phie);thetae(:,k)=thetae_1+K*(y(k)-phie'*thetae_1);P=(eye(na+nb+1+nc)-K*phie')*P;xie=y(k)-phie'*thetae(:,k);%白噪声的估值%提取辨识参数ae=1thetae(1:na,k)'be=thetae(na+1:na+nb+1,k)'ce=1e,f,g=sindiophantine(ae,be,ce,d);%求解单步Diophantine方程u(k)

20、=(-f(2:nf+1)*uk(1:nf)+ce*yr(k+d:-1:k+d-min(d,nc);yrk(1:nc-d)-g*y(k);yk(1:na-1)/f(1);%更新数据thetae_1=thetae(:,k);fori=d+nb:-1:2uk(i)=uk(i-1);enduk(1)=u(k);fori=na:-1:2yk(i)=yk(i-1);endyk(1)=y(k);endfigure(1);subplot(2,1,1);plot(time,yr(1:L),'r:',time,y);xlabel('k');ylabel('y_r(k)?0y(k)');legend('y_r(k)','y(k)');axis(0L-55)