方 宇，李曉斌，金子博，高 偉

(上海工程技術(shù)大學(xué) 城市軌道交通學(xué)院，上海 201620)

懸掛系統(tǒng)的狀態(tài)直接影響軌道車輛運(yùn)行的安全性、平穩(wěn)性和舒適性。據(jù)統(tǒng)計(jì)，在軌道車輛產(chǎn)生的眾多故障中，轉(zhuǎn)向架懸掛系統(tǒng)故障頻率最高，占總數(shù)的1/3左右[1]，是車輛故障中比例最大的一部分。目前懸掛系統(tǒng)在線監(jiān)測(cè)主要依賴于信號(hào)分析方法，需在軌道車輛上布置數(shù)目較多的傳感器，主要遵循研究關(guān)鍵部件檢測(cè)信號(hào)特性的突變，包括幅值、頻率特性以及對(duì)動(dòng)力學(xué)參數(shù)的統(tǒng)計(jì)特征進(jìn)行分析，檢測(cè)結(jié)果可靠性較低，具有一定的局限性[2]。

參數(shù)估計(jì)是近年來提出的一種全新的車輛狀態(tài)監(jiān)測(cè)方法，在實(shí)際應(yīng)用時(shí)只需利用少量的傳感器獲取車輛振動(dòng)狀態(tài)信息，就能夠估算出車輛關(guān)鍵部件實(shí)際參數(shù)值，因此相對(duì)信號(hào)分析方法具備一定的優(yōu)勢(shì)。文獻(xiàn)[3，4]提出了一種基于卡爾曼濾波和粒子濾波的參數(shù)估計(jì)算法，利用該算法估計(jì)了軌道車輛正常狀態(tài)下二系橫向減震器阻尼、抗蛇行減震器阻尼及輪對(duì)等效錐度的參數(shù)值，獲得較為理想的效果，但是沒有深入研究車輛故障狀態(tài)下這些部件的參數(shù)估計(jì)算法。文獻(xiàn)[5]提出一種新型的軌道車輛懸掛系統(tǒng)狀態(tài)監(jiān)測(cè)方法，引入了再次均勻采樣策略，有效克服了傳統(tǒng)參數(shù)估計(jì)算法因粒子枯竭而無法對(duì)懸掛系統(tǒng)故障進(jìn)行監(jiān)測(cè)的缺陷，但是算法中所需的觀測(cè)值是通過車輛系統(tǒng)仿真模擬得到，沒有深入研究從軌道車輛SIMPACK整車模型或?qū)嵻囋囼?yàn)中獲得觀測(cè)信息。

本文主要研究RBPF參數(shù)估計(jì)方法，利用多體動(dòng)力學(xué)仿真和實(shí)車試驗(yàn)獲取振動(dòng)數(shù)據(jù)，結(jié)合RBPF算法估計(jì)軌道車輛在不同工況下的二系橫向減震器實(shí)際參數(shù)值，估計(jì)結(jié)果驗(yàn)證了參數(shù)估計(jì)方法的有效性和準(zhǔn)確性。

1 軌道車輛SIMPACK整車模型的建立

以國內(nèi)某典型軌道車輛為研究對(duì)象，車輛部分參數(shù)見表1，利用SIMPACK軟件建立軌道車輛整車模型。在多體系統(tǒng)中，拓?fù)浣Y(jié)構(gòu)建立是系統(tǒng)建模的一個(gè)十分關(guān)鍵的前期工作，多體系統(tǒng)拓?fù)鋱D建立的好壞直接和模型建立的正確與否密切相關(guān)。拓?fù)潢P(guān)系主要由車體、構(gòu)架和輪對(duì)三大部件組成[6]，通過約束形式、鉸接形式以及力元形式來表示不同剛體之間的連接關(guān)系，軌道車輛系統(tǒng)拓?fù)鋱D如圖1所示。

表1 車輛參數(shù)

利用SIMPACK軟件建立的軌道車輛整車模型如圖2所示，文中主要利用建立的整車模型獲得車體、構(gòu)架、輪對(duì)的振動(dòng)加速度、角加速度等后處理信息，作為RBPF算法中的實(shí)際觀測(cè)值，對(duì)軌道車輛懸掛系統(tǒng)關(guān)鍵部件參數(shù)進(jìn)行估計(jì)。

圖1 軌道車輛系統(tǒng)拓?fù)鋱D

