CN107246973A - 基于非线性滤波的抗蛇行减震器性能参数及故障辨识方法 - Google Patents

基于非线性滤波的抗蛇行减震器性能参数及故障辨识方法 Download PDF

Info

Publication number
CN107246973A
CN107246973A CN201710371939.2A CN201710371939A CN107246973A CN 107246973 A CN107246973 A CN 107246973A CN 201710371939 A CN201710371939 A CN 201710371939A CN 107246973 A CN107246973 A CN 107246973A
Authority
CN
China
Prior art keywords
mrow
msub
sampled point
parameter
suspension
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Granted
Application number
CN201710371939.2A
Other languages
English (en)
Other versions
CN107246973B (zh
Inventor
魏秀琨
张晓中
尹贤贤
魏德华
贾利民
李岩
刘玉鑫
王腾腾
闫冬
张靖林
江思阳
杨子明
李赛
孟鸿飞
赵利瑞
王熙楠
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Beijing Jiaotong University
Original Assignee
Beijing Jiaotong University
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by Beijing Jiaotong University filed Critical Beijing Jiaotong University
Priority to CN201710371939.2A priority Critical patent/CN107246973B/zh
Publication of CN107246973A publication Critical patent/CN107246973A/zh
Application granted granted Critical
Publication of CN107246973B publication Critical patent/CN107246973B/zh
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Classifications

    • GPHYSICS
    • G01MEASURING; TESTING
    • G01MTESTING STATIC OR DYNAMIC BALANCE OF MACHINES OR STRUCTURES; TESTING OF STRUCTURES OR APPARATUS, NOT OTHERWISE PROVIDED FOR
    • G01M17/00Testing of vehicles
    • G01M17/08Railway vehicles
    • G01M17/10Suspensions, axles or wheels

Abstract

本发明公开一种基于非线性滤波的抗蛇行减震器性能参数及故障辨识方法,包括如下步骤:S1、建立车辆悬挂系统横向动力学模型;S2、根据车辆悬挂系统横向动力学模型,建立悬挂系统的离散状态方程和离散观测方程;S3、根据悬挂系统的离散状态方程和离散观测方程,基于边缘化粒子滤波算法对悬挂系统的抗蛇行减震器进行性能参数辨识和故障辨识。本发明缩小了模型与实际运营的高速列车实际参数间的差异,提高了辨识结果的准确性。

Description

基于非线性滤波的抗蛇行减震器性能参数及故障辨识方法
技术领域
本发明涉及高速列车悬挂系统故障辨识领域。更具体地,涉及一种基于非线性滤波的抗蛇行减震器性能参数及故障辨识方法。
背景技术
铁路是国民经济发展的大动脉,是实现国家现代化建设的重要支撑。高速铁路作为一种安全、高速、舒适的交通运输方式,对国家经济社会的发展起到巨大的推动作用。2016年,调整后的《中长期铁路网规划》明确提出建设以“八纵八横”主通道为骨架、区域连接线衔接、城际铁路补充的高速铁路网。高速铁路运输在我国交通运输尤其是铁路运输中的地位变得越来越重要。安全性是高速铁路的核心竞争力,保障高速列车的安全运行不仅关系到乘客的生命安全,更会对我国高速铁路的发展和走出国门产生深远的影响。此外,随着我国高速铁路的快速发展,大量高速列车投入运营,由此产生的一系列维修保养问题亟待解决。
作为保障高速列车安全运行的关键设备,抗蛇行减震器能够有效抑制车辆的蛇行运动,对车辆运行的横向稳定性非常重要。现有的抗蛇行减震器性能参数及故障辨识方法仅通过间接分析车辆振动数据判断零部件的故障状态,容易受到采集设备、数据处理方法和耦合故障的影响,辨识效率较低。
因此,需要提供一种直接对抗蛇行减震器的性能参数进行辨识并基于性能参数辨识结果对抗蛇行减震器的故障进行辨识的基于非线性滤波的抗蛇行减震器性能参数及故障辨识方法。
发明内容
本发明的目的在于提供一种基于非线性滤波的抗蛇行减震器性能参数及故障辨识方法,解决现有的高速列车悬挂系统性能参数及故障辨识普遍采用半车模型进行算法研究,模型简化程度较高,与实际情况差别较大的问题,缩小模型与实际运营的高速列车间的差异,提高辨识结果的准确性。
为达到上述目的,本发明采用下述技术方案:
一种基于非线性滤波的抗蛇行减震器性能参数及故障辨识方法,包括如下步骤:
S1、建立车辆悬挂系统横向动力学模型;
S2、根据车辆悬挂系统横向动力学模型,建立悬挂系统的离散状态方程和离散观测方程;
S3、根据悬挂系统的离散状态方程和离散观测方程,基于边缘化粒子滤波算法对悬挂系统的抗蛇行减震器进行性能参数辨识和故障辨识。
优选地,步骤S2中,
悬挂系统的离散状态方程为:
xk+1=Axk+Buk+Ewk
其中,xk为第k个采样点的状态变量,uk为第k个采样点的系统输入,wk为第k个采样点的过程噪声,离散状态方程中悬挂系统矩阵 离散状态方程中过程噪声矩阵 为连续状态方程中悬挂系统矩阵,为连续状态方程中过程噪声矩阵;
悬挂系统的离散观测方程为:
yk=Cxk+Duk+Fvk
其中,离散状态方程中悬挂系统矩阵 为连续状态方程中悬挂系统矩阵,离散状态方程中观测噪声矩阵 为连续状态方程中观测噪声矩阵,vk为第k个采样点的观测噪声。
优选地,步骤S3的具体过程为:
S3.1、建立悬挂系统的非线性系统模型:
θk~p(θkk-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=[xk T θk T]T
由贝叶斯定理可得:
p(ξk|Zk)=p(xkk|Zk)=p(xkk,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);
设系统初始状态值为对应协方差矩阵为P0,对进行采样,获得系统状态的预测初值x1|0(i),并根据下式计算对应的协方差矩阵P1|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的偏差计算粒子权重并进行归一化,得到αk(i)
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)CTk|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)CTk(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;为第k个采样点的各参数粒子的蒙特卡洛均值, 为噪声,且有Vk为第k个采样点的参数粒子的方差矩阵,
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))Pk T(i)ATk+1|k(i))+E(θk+1|k(i))QwETk+1|k(i))
其中,Qw为系统过程噪声的协方差矩阵;
S3.4、在各采样点,根据参数辨识结果跟踪抗蛇行减震器的性能参数衰减,进而实现故障辨识。
优选地,步骤S3.4的具体方法为:设置故障阈值,将各采样点对应的性能参数值与故障阈值比较,若采样点对应的性能参数值小于等于故障阈值则判断抗蛇行减震器在该采样点的采样时刻出现故障,其中,设性能参数的正常值为θnormal,则故障阈值θfault设置为θfault=0.5θnormal
优选地,步骤S3.3还包括:将系统观测输出等间隔分成若干段,在完成基于第k`-1段观测数据的参数辨识后,以第k`-1段观测数据的参数辨识中最后一次参数辨识的结果作为参考值,在范围内均匀采样重新生成参数粒子,作为第k`段参数辨识过程的参数粒子初始值。
本发明的有益效果如下:
本发明所述技术方案将基于非线性滤波的参数辨识算法应用于高速列车抗蛇行减震器性能参数及故障的辨识中。通过建立高速列车悬挂系统的整车模型,克服了模型简化程度过高的问题,使得参数辨识过程更加贴近高速列车运行的实际过程。在此基础上,本发明将再次均匀采样策略引入边缘化粒子滤波算法,使参数粒子重新获得多样性,有效克服了算法收敛后的粒子贫化问题。通过对突发故障情况下的抗蛇行减震器性能参数进行辨识,可以实现性能参数的跟踪监测,从而为车辆维修部门提供运维保障的决策支持。
附图说明
下面结合附图对本发明的具体实施方式作进一步详细的说明;
图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=[ua1 ua2 ua3 ua4]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
其中,离散状态方程中悬挂系统矩阵 为连续状态方程中悬挂系统矩阵,离散状态方程中观测噪声矩阵 为连续状态方程中观测噪声矩阵,vk为第k个采样点的观测噪声。
步骤S3的具体过程为:
S3.1、建立悬挂系统的非线性系统模型:
θk~p(θkk-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=[xk Tθk T]T
由贝叶斯定理可得
p(ξk|Zk)=p(xkk|Zk)=p(xkk,Zk)p(θk|Zk)
其中,Zk为第k个采样点的系统观测值,p(xkk,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);
设系统初始状态值为对应协方差矩阵为P0,对进行采样,获得系统状态的预测初值x1|0(i),并根据下式计算对应的协方差矩阵P1|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的偏差计算粒子权重并进行归一化,得到αk(i)
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)CTk|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)CTk(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;为第k个采样点的各参数粒子的蒙特卡洛均值, 为噪声,且有Vk为第k个采样点的参数粒子的方差矩阵,
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段观测数据的参数辨识中最后一次参数辨识的结果作为参考值(即,以第k`-1组的参数辨识中最后一次参数辨识的结果作为参考值),在范围内均匀采样重新生成参数粒子,作为第k`段参数辨识过程的参数粒子初始值。
下面代入具体设置作进一步说明:
设置系统观测输出分为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段观测数据的参数辨识后,以参数辨识结果作为参考值,并在范围内均匀采样重新生成参数粒子,作为第k`段参数辨识过程的初始参数粒子。由图12可知,算法对突发故障的反应时间由两方面因素决定:一是算法本身的辨识速度,二是系统观测输出的分段长度。当算法辨识速度较快时,需要的迭代次数较少,因此可以将系统观测输出的分段长度设置得尽量短,使算法可以很快辨识出突变后的参数。
各阶段的具体情况如表二所示:
表二故障辨识各阶段的设置情况
具体说明如下:
(1)阶段2在第126个采样点处发生突发故障,参数由正常值θnormal跳变为故障阈值θfault
(2)阶段2开始于正常工况,结束于故障工况,如果辨识结果在参数突变前已经收敛,则粒子贫化现象将导致算法不能及时对突发故障做出反应,因此阶段2的参数辨识结果理论上应收敛于正常值θnormal
(3)阶段3均匀采样的参考值理论上应接近参数的正常值θnormal,并且如果算法有效,参数辨识结果理论上应该收敛于参数的故障值θfault
如图13所示,对抗蛇行减震器进行故障仿真,设置车体前侧的两个抗蛇行减震器发生突发故障,车体后侧的两个抗蛇行减震器处于正常状态。同时对故障和正常的抗蛇行减震器进行参数辨识,结果如图14所示。由图可知,对正常的抗蛇行减震器而言,参数辨识结果能够收敛在正常值附近;对故障的抗蛇行减震器而言,参数辨识结果能够跟踪参数跳变,并且辨识值能够收敛在真实值附近。综上所述,基于再次均匀采样策略的边缘化粒子滤波算法能够实现突发故障情况下性能参数的跟踪监测,从而实现故障辨识。
显然,本发明的上述实施例仅仅是为清楚地说明本发明所作的举例,而并非是对本发明的实施方式的限定,对于所属领域的普通技术人员来说,在上述说明的基础上还可以做出其它不同形式的变化或变动,这里无法对所有的实施方式予以穷举,凡是属于本发明的技术方案所引伸出的显而易见的变化或变动仍处于本发明的保护范围之列。

