本发明涉及高速列车悬挂系统故障辨识领域。更具体地,涉及一种基于非线性滤波的抗蛇行减震器性能参数及故障辨识方法。
背景技术:
铁路是国民经济发展的大动脉,是实现国家现代化建设的重要支撑。高速铁路作为一种安全、高速、舒适的交通运输方式,对国家经济社会的发展起到巨大的推动作用。2016年,调整后的《中长期铁路网规划》明确提出建设以“八纵八横”主通道为骨架、区域连接线衔接、城际铁路补充的高速铁路网。高速铁路运输在我国交通运输尤其是铁路运输中的地位变得越来越重要。安全性是高速铁路的核心竞争力,保障高速列车的安全运行不仅关系到乘客的生命安全,更会对我国高速铁路的发展和走出国门产生深远的影响。此外,随着我国高速铁路的快速发展,大量高速列车投入运营,由此产生的一系列维修保养问题亟待解决。
作为保障高速列车安全运行的关键设备,抗蛇行减震器能够有效抑制车辆的蛇行运动,对车辆运行的横向稳定性非常重要。现有的抗蛇行减震器性能参数及故障辨识方法仅通过间接分析车辆振动数据判断零部件的故障状态,容易受到采集设备、数据处理方法和耦合故障的影响,辨识效率较低。
因此,需要提供一种直接对抗蛇行减震器的性能参数进行辨识并基于性能参数辨识结果对抗蛇行减震器的故障进行辨识的基于非线性滤波的抗蛇行减震器性能参数及故障辨识方法。
技术实现要素:
本发明的目的在于提供一种基于非线性滤波的抗蛇行减震器性能参数及故障辨识方法,解决现有的高速列车悬挂系统性能参数及故障辨识普遍采用半车模型进行算法研究,模型简化程度较高,与实际情况差别较大的问题,缩小模型与实际运营的高速列车间的差异,提高辨识结果的准确性。
为达到上述目的,本发明采用下述技术方案:
一种基于非线性滤波的抗蛇行减震器性能参数及故障辨识方法,包括如下步骤:
s1、建立车辆悬挂系统横向动力学模型;
s2、根据车辆悬挂系统横向动力学模型,建立悬挂系统的离散状态方程和离散观测方程;
s3、根据悬挂系统的离散状态方程和离散观测方程,基于边缘化粒子滤波算法对悬挂系统的抗蛇行减震器进行性能参数辨识和故障辨识。
优选地,步骤s2中,
悬挂系统的离散状态方程为:
xk+1=axk+buk+ewk
其中,xk为第k个采样点的状态变量,uk为第k个采样点的系统输入,wk为第k个采样点的过程噪声,离散状态方程中悬挂系统矩阵
悬挂系统的离散观测方程为:
yk=cxk+duk+fvk
其中,离散状态方程中悬挂系统矩阵
优选地,步骤s3的具体过程为:
s3.1、建立悬挂系统的非线性系统模型:
θk~p(θk|θk-1)
xk=f(xk-1,uk-1,wk-1)=a(θk-1)xk-1+b(θk-1)uk-1+e(θk-1)wk-1
yk=h(xk,uk,vk)=c(θk)xk+d(θk)uk+f(θk)vk
其中,θk为第k个采样点的悬挂系统的抗蛇行减震器性能参数,p(·|θk-1)表示已知θk-1时的概率密度函数,xk为第k个采样点的状态变量,yk为第k个采样点的观测变量,uk为第k个采样点的系统输入,wk和vk分别为第k个采样点的过程噪声和观测噪声,f(·)为非线性状态转移方程,h(·)为非线性观测方程;
系统待辨识变量ξk为:
ξk=[xktθkt]t
由贝叶斯定理可得:
p(ξk|zk)=p(xk,θk|zk)=p(xk|θk,zk)p(θk|zk)
其中,zk为第k个采样点的系统观测值;
s3.2、初始化第1个采样点的悬挂系统的抗蛇行减震器性能参数θ1|0(i):
如果已知初始概率密度函数为p(θ0|z0),对其进行采样可以得到作为参数粒子初始值的第1个采样点的悬挂系统的抗蛇行减震器性能参数θ1|0(i);如果没有先验知识作为依据获得p(θ0|z0),则通过在参数取值范围[θmin,θmax]内均匀采样获得θ1|0(i);
设系统初始状态值为
s3.3、依次对k=2,k=3,…,k=t时,第k个采样点的悬挂系统的抗蛇行减震器性能参数进行辨识,具体包括:
s3.3.1、权值更新及归一化:
根据第k个采样点的状态预测粒子xk|k-1(i)计算对应的观测值yk|k-1(i),根据yk|k-1(i)与实际观测值yk的偏差计算粒子权重
yk|k-1(i)=c(θk|k-1(i))xk|k-1(i)+d(θk|k-1(i))uk
rk(i)=c(θk|k-1(i))pk|k-1(i)ct(θk|k-1(i))+qv
其中,qv为系统观测噪声的协方差矩阵;
s3.3.2、参数辨识:
计算得到第k个采样点的参数辨识值为
s3.3.3、重采样:
对粒子{θk|k-1(i),xk|k-1(i),pk|k-1(i):i=1,2,…,n}进行重采样,得到对应的采样结果为{θk(j),xk|k-1(j),pk|k-1(j):j=1,2,…,n},满足pr{θk(j)=θk|k-1(i)}=αk(i);
s3.3.4、卡尔曼滤波测量更新
计算第k个采样点的卡尔曼滤波增益矩阵kk(i)、状态变量的辨识结果xk(i)和均方误差矩阵的辨识结果pk(i):
rk(i)=c(θk(i))pk|k-1(i)ct(θk(i))+qv
xk(i)=xk|k-1(i)+kk(i)(yk-c(θk(i))xk|k-1(i))
pk(i)=pk|k-1(i)-kk(i)c(θk(i))pk|k-1(i);
s3.3.5、粒子滤波预测更新:
粒子滤波预测更新过程应用参数演化中心平滑模型实现,具体为
其中,a=(3δ-1)/2δ;δ为折扣因子,取值范围为(0,1],一般情况下取值为0.95~0.99;
s3.3.6、卡尔曼滤波预测更新:
根据卡尔曼滤波预测更新的结果以及粒子滤波预测更新的结果计算下一时刻的状态预测值xk+1|k(i)以及均方误差矩阵的预测值pk+1|k(i):
xk+1|k(i)=a(θk+1|k(i))xk(i)+b(θk+1|k(i))uk
pk+1|k(i)=a(θk+1|k(i))pkt(i)at(θk+1|k(i))+e(θk+1|k(i))qwet(θk+1|k(i))
其中,qw为系统过程噪声的协方差矩阵;
s3.4、在各采样点,根据参数辨识结果跟踪抗蛇行减震器的性能参数衰减,进而实现故障辨识。
优选地,步骤s3.4的具体方法为:设置故障阈值,将各采样点对应的性能参数值与故障阈值比较,若采样点对应的性能参数值小于等于故障阈值则判断抗蛇行减震器在该采样点的采样时刻出现故障,其中,设性能参数的正常值为θnormal,则故障阈值θfault设置为θfault=0.5θnormal。
优选地,步骤s3.3还包括:将系统观测输出等间隔分成若干段,在完成基于第k`-1段观测数据的参数辨识后,以第k`-1段观测数据的参数辨识中最后一次参数辨识的结果
本发明的有益效果如下:
本发明所述技术方案将基于非线性滤波的参数辨识算法应用于高速列车抗蛇行减震器性能参数及故障的辨识中。通过建立高速列车悬挂系统的整车模型,克服了模型简化程度过高的问题,使得参数辨识过程更加贴近高速列车运行的实际过程。在此基础上,本发明将再次均匀采样策略引入边缘化粒子滤波算法,使参数粒子重新获得多样性,有效克服了算法收敛后的粒子贫化问题。通过对突发故障情况下的抗蛇行减震器性能参数进行辨识,可以实现性能参数的跟踪监测,从而为车辆维修部门提供运维保障的决策支持。
附图说明
下面结合附图对本发明的具体实施方式作进一步详细的说明;
图1示出基于非线性滤波的抗蛇行减震器性能参数及故障辨识方法的流程图。
图2示出高速列车悬挂系统结构示意图,其中,①为轴箱弹簧、②为一系垂向减震器、③为空气弹簧、④为抗蛇行减震器、⑤为二系横向减震器。
图3示出高速列车悬挂系统横向动力学模型俯视图。
图4示出高速列车悬挂系统横向动力学模型俯视图。
图5示出高速列车悬挂系统横向动力学模型各自由度运动的正方向示意图。
图6示出高速列车悬挂系统横向动力学模型的仿真模型示意图。
图7示出高速列车悬挂系统横向动力学模型零极点分布示意图。
图8示出高速列车悬挂系统横向动力学模型车体子系统的零极点分布示意图。
图9示出高速列车悬挂系统横向动力学模型转向架子系统的零极点分布示意图。
图10示出高速列车悬挂系统横向动力学模型转向架子系统的零极点分布示意图。
图11示出本发明中的边缘化粒子滤波算法流程流程图。
图12示出再次均匀采样策略原理示意图。
图13示出抗蛇行减震器故障设置图。
图14示出抗蛇行减震器阻尼突发故障情况下的性能参数辨识结果示意图。
具体实施方式
为了更清楚地说明本发明,下面结合优选实施例和附图对本发明做进一步的说明。附图中相似的部件以相同的附图标记进行表示。本领域技术人员应当理解,下面所具体描述的内容是说明性的而非限制性的,不应以此限制本发明的保护范围。
本发明公开的基于非线性滤波的抗蛇行减震器性能参数及故障辨识方法(以下将抗蛇行减震器性能参数简述为参数)采用基于非线性滤波的参数辨识算法和再次均匀采样策略,实现对抗蛇行减震器进行性能参数辨识,进而实现对抗蛇行减震器进行故障辨识。如图1所示,该方法包括如下步骤:
s1、建立车辆悬挂系统横向动力学模型;
s2、根据车辆悬挂系统横向动力学模型,建立悬挂系统的离散状态方程和离散观测方程;
s3、根据悬挂系统的离散状态方程和离散观测方程,基于边缘化粒子滤波算法对悬挂系统的抗蛇行减震器进行性能参数辨识和故障辨识。
其中,步骤s1的具体过程为:
如图2所示,车辆悬挂系统是指位于车体与转向架及转向架与轮对之间,起到支撑、缓冲和引导作用的一系列弹簧和阻尼元件,一般包括两系悬挂系统。一系悬挂系统位于转向架构架与轮对之间,起到缓冲轮轨冲击、支撑转向架和引导车辆运行的作用,主要包括轴箱、轴箱弹簧和垂向减震器等;二系悬挂系统位于车体与转向架构架之间,起到缓冲高频振动、支撑车体和引导车辆运行的作用,主要包括空气弹簧、横向减震器和抗蛇行减震器等。
本发明对二系悬挂系统抗蛇行减震器的故障辨识问题开展研究,首先要建立悬挂系统的整体模型。车辆系统是复杂的非线性多体系统,车体、转向架构架和轮对的刚度要比悬挂系统大很多,因此可以不考虑其弹性,并将车辆系统简化为多刚体系统。分别考虑车体和转向架构架的横移、摇头和侧滚运动以及轮对的横移和摇头运动,建立车辆悬挂系统横向动力学模型:
首先,建立车辆悬挂系统状态空间模型:
y=cx+du
其中,
状态变量x为
观测变量y为
其中,yc、φc和ρc分别表示车体的横向运动位移、摇头角和侧滚角,yb1、φb1和ρb1分别表示前转向架的横向运动位移、摇头角和侧滚角,yb2、φb2和ρb2分别表示后转向架的横向运动位移、摇头角和侧滚角,yw1和φw1分别表示前转向架前侧轮对的横向运动位移和摇头角,yw2和φw2分别表示前转向架后侧轮对的横向运动位移和摇头角,yw3和φw3分别表示后转向架前侧轮对的横向运动位移和摇头角,yw4和φw4分别表示后转向架后侧轮对的横向运动位移和摇头角,
u为轨道横向不平顺激励,作用于车辆的四组轮对,u=[ua1ua2ua3ua4]t,ua1表示前转向架前侧轮对的轨道横向不平顺激励,ua2表示前转向架后侧轮对的轨道横向不平顺激励,ua3表示后转向架前侧轮对的轨道横向不平顺激励,ua4表示后转向架后侧轮对的轨道横向不平顺激励;
a、b、c和d分别为悬挂系统矩阵,可由运动体各自由度运动的微分方程得到;
引入系统过程噪声w和观测噪声v,对车辆悬挂系统状态空间模型进行变换,得到车辆悬挂系统横向动力学模型:
y=cx+du+fv=cx+d(u+b-1ew)+fv-db-1ew
其中,e和f分别为过程噪声矩阵和观测噪声矩阵,u+b-1ew作为新的系统输入,实现过程噪声和观测噪声的引入。
图3和图4分别示出了悬挂系统横向动力学模型的俯视图和后视图。为方便建立模型,指定运动体各自由度运动的正方向,如图5所示。
如图6所示,可在simulink中搭建高速列车悬挂系统横向动力学模型的仿真模型,并对整车系统、车体子系统、前/后转向架构架子系统和各位轮对子系统的零极点分布进行分析。
步骤s2的具体过程为:
基于非线性滤波的参数辨识是在离散域进行的,因此在进行参数辨识前需要先对连续系统进行离散化。
悬挂系统的连续状态方程为
其中,x(t)为状态变量,y(t)为输出变量,u(t)为系统输入,w(t)和v(t)分别为过程噪声和观测噪声,
连续系统离散化主要是对描述系统动态特性的状态方程而言,输出方程为静态的代数方程,离散化后应保持不变,即
基于连续系统状态方程的解,可得到悬挂系统的离散状态方程为
xk+1=axk+buk+ewk
其中,xk为第k个采样点的状态变量,uk为第k个采样点的系统输入,wk为第k个采样点的过程噪声,离散状态方程中悬挂系统矩阵
悬挂系统的离散观测方程为:
yk=cxk+duk+fvk
其中,离散状态方程中悬挂系统矩阵
步骤s3的具体过程为:
s3.1、建立悬挂系统的非线性系统模型:
θk~p(θk|θk-1)
xk=f(xk-1,uk-1,wk-1)=a(θk-1)xk-1+b(θk-1)uk-1+e(θk-1)wk-1
yk=h(xk,uk,vk)=c(θk)xk+d(θk)uk+f(θk)vk
其中,θk为第k个采样点的悬挂系统的抗蛇行减震器性能参数,p(·|θk-1)表示已知θk-1时的概率密度函数,xk为第k个采样点的状态变量,yk为第k个采样点的观测变量,uk为第k个采样点的系统输入,wk和vk分别为第k个采样点的过程噪声和观测噪声,f(·)为非线性状态转移方程,h(·)为非线性观测方程。
将系统待辨识变量ξk分为两部分:
ξk=[xktθkt]t
由贝叶斯定理可得
p(ξk|zk)=p(xk,θk|zk)=p(xk|θk,zk)p(θk|zk)
其中,zk为第k个采样点的系统观测值,p(xk|θk,zk)易于求出解析解,可通过卡尔曼滤波器实现。而p(θk|zk)无法直接求得解析解,本发明中通过粒子滤波器实现。则系统待辨识变量ξk的最小均方误差估计为
s3.2、初始化第1个采样点的悬挂系统的抗蛇行减震器性能参数θ1|0(i):
如果已知初始概率密度函数为p(θ0|z0),对其进行采样可以得到作为参数粒子初始值的第1个采样点的悬挂系统的抗蛇行减震器性能参数θ1|0(i);如果没有先验知识作为依据获得p(θ0|z0),则通过在参数取值范围[θmin,θmax]内均匀采样获得θ1|0(i);
设系统初始状态值为
s3.3、依次对k=2,k=3,…,k=t时,第k个采样点的悬挂系统的抗蛇行减震器性能参数进行辨识,具体包括:
s3.3.1、权值更新及归一化:
根据第k个采样点的状态预测粒子xk|k-1(i)计算对应的观测值yk|k-1(i),根据yk|k-1(i)与实际观测值yk的偏差计算粒子权重
yk|k-1(i)=c(θk|k-1(i))xk|k-1(i)+d(θk|k-1(i))uk
rk(i)=c(θk|k-1(i))pk|k-1(i)ct(θk|k-1(i))+qv
其中,qv为系统观测噪声的协方差矩阵;
s3.3.2、参数辨识:
计算得到第k个采样点的参数辨识值为
s3.3.3、重采样(粒子滤波测量更新):
对粒子{θk|k-1(i),xk|k-1(i),pk|k-1(i):i=1,2,…,n}进行重采样,得到对应的采样结果为{θk(j),xk|k-1(j),pk|k-1(j):j=1,2,…,n},满足pr{θk(j)=θk|k-1(i)}=αk(i);
s3.3.4、卡尔曼滤波测量更新
计算第k个采样点的卡尔曼滤波增益矩阵kk(i)、状态变量的辨识结果xk(i)和均方误差矩阵的辨识结果pk(i):
rk(i)=c(θk(i))pk|k-1(i)ct(θk(i))+qv
xk(i)=xk|k-1(i)+kk(i)(yk-c(θk(i))xk|k-1(i))
pk(i)=pk|k-1(i)-kk(i)c(θk(i))pk|k-1(i);
s3.3.5、粒子滤波预测更新:
粒子滤波预测更新过程应用参数演化中心平滑模型实现,具体为
其中,a=(3δ-1)/2δ;δ为折扣因子,取值范围为(0,1],一般情况下取值为0.95~0.99;
s3.3.6、卡尔曼滤波预测更新:
根据卡尔曼滤波预测更新的结果以及粒子滤波预测更新的结果计算下一时刻的状态预测值xk+1|k(i)以及均方误差矩阵的预测值pk+1|k(i):
xk+1|k(i)=a(θk+1|k(i))xk(i)+b(θk+1|k(i))uk
其中,qw为系统过程噪声的协方差矩阵。
s3.4、在各采样点,根据参数辨识结果跟踪抗蛇行减震器的性能参数衰减,进而实现故障辨识,具体方法为:设置故障阈值,将各采样点对应的性能参数值与故障阈值比较,若采样点对应的性能参数值小于等于故障阈值则判断抗蛇行减震器在该采样点的采样时刻出现故障,其中,设性能参数的正常值为θnormal,则故障阈值θfault可设置为θfault=0.5θnormal。
本发明采用的非线性滤波算法为边缘化粒子滤波算法,将状态辨识过程和参数辨识过程进行有效分离,分别采用卡尔曼滤波和粒子滤波进行辨识。当边缘化粒子滤波算法的参数辨识结果收敛后,所有参数粒子的权重基本相同,并且粒子值接近参数的正常值。此时即使抗蛇行减震器发生突发故障,由于参数粒子已经贫化,无法继续逼近参数的故障值,导致算法不能及时对突发故障做出反应。因此,实现故障辨识的关键是使参数粒子重新获得多样性,并在系统观测输出的修正下持续逼近故障值。因此本发明引入再次均匀采样策略,将参数辨识过程分割为若干段,每一段参数辨识过程结束后,对参数粒子进行均匀采样,使参数粒子重新获得多样性,以解决粒子贫化问题,步骤s3.3还包括:将系统观测输出等间隔分成若干段(即,将总采样点分为多组),在完成基于第k`-1段观测数据的参数辨识后(即,在完成第k`-1组的参数辨识后),以第k`-1段观测数据的参数辨识中最后一次参数辨识的结果
下面代入具体设置作进一步说明:
设置系统观测输出分为4段,分别为阶段1、阶段2、阶段3和阶段4,系统观测输出的分段长度t=250,粒子数ns=1000,采样间隔ts=0.5ms,被辨识参数的正常值为θnormal,故障值为θfault=0.5θnormal,均匀采样的参考值为
在simulink中搭建悬挂系统横向动力学模型的仿真模型后,设置车辆运行速度为80m/s,系统采样频率为2000hz,系统采样间隔为0.5ms。对悬挂系统横向动力学模型的零极点分布进行分析,整车系统、车体子系统、前/后转向架构架子系统和各位轮对子系统的零极点分布分别如图7至图10所示。如果各级系统的零极点均位于虚轴左侧,说明系统是稳定的。
本发明采用的非线性滤波算法是边缘化粒子滤波算法,该算法的基本流程如图11所示。在此基础上,引入在此均匀采样策略,基本原理如图12所示。
再次均匀采样策略的基本过程是:将系统观测输出等间隔分成若干段,在完成基于第k`-1段观测数据的参数辨识后,以参数辨识结果
各阶段的具体情况如表二所示:
表二故障辨识各阶段的设置情况
具体说明如下:
(1)阶段2在第126个采样点处发生突发故障,参数由正常值θnormal跳变为故障阈值θfault;
(2)阶段2开始于正常工况,结束于故障工况,如果辨识结果在参数突变前已经收敛,则粒子贫化现象将导致算法不能及时对突发故障做出反应,因此阶段2的参数辨识结果理论上应收敛于正常值θnormal;
(3)阶段3均匀采样的参考值
如图13所示,对抗蛇行减震器进行故障仿真,设置车体前侧的两个抗蛇行减震器发生突发故障,车体后侧的两个抗蛇行减震器处于正常状态。同时对故障和正常的抗蛇行减震器进行参数辨识,结果如图14所示。由图可知,对正常的抗蛇行减震器而言,参数辨识结果能够收敛在正常值附近;对故障的抗蛇行减震器而言,参数辨识结果能够跟踪参数跳变,并且辨识值能够收敛在真实值附近。综上所述,基于再次均匀采样策略的边缘化粒子滤波算法能够实现突发故障情况下性能参数的跟踪监测,从而实现故障辨识。
显然,本发明的上述实施例仅仅是为清楚地说明本发明所作的举例,而并非是对本发明的实施方式的限定,对于所属领域的普通技术人员来说,在上述说明的基础上还可以做出其它不同形式的变化或变动,这里无法对所有的实施方式予以穷举,凡是属于本发明的技术方案所引伸出的显而易见的变化或变动仍处于本发明的保护范围之列。