圖2 軌道車輛整車模型

2 軌道車輛系統(tǒng)橫向動(dòng)力學(xué)半車模型和狀態(tài)空間模型建立

軌道車輛后半車模型平面圖如圖3所示，對(duì)應(yīng)于圖2中輪對(duì)1、輪對(duì)2、后轉(zhuǎn)向架構(gòu)架及半個(gè)車體，由于只建立半車模型，因此車體只考慮1個(gè)橫向運(yùn)動(dòng)的自由度，而輪對(duì)和構(gòu)架均有橫向和搖頭運(yùn)動(dòng)2個(gè)自由度。圖中：yw1和yw2分別為輪對(duì)1、輪對(duì)2的橫移量；ψw1和ψw2分別為輪對(duì)1、輪對(duì)2的搖頭角；yb和ybd分別為構(gòu)架和車體的橫向位移；ψb為構(gòu)架的搖頭角。同理，也可建立由輪對(duì)3、輪對(duì)4、前轉(zhuǎn)向架構(gòu)架及半個(gè)車體構(gòu)成的前半車模型。

圖3 軌道車輛后半車模型平面圖

軌道車輛系統(tǒng)橫向動(dòng)力學(xué)半車模型振動(dòng)方程為

(1)

(2)

(3)

(4)

(5)

(6)

(7)

(8)

(9)

(10)

式中：yt1和yt2分別為作用在2個(gè)輪對(duì)上的軌道橫向不平順；yd1和yd2分別為2個(gè)輪對(duì)在輪軌接觸點(diǎn)處的橫向位移；Iw1和Iw2分別為2個(gè)輪對(duì)的轉(zhuǎn)動(dòng)慣量；Ib為構(gòu)架的轉(zhuǎn)動(dòng)慣量；f11和f22分別為縱向和橫向蠕滑系數(shù)；λ為輪對(duì)等效錐度；v為車輛運(yùn)行速度。將系統(tǒng)橫向動(dòng)力學(xué)模型變換為狀態(tài)空間模型如下

(11)

y(t)=Hx(t)+v(t)

(12)

式(11)為時(shí)間t的連續(xù)函數(shù)，將其離散化[7]并結(jié)合式(12)可得到離散化后的系統(tǒng)橫向動(dòng)態(tài)空間模型為

xk+1=A(θ)xk+B(θ)wk

(13)

yk+1=H(θ)xk+1+vk+1

(14)

式中：w和v分別為相互獨(dú)立的狀態(tài)噪聲和觀測(cè)噪聲向量；A(θ)，B(θ)和H(θ)是由θ決定的矩陣。

3 RBPF算法在軌道車輛懸掛系統(tǒng)參數(shù)估計(jì)中的應(yīng)用

圖4 RBPF參數(shù)估計(jì)算法流程圖

此時(shí)，系統(tǒng)動(dòng)態(tài)空間模型可以表示為

xk=A(θ)xk-1+B(θ)wk-1

(15)

yk=H(θ)xk+vk

(16)

參數(shù)估計(jì)過程中狀態(tài)噪聲w為軌道橫向不平順，方差陣Q[9]為

式中：Ar為軌道粗糙度因素，取Ar=2×0.33×10-3；Ts為采樣時(shí)間間隔，取Ts=0.001 s。

觀測(cè)噪聲協(xié)方差矩陣設(shè)為均值為零的定長平穩(wěn)隨機(jī)白噪聲序列，方差陣R為[10]

R=0.001|max(yk(n))-min(yk(n))|

n=1，2，3，4，5，6，7

式中：max(yk(n))和min(yk(n))分別是第n個(gè)觀測(cè)量在整個(gè)時(shí)間歷程中的最大值和最小值。

假設(shè)初始狀態(tài)估計(jì)值為x0，初始的狀態(tài)估計(jì)方差P0是17階稀疏矩陣，一般方差P0設(shè)置較小，這樣可避免因?qū)V波值造成過大的擾動(dòng)而導(dǎo)致失真。未知參數(shù)向量θ的變化范圍為[θmin，θmax]，粒子數(shù)量為M。RBPF算法步驟如下：

步驟1初始化。在區(qū)間[θmin，θmax]內(nèi)均勻采樣，組成粒子集θ1|0(i)，初始狀態(tài)x1|0(i)=x0，初始狀態(tài)方差P1|0(i)=P0，i=1，2，3，…，M。