Claims (5)

1.一种基于非线性滤波的抗蛇行减震器性能参数及故障辨识方法,其特征在于,该方法包括如下步骤:
S1、建立车辆悬挂系统横向动力学模型;
S2、根据车辆悬挂系统横向动力学模型,建立悬挂系统的离散状态方程和离散观测方程;
S3、根据悬挂系统的离散状态方程和离散观测方程,基于边缘化粒子滤波算法对悬挂系统的抗蛇行减震器进行性能参数辨识和故障辨识。
2.根据权利要求1所述的基于非线性滤波的抗蛇行减震器性能参数及故障辨识方法,其特征在于,步骤S2中,
悬挂系统的离散状态方程为:
xk+1=Axk+Buk+Ewk
其中,xk为第k个采样点的状态变量,uk为第k个采样点的系统输入,wk为第k个采样点的过程噪声,离散状态方程中悬挂系统矩阵 离散状态方程中过程噪声矩阵 为连续状态方程中悬挂系统矩阵,为连续状态方程中过程噪声矩阵;
悬挂系统的离散观测方程为:
yk=Cxk+Duk+Fvk
其中,离散状态方程中悬挂系统矩阵 为连续状态方程中悬挂系统矩阵,离散状态方程中观测噪声矩阵 为连续状态方程中观测噪声矩阵,vk为第k个采样点的观测噪声。
3.根据权利要求2所述的基于非线性滤波的抗蛇行减震器性能参数及故障辨识方法,其特征在于,步骤S3的具体过程为:
S3.1、建立悬挂系统的非线性系统模型:
θk~p(θkk-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=[xk Tθk T]T
由贝叶斯定理可得:
p(ξk|Zk)=p(xkk|Zk)=p(xkk,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);
设系统初始状态值为对应协方差矩阵为P0,对进行采样,获得系统状态的预测初值x1|0(i),并根据下式计算对应的协方差矩阵P1|0(i):
<mfenced open = "{" close = ""> <mtable> <mtr> <mtd> <mrow> <msub> <mover> <mi>x</mi> <mo>&amp;OverBar;</mo> </mover> <mrow> <mn>1</mn> <mo>|</mo> <mn>0</mn> </mrow> </msub> <mo>=</mo> <msub> <mover> <mi>x</mi> <mo>^</mo> </mover> <mn>0</mn> </msub> </mrow> </mtd> </mtr> <mtr> <mtd> <mrow> <msub> <mi>P</mi> <mrow> <mn>1</mn> <mo>|</mo> <mn>0</mn> </mrow> </msub> <mrow> <mo>(</mo> <mi>i</mi> <mo>)</mo> </mrow> <mo>=</mo> <mi>E</mi> <mo>&amp;lsqb;</mo> <mrow> <mo>(</mo> <msub> <mi>x</mi> <mrow> <mn>1</mn> <mo>|</mo> <mn>0</mn> </mrow> </msub> <mo>(</mo> <mi>i</mi> <mo>)</mo> <mo>-</mo> <msub> <mover> <mi>x</mi> <mo>&amp;OverBar;</mo> </mover> <mrow> <mn>1</mn> <mo>|</mo> <mn>0</mn> </mrow> </msub> <mo>)</mo> </mrow> <msup> <mrow> <mo>(</mo> <msub> <mi>x</mi> <mrow> <mn>1</mn> <mo>|</mo> <mn>0</mn> </mrow> </msub> <mo>(</mo> <mi>i</mi> <mo>)</mo> <mo>-</mo> <msub> <mover> <mi>x</mi> <mo>&amp;OverBar;</mo> </mover> <mrow> <mn>1</mn> <mo>|</mo> <mn>0</mn> </mrow> </msub> <mo>)</mo> </mrow> <mi>T</mi> </msup> <mo>&amp;rsqb;</mo> </mrow> </mtd> </mtr> </mtable> </mfenced>
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的偏差计算粒子权重并进行归一化,得到αk(i)
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)CTk|k-1(i))+Qv
<mrow> <msub> <mover> <mi>&amp;alpha;</mi> <mo>~</mo> </mover> <mi>k</mi> </msub> <mrow> <mo>(</mo> <mi>i</mi> <mo>)</mo> </mrow> <mo>=</mo> <mi>p</mi> <mrow> <mo>(</mo> <msub> <mi>y</mi> <mi>k</mi> </msub> <mo>|</mo> <msub> <mi>Z</mi> <mrow> <mi>k</mi> <mo>-</mo> <mn>1</mn> </mrow> </msub> <mo>,</mo> <msub> <mi>&amp;theta;</mi> <mrow> <mi>k</mi> <mo>|</mo> <mi>k</mi> <mo>-</mo> <mn>1</mn> </mrow> </msub> <mo>(</mo> <mi>i</mi> <mo>)</mo> <mo>)</mo> </mrow> <mo>~</mo> <mi>N</mi> <mrow> <mo>(</mo> <msub> <mi>y</mi> <mrow> <mi>k</mi> <mo>|</mo> <mi>k</mi> <mo>-</mo> <mn>1</mn> </mrow> </msub> <mo>(</mo> <mi>i</mi> <mo>)</mo> <mo>,</mo> <msub> <mi>R</mi> <mi>k</mi> </msub> <mo>(</mo> <mi>i</mi> <mo>)</mo> <mo>)</mo> </mrow> </mrow>
<mrow> <msub> <mi>&amp;alpha;</mi> <mi>k</mi> </msub> <mrow> <mo>(</mo> <mi>i</mi> <mo>)</mo> </mrow> <mo>=</mo> <mfrac> <mrow> <msub> <mover> <mi>&amp;alpha;</mi> <mo>~</mo> </mover> <mi>k</mi> </msub> <mrow> <mo>(</mo> <mi>i</mi> <mo>)</mo> </mrow> </mrow> <mrow> <msubsup> <mi>&amp;Sigma;</mi> <mrow> <mi>j</mi> <mo>=</mo> <mn>1</mn> </mrow> <mi>N</mi> </msubsup> <msub> <mover> <mi>&amp;alpha;</mi> <mo>~</mo> </mover> <mi>k</mi> </msub> <mrow> <mo>(</mo> <mi>i</mi> <mo>)</mo> </mrow> </mrow> </mfrac> </mrow>
其中,Qv为系统观测噪声的协方差矩阵;
S3.3.2、参数辨识:
计算得到第k个采样点的参数辨识值为
<mrow> <msub> <mover> <mi>&amp;theta;</mi> <mo>^</mo> </mover> <mi>k</mi> </msub> <mo>=</mo> <munderover> <mo>&amp;Sigma;</mo> <mrow> <mi>i</mi> <mo>=</mo> <mn>1</mn> </mrow> <mi>N</mi> </munderover> <msub> <mi>&amp;alpha;</mi> <mi>k</mi> </msub> <mrow> <mo>(</mo> <mi>i</mi> <mo>)</mo> </mrow> <msub> <mi>&amp;theta;</mi> <mrow> <mi>k</mi> <mo>|</mo> <mi>k</mi> <mo>-</mo> <mn>1</mn> </mrow> </msub> <mrow> <mo>(</mo> <mi>i</mi> <mo>)</mo> </mrow> </mrow>
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):
<mrow> <msub> <mi>K</mi> <mi>k</mi> </msub> <mrow> <mo>(</mo> <mi>i</mi> <mo>)</mo> </mrow> <mo>=</mo> <msub> <mi>P</mi> <mrow> <mi>k</mi> <mo>|</mo> <mi>k</mi> <mo>-</mo> <mn>1</mn> </mrow> </msub> <mrow> <mo>(</mo> <mi>i</mi> <mo>)</mo> </mrow> <msup> <mi>C</mi> <mi>T</mi> </msup> <mrow> <mo>(</mo> <msub> <mi>&amp;theta;</mi> <mi>k</mi> </msub> <mrow> <mo>(</mo> <mi>i</mi> <mo>)</mo> </mrow> <mo>)</mo> </mrow> <msubsup> <mi>R</mi> <mi>k</mi> <mrow> <mo>-</mo> <mn>1</mn> </mrow> </msubsup> <mrow> <mo>(</mo> <mi>i</mi> <mo>)</mo> </mrow> </mrow>
Rk(i)=C(θk(i))Pk|k-1(i)CTk(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、粒子滤波预测更新:
粒子滤波预测更新过程应用参数演化中心平滑模型实现,具体为
<mrow> <msub> <mi>&amp;theta;</mi> <mrow> <mi>k</mi> <mo>|</mo> <mi>k</mi> <mo>+</mo> <mn>1</mn> </mrow> </msub> <mrow> <mo>(</mo> <mi>i</mi> <mo>)</mo> </mrow> <mo>=</mo> <msub> <mi>a&amp;theta;</mi> <mi>k</mi> </msub> <mrow> <mo>(</mo> <mi>i</mi> <mo>)</mo> </mrow> <mo>+</mo> <mrow> <mo>(</mo> <mn>1</mn> <mo>-</mo> <mi>a</mi> <mo>)</mo> </mrow> <msub> <mover> <mi>&amp;theta;</mi> <mo>&amp;OverBar;</mo> </mover> <mi>k</mi> </msub> <mo>+</mo> <msubsup> <mi>W</mi> <mi>k</mi> <mi>&amp;theta;</mi> </msubsup> </mrow>
其中,a=(3δ-1)/2δ;δ为折扣因子,取值范围为(0,1],一般情况下取值为0.95~0.99;为第k个采样点的各参数粒子的蒙特卡洛均值, 为噪声,且有Vk为第k个采样点的参数粒子的方差矩阵,
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
<mrow> <msub> <mi>P</mi> <mrow> <mi>k</mi> <mo>+</mo> <mn>1</mn> <mo>|</mo> <mi>k</mi> </mrow> </msub> <mrow> <mo>(</mo> <mi>i</mi> <mo>)</mo> </mrow> <mo>=</mo> <mi>A</mi> <mrow> <mo>(</mo> <msub> <mi>&amp;theta;</mi> <mrow> <mi>k</mi> <mo>+</mo> <mn>1</mn> <mo>|</mo> <mi>k</mi> </mrow> </msub> <mo>(</mo> <mi>i</mi> <mo>)</mo> <mo>)</mo> </mrow> <msubsup> <mi>P</mi> <mi>k</mi> <mi>T</mi> </msubsup> <mrow> <mo>(</mo> <mi>i</mi> <mo>)</mo> </mrow> <msup> <mi>A</mi> <mi>T</mi> </msup> <mrow> <mo>(</mo> <msub> <mi>&amp;theta;</mi> <mrow> <mi>k</mi> <mo>+</mo> <mn>1</mn> <mo>|</mo> <mi>k</mi> </mrow> </msub> <mo>(</mo> <mi>i</mi> <mo>)</mo> <mo>)</mo> </mrow> <mo>+</mo> <mi>E</mi> <mrow> <mo>(</mo> <msub> <mi>&amp;theta;</mi> <mrow> <mi>k</mi> <mo>+</mo> <mn>1</mn> <mo>|</mo> <mi>k</mi> </mrow> </msub> <mo>(</mo> <mi>i</mi> <mo>)</mo> <mo>)</mo> </mrow> <msub> <mi>Q</mi> <mi>w</mi> </msub> <msup> <mi>E</mi> <mi>T</mi> </msup> <mrow> <mo>(</mo> <msub> <mi>&amp;theta;</mi> <mrow> <mi>k</mi> <mo>+</mo> <mn>1</mn> <mo>|</mo> <mi>k</mi> </mrow> </msub> <mo>(</mo> <mi>i</mi> <mo>)</mo> <mo>)</mo> </mrow> </mrow>
其中,Qw为系统过程噪声的协方差矩阵;
S3.4、在各采样点,根据参数辨识结果跟踪抗蛇行减震器的性能参数衰减,进而实现故障辨识。
4.根据权利要求3所述的基于非线性滤波的抗蛇行减震器性能参数及故障辨识方法,其特征在于,步骤S3.4的具体方法为:设置故障阈值,将各采样点对应的性能参数值与故障阈值比较,若采样点对应的性能参数值小于等于故障阈值则判断抗蛇行减震器在该采样点的采样时刻出现故障,其中,设性能参数的正常值为θnormal,则故障阈值θfault设置为θfault=0.5θnormal
5.根据权利要求3所述的基于非线性滤波的抗蛇行减震器性能参数及故障辨识方法,其特征在于,步骤S3.3还包括:将系统观测输出等间隔分成若干段,在完成基于第k`-1段观测数据的参数辨识后,以第k`-1段观测数据的参数辨识中最后一次参数辨识的结果作为参考值,在范围内均匀采样重新生成参数粒子,作为第k`段参数辨识过程的参数粒子初始值。
CN201710371939.2A 2017-05-24 2017-05-24 基于非线性滤波的抗蛇行减震器性能参数及故障辨识方法 Active CN107246973B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN201710371939.2A CN107246973B (zh) 2017-05-24 2017-05-24 基于非线性滤波的抗蛇行减震器性能参数及故障辨识方法

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN201710371939.2A CN107246973B (zh) 2017-05-24 2017-05-24 基于非线性滤波的抗蛇行减震器性能参数及故障辨识方法

Publications (2)

Publication Number Publication Date
CN107246973A true CN107246973A (zh) 2017-10-13
CN107246973B CN107246973B (zh) 2019-06-14

Family

ID=60017637

Family Applications (1)

Application Number Title Priority Date Filing Date
CN201710371939.2A Active CN107246973B (zh) 2017-05-24 2017-05-24 基于非线性滤波的抗蛇行减震器性能参数及故障辨识方法

Country Status (1)

Country Link
CN (1) CN107246973B (zh)

Cited By (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN109472887A (zh) * 2018-11-23 2019-03-15 中国铁路总公司 车辆悬挂运行状态识别方法及装置
CN110823542A (zh) * 2019-11-06 2020-02-21 中车青岛四方机车车辆股份有限公司 减振器测试装置及减振器测试方法

Citations (5)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN102607867A (zh) * 2012-02-22 2012-07-25 北京交通大学 基于glrt列车悬挂系在途故障检测系统及检测方法
CN103310044A (zh) * 2013-05-27 2013-09-18 上海工程技术大学 基于改进粒子滤波算法的轨道车辆悬挂系统参数估计方法
CN103471865A (zh) * 2013-09-12 2013-12-25 北京交通大学 基于线性判别法的列车悬挂系统故障分离方法
CN104360679A (zh) * 2014-10-21 2015-02-18 南京航空航天大学 基于动态执行器的列车悬挂系统故障诊断与容错控制方法
CN105370783A (zh) * 2015-12-01 2016-03-02 唐山轨道客车有限责任公司 抗蛇形减振器、非动力转向架及轨道车辆

Patent Citations (5)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN102607867A (zh) * 2012-02-22 2012-07-25 北京交通大学 基于glrt列车悬挂系在途故障检测系统及检测方法
CN103310044A (zh) * 2013-05-27 2013-09-18 上海工程技术大学 基于改进粒子滤波算法的轨道车辆悬挂系统参数估计方法
CN103471865A (zh) * 2013-09-12 2013-12-25 北京交通大学 基于线性判别法的列车悬挂系统故障分离方法
CN104360679A (zh) * 2014-10-21 2015-02-18 南京航空航天大学 基于动态执行器的列车悬挂系统故障诊断与容错控制方法
CN105370783A (zh) * 2015-12-01 2016-03-02 唐山轨道客车有限责任公司 抗蛇形减振器、非动力转向架及轨道车辆

Cited By (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN109472887A (zh) * 2018-11-23 2019-03-15 中国铁路总公司 车辆悬挂运行状态识别方法及装置
CN110823542A (zh) * 2019-11-06 2020-02-21 中车青岛四方机车车辆股份有限公司 减振器测试装置及减振器测试方法
CN110823542B (zh) * 2019-11-06 2021-08-20 中车青岛四方机车车辆股份有限公司 减振器测试装置及减振器测试方法

Also Published As

Publication number Publication date
CN107246973B (zh) 2019-06-14

Similar Documents

Publication Publication Date Title
CN110334371A (zh) 一种基于有限元模型的车-桥耦合系统振动计算方法
CN110155101A (zh) 横向全主动控制减振系统及其中控制器的控制方法
CN108268711A (zh) 一种风车轨桥耦合模型及桥上抗风行车准则制定方法
Yang et al. An iterative interacting method for dynamic analysis of the maglev train–guideway/foundation–soil system
CN105973457A (zh) 动车组车载平稳性监控装置及方法
Li et al. Natural frequency of railway girder bridges under vehicle loads
Guo et al. Running safety analysis of a train on the Tsing Ma Bridge under turbulent winds
CN111859580B (zh) 一种铁路线路线型动态分析与设计方法
CN107246973B (zh) 基于非线性滤波的抗蛇行减震器性能参数及故障辨识方法
CN109766635A (zh) 一种机车机械部件状态感知传感器优化布局方法
Xiang et al. Dynamic interaction analysis of high-speed maglev train and guideway with a control loop failure
Shi et al. An efficient non-iterative hybrid method for analyzing train–rail–bridge interaction problems
Min et al. Gust wind effects on stability and ride quality of actively controlled maglev guideway systems
Abood et al. Railway carriage simulation model to study the influence of vertical secondary suspension stiffness on ride comfort of railway carbody
Ren et al. Analysis of vibration and frequency transmission of high speed EMU with flexible model
Guo et al. Coupling vibration analysis of high-speed maglev train-viaduct systems with control loop failure
Cheng et al. Dynamics analysis of high-speed railway vehicles excited by wind loads
Yau Wave passage effects on the seismic response of a maglev vehicle moving on multi-span guideway
Peng et al. Research on control strategy of magnetorheological semi-active suspension for high-speed train
Yıldız et al. Improving curving performance of a straddle type monorail vehicle by using semi-active devices
CN105117554B (zh) 高速轨道车辆一系垂向悬架最优阻尼比的设计方法
Zhu et al. Dynamic performance influences on Hopf bifurcation characteristics for vehicles
Wang et al. Dynamic performance of low-and medium-speed maglev train running on the turnout
Hasnan et al. Development of a simulation model of semi-active suspension for monorail
Wu et al. Train-Bridge Dynamic Behaviour of Long-Span Asymmetrical-Stiffness Cable-Stayed Bridge

Legal Events

Date Code Title Description
PB01 Publication
PB01 Publication
SE01 Entry into force of request for substantive examination
SE01 Entry into force of request for substantive examination
GR01 Patent grant
GR01 Patent grant