步驟2在每個(gè)時(shí)刻k(k=1，2，…)，重復(fù)以下步驟：

yk|k-1(i)=H(θk|k-1(i))xk|k-1(i)

(17)

Rk(i)=H(θk|k-1(i))Pk|k-1(i)HΤ(θk|k-1(i))+R

(18)

(19)

(20)

式中：各時(shí)刻觀測(cè)向量yk可從SIMPACK整車模型和實(shí)車試驗(yàn)中獲得，各個(gè)粒子的權(quán)值服從多維正態(tài)分布函數(shù)。

(21)

(3)重采樣。首先產(chǎn)生M個(gè)在[0，1]區(qū)間均勻分布的隨機(jī)數(shù)up(p=1，2，…，M)，通過搜索算法找到滿足以下條件的整數(shù)l。

(22)

l=1，2，3，…，M

記錄第l個(gè)粒子，即對(duì)該粒子進(jìn)行了M重復(fù)采樣，同時(shí)保留相應(yīng)的θk(i)，xk|k-1(l)，Pk|k-1(l)，得到新的參數(shù)、狀態(tài)和方差組合。

(4)卡爾曼濾波測(cè)量更新。獲得k時(shí)刻卡爾曼濾波增益Kk(i)、狀態(tài)估計(jì)值xk(i)估計(jì)均方誤差Pk(i)分別為

Rk(i)=H(θk(i))Pk|k-1(i)HT(θk(i))+R

(23)

(24)

xk(i)=xk|k-1(i)+Kk(i)(yk-H(θk(i))xk|k-1(i))

(25)

Pk(i)=Pk|k-1(i)-Kk(i)H(θk(i))Pk|k-1(i)

(26)

(5)粒子濾波時(shí)間預(yù)測(cè)更新。

(27)

(6)卡爾曼濾波測(cè)量更新，得到k+1時(shí)刻估計(jì)的狀態(tài)xk+1|k(i)和估計(jì)方差Pk+1|k(i)一步預(yù)測(cè)更新。

xk+1|k(i)=A(θk+1|k(i))xk(i)

(28)

Pk+1|k(i)=A(θk+1|k(i))Pk(i)AT(θk+1|k(i))+

B(θk+1|k(i))QBT(θk+1|k(i))

(29)

4 正常情況下算法驗(yàn)證

4.1 試驗(yàn)?zāi)康募皸l件

選取既定的軌道車輛進(jìn)行試驗(yàn)，目前該車轉(zhuǎn)向架二系橫向減震器處于工作正常狀態(tài)，利用真實(shí)的振動(dòng)測(cè)試數(shù)據(jù)作為RBPF參數(shù)估計(jì)算法中觀測(cè)值，估計(jì)前、后轉(zhuǎn)向架二系橫向減震器實(shí)際阻尼值，從而驗(yàn)證車輛關(guān)鍵部件正常情況下參數(shù)估計(jì)算法的正確性。

4.1.1 試驗(yàn)平臺(tái)和系統(tǒng)組成

整個(gè)車載系統(tǒng)主要由數(shù)據(jù)采集器、工業(yè)交換機(jī)、車載服務(wù)器和GPS(測(cè)速模塊)構(gòu)成，車載服務(wù)器通過車內(nèi)局域網(wǎng)對(duì)整個(gè)采集系統(tǒng)進(jìn)行控制。車載服務(wù)器系統(tǒng)采集各部件的振動(dòng)狀態(tài)信息。速度信息由GPS提供，并獲取當(dāng)前車輛運(yùn)行工況。在整個(gè)試驗(yàn)過程中，車輛檢測(cè)系統(tǒng)通過存儲(chǔ)卡將原始信號(hào)數(shù)據(jù)全部本地存儲(chǔ)，待停止試驗(yàn)后取出振動(dòng)數(shù)據(jù)。

4.1.2 車輛參數(shù)和傳感器布置

圖5 車體加速度計(jì)安裝實(shí)物圖

圖6 構(gòu)架加速度計(jì)安裝實(shí)物圖

圖7 軸箱加速度計(jì)安裝實(shí)物圖

圖8 構(gòu)架陀螺儀安裝實(shí)物圖

4.2 試驗(yàn)結(jié)果分析

文中主要以后半車模型懸掛系統(tǒng)為例，估計(jì)二系橫向減震器阻尼值。整個(gè)計(jì)算過程中，選取一段時(shí)間的數(shù)據(jù)，其速度-時(shí)間關(guān)系曲線如圖9所示。

圖9 速度-時(shí)間關(guān)系曲線

對(duì)現(xiàn)場(chǎng)數(shù)據(jù)采集器采集的實(shí)際振動(dòng)測(cè)試數(shù)據(jù)，采用巴特沃斯低通濾波器進(jìn)行預(yù)處理，利用MATLAB軟件編寫程序來驗(yàn)證正常情況下RBPF算法的正確性，選取粒子數(shù)量M=200，圖10分別給出了車速為v=5 km/h、v=10 km/h、v=15 km/h、v=20 km/h四種不同工況下的參數(shù)估計(jì)結(jié)果。

(a)v=5 km/h

(b)v=10 km/h

(c)v=15 km/h

(d)v=20 km/h圖10 不同速度下cyb估計(jì)結(jié)果

圖10中，二系橫向減震器cyb估計(jì)值已進(jìn)行歸一化處理，即圖中縱坐標(biāo)數(shù)值為“1”代表二系橫向減震器部件正常(cyb=29 500 N·s·m-1)；縱坐標(biāo)數(shù)值為“0.5”代表二系橫向減震器阻尼參數(shù)值縮小1倍(cyb=14 750 N·s·m-1)；縱坐標(biāo)數(shù)值為“0”代表二系橫向減震器已完全失效，以此類推。由圖10可知：實(shí)車試驗(yàn)獲得的觀測(cè)值雖然少了輪對(duì)搖頭角加速度，估計(jì)值仍能夠很快的收斂到正常值(這里指cyb無故障時(shí)實(shí)際參數(shù)值)附近，并且車速快慢對(duì)參數(shù)估計(jì)結(jié)果影響不大。雖然在整個(gè)試驗(yàn)過程中可能會(huì)因?yàn)閭鞲衅鞯臏y(cè)量誤差及安裝位置會(huì)造成估計(jì)值與正常值有一定的偏差，但參數(shù)估計(jì)值整體趨勢(shì)基本與正常值相吻合，參數(shù)估計(jì)結(jié)果與實(shí)際情況相符，驗(yàn)證了正常情況下參數(shù)估計(jì)算法的正確性，用同樣的方法估計(jì)得到前轉(zhuǎn)向架二系橫向減震器目前同樣處于正常狀態(tài)，估計(jì)結(jié)果與實(shí)際情況也相符。

5 模擬故障情況下算法驗(yàn)證

(a)前轉(zhuǎn)向架cyb縮小至0.5

(b)后轉(zhuǎn)向架cyb正常

(c)前轉(zhuǎn)向架cyb正常

(d)后轉(zhuǎn)向架cyb縮小至0.5

(e)前轉(zhuǎn)向架cyb縮小至0.7

(f)后轉(zhuǎn)向架cyb縮小至0.3

(h)前轉(zhuǎn)向架cyb正常

(g)后轉(zhuǎn)向架cyb正常圖11 前后轉(zhuǎn)向架參數(shù)估計(jì)結(jié)果

圖11是軌道車輛前、后轉(zhuǎn)向架二系橫向減震器不同故障下進(jìn)行的參數(shù)估計(jì)，分為前轉(zhuǎn)向架cyb故障，后轉(zhuǎn)向架cyb正常；前轉(zhuǎn)向架cyb正常，后轉(zhuǎn)向架cyb故障；前、后轉(zhuǎn)向架cyb均故障；前、后轉(zhuǎn)向架cyb均正常四種工況。以前、后轉(zhuǎn)向架cyb均故障為例，圖11(e)表示前轉(zhuǎn)向架cyb阻尼值縮小至正常值0.7倍左右，而圖11(f)表示后轉(zhuǎn)向架cyb阻尼值縮小至正常值0.3倍左右，相比得到后轉(zhuǎn)向架二系橫向減震器故障程度比前轉(zhuǎn)向架更嚴(yán)重。從圖11可以看出，利用SIMPACK整車模型輸出的觀測(cè)向量yk，結(jié)合RBPF參數(shù)估計(jì)算法都能夠有效對(duì)二系橫向減震器cyb故障情況下實(shí)際參數(shù)值進(jìn)行估計(jì)。計(jì)算步數(shù)達(dá)到約500時(shí)估計(jì)值就能收斂到真實(shí)值附近，并且估計(jì)值非常接近參數(shù)的真實(shí)值，故障情況下算法的正確性得到驗(yàn)證，監(jiān)測(cè)實(shí)時(shí)性也得到保障。

6 結(jié)束語

本文建立了國內(nèi)某型軌道車輛橫向動(dòng)力學(xué)半車模型和狀態(tài)空間模型，引入RBPF算法，結(jié)合多體動(dòng)力學(xué)仿真和實(shí)車試驗(yàn)的振動(dòng)數(shù)據(jù)，對(duì)軌道車輛懸掛系統(tǒng)關(guān)鍵部件之一的二系橫向減震器參數(shù)進(jìn)行估計(jì)，結(jié)果表明：該方法能準(zhǔn)確有效地估計(jì)得到懸掛系統(tǒng)關(guān)鍵部件正常或故障工況下實(shí)際參數(shù)值，為車輛關(guān)鍵部件在線狀態(tài)監(jiān)測(cè)和安全預(yù)警提供了重要的理論依據(jù)。

參考文獻(xiàn)：

[1]陳龍.基于粒子濾波的城軌車輛狀態(tài)監(jiān)測(cè)方法研究[D].上海:上海工程技術(shù)大學(xué),2013:56-58.

[2]漆瑾.高速機(jī)車走行部實(shí)時(shí)在線監(jiān)測(cè)與診斷方法的研究[D].成都：西南交通大學(xué),2006:4-6.

[3] LI P,GOODALL R, WESTON P.Estimation of Railway Vehicle Suspension Parameters for Condition Monitoring[J].Control Engineering Practice,2007,15(1):43-55.

[4]LI P,GOODALL R, KADIRKAMANATHAN V. Estimation of Parameters in a Linear State Space Model Using a Rao Blackwellised Particle Filter[J].Control Theory Application,2004,151(6):727-738.

[5]方宇，陳龍，李曉斌，等.基于再次均勻采樣策略改進(jìn)軌道車輛二系懸掛系統(tǒng)參數(shù)估計(jì)方法[J].中國鐵道科學(xué),2013,34(3)：72-78.

FANG Yu,CHEN Long,LI Xiaobin,et al.Estimation Method for the Secondary Suspension System Parameters of Railway Vehicle Improved by Repeated Uniform Sampling Strategy[J].China Railway Science,2013,34(3):72-78.

[6]繆炳榮，羅仁，王哲，等.SIMPACK動(dòng)力學(xué)分析高級(jí)教程[M].成都：西南交通大學(xué)出版社，2010.

[7]方宇，陳龍，鄭樹彬，等.基于參數(shù)估計(jì)的軌道車輛懸掛系統(tǒng)狀態(tài)監(jiān)測(cè)方法[J].鐵道學(xué)報(bào)，2013，35(5)：15-20.

FANG Yu,CHEN Long,ZHENG Shubin,et al.Condition Monitoring of Rail Vehicle Suspension System Based on Parameter Estimation[J].Journal of the China Railway Society,2013,35(5):15-20.

[8]DOUCET A,GODSILL S,ANDRIEU C.On Sequential Monte Carlo Sampling Methods for Bayesian Filtering[J].Statistics and Computing,2000,10(3):197-208.

[9]LI H.Measuring Systems for Active Steering of Railway Vehicles[D]. Loughborough: Loughborough University, 2001.

[10]林棻，趙又群，徐朔南.基于粒子濾波算法的汽車狀態(tài)估計(jì)技術(shù)[J].農(nóng)業(yè)機(jī)械學(xué)報(bào),2011,42(2)：23-27.

LIN Fen,ZHAO Youqun,XU Shuonan.Vehicle States Estimation Technology Based on Particle Filter Algorithm[J].Transactions of the Chinese Society for Agricultural Machinery,2011,42(2)：23-27.

[11]丁建明.車輛動(dòng)力學(xué)性能參數(shù)估計(jì)方法研究[D].成都:西南交通大學(xué)，2011:38-39.

[12]劉春生，胡壽松.一類基于狀態(tài)估計(jì)的非線性系統(tǒng)的智能故障診斷[J].控制與決策,2005,20(5):557-561.

LIU Chunsheng,HU Shousong.Intelligent Nonlinear Fault Diagnosis Based on State Estimator[J].Control and Decision,2005,20(5):557-561.