楊 濤,張鈺奇,付春健,趙華東

(1.鄭州大學(xué) 機(jī)械與動(dòng)力工程學(xué)院,鄭州 450001; 2.河南省智能制造研究院,鄭州 450001)

0 引 言

水工閘門由于工作環(huán)境惡劣,在運(yùn)行過(guò)程中會(huì)不可避免地出現(xiàn)性能退化,最終導(dǎo)致破壞。一旦閘門發(fā)生破壞,會(huì)造成巨大的經(jīng)濟(jì)損失。因此準(zhǔn)確預(yù)測(cè)閘門剩余壽命(Remaining Useful Life,RUL)對(duì)于實(shí)時(shí)評(píng)估閘門健康狀態(tài),確定合適的維修策略以保證閘門運(yùn)行安全具有重要意義。

目前研究主要認(rèn)為銹蝕是影響閘門壽命的最大因素[1]。Kuz’Mitskii等[2]基于銹蝕分階段模型預(yù)測(cè)了閘門剩余壽命;李偉康等[3]考慮到銹蝕會(huì)影響材料性能,基于D-H曲線對(duì)閘門進(jìn)行RUL預(yù)測(cè);然而,目前的研究存在著預(yù)測(cè)指標(biāo)單一、精度不高,難以給出剩余壽命的概率分布等問(wèn)題。

目前,閘門上的銹蝕模型一般為確定性函數(shù),難以體現(xiàn)銹蝕過(guò)程的不確定性。而Gamma過(guò)程可以很好地模擬銹蝕演化過(guò)程[4],體現(xiàn)閘門銹蝕過(guò)程的不確定性。因此,考慮采用Gamma過(guò)程進(jìn)行銹蝕模擬,通過(guò)數(shù)值模擬方法得到閘門退化數(shù)據(jù)。在獲取閘門的退化數(shù)據(jù)后,同時(shí)考慮到工程實(shí)際中設(shè)備運(yùn)行工況復(fù)雜,反映設(shè)備退化的指標(biāo)往往并不唯一。因此,如何將多種信息進(jìn)行融合來(lái)共同反映設(shè)備退化狀態(tài),已成為剩余壽命預(yù)測(cè)的關(guān)鍵。目前,許多研究通過(guò)將多種特征融合成健康因子(Health Index,HI)來(lái)解決這一問(wèn)題。任子強(qiáng)等[5]將多維數(shù)據(jù)融合成一個(gè)健康指標(biāo)來(lái)表征發(fā)動(dòng)機(jī)的退化過(guò)程,提高了RUL預(yù)測(cè)的精度;趙廣社等[6]基于歐式距離構(gòu)建了健康因子,并對(duì)發(fā)動(dòng)機(jī)進(jìn)行了RUL預(yù)測(cè)。

RUL預(yù)測(cè)的方法可以大致分為基于數(shù)據(jù)驅(qū)動(dòng)和退化模型兩種[7]。其中,基于退化模型的方法不依賴于大量同類型設(shè)備的全壽命周期數(shù)據(jù),而是根據(jù)歷史數(shù)據(jù)建立退化模型來(lái)描述設(shè)備的退化趨勢(shì),因此得到了廣泛的運(yùn)用[8]。LI等[9]利用維納過(guò)程刻畫了渦輪發(fā)動(dòng)機(jī)的退化過(guò)程;張英波等[10]基于Gamma過(guò)程建立了行星架RUL預(yù)測(cè)模型。

粒子濾波是基于退化模型實(shí)現(xiàn)RUL預(yù)測(cè)的有效手段,在非線性、非高斯系統(tǒng)中具有明顯優(yōu)勢(shì)?；诹Ｗ訛V波的RUL預(yù)測(cè)方法可以給出剩余壽命的概率分布,改善長(zhǎng)期預(yù)測(cè)的不確定性問(wèn)題,相較于單一的預(yù)測(cè)值,更具有參考價(jià)值。李玥鋅等[11]基于維納過(guò)程對(duì)電池退化過(guò)程進(jìn)行建模,并基于粒子濾波算法得到了其剩余壽命概率分布;王璽等[12]基于粒子濾波算法對(duì)陀螺儀進(jìn)行了RUL預(yù)測(cè),取得較高的精度。

基于上述分析,本文提出了一種多特征信息融合預(yù)測(cè)閘門剩余壽命的方法,利用非線性維納過(guò)程對(duì)閘門進(jìn)行退化建模,并基于粒子濾波算法實(shí)時(shí)跟蹤閘門退化曲線并進(jìn)行RUL預(yù)測(cè),同時(shí)給出了閘門剩余壽命的概率分布,可以更好地指導(dǎo)工程實(shí)際。

1 閘門銹蝕過(guò)程模擬

由于閘門運(yùn)行工況復(fù)雜,影響其銹蝕的因素(如水流、干濕循環(huán)等)具有不確定性,因此閘門的銹蝕過(guò)程也伴隨著隨機(jī)性。此外,銹蝕過(guò)程還具有累積性。為了表征銹蝕隨時(shí)間發(fā)展的不確定性和累積性,考慮采用Gamma隨機(jī)過(guò)程模擬銹蝕深度的發(fā)展過(guò)程。

Gamma過(guò)程是一種具有非負(fù)增量的隨機(jī)過(guò)程,特別適合模擬具有微小增量的累計(jì)漸變過(guò)程,如銹蝕、疲勞、蠕變等[13]。假設(shè)t時(shí)刻的銹蝕量d(t)滿足Gamma分布Ga(α(t),β),其概率密度函數(shù)為

(1)

式中:d為銹蝕量;α(t)和β分別為形狀參數(shù)和尺度參數(shù);Γ(α(t))為Gamma函數(shù)。

根據(jù)Gamma過(guò)程的性質(zhì),t時(shí)刻銹蝕量的均值可以表示為

E[d(t)] =α(t)β。

(2)

當(dāng)采用Gamma過(guò)程描述退化過(guò)程時(shí),一般認(rèn)為尺度參數(shù)β為常數(shù),形狀參數(shù)α(t)與時(shí)間的冪成正比[14],即

α(t)=ctb;c>0,b>0 。

(3)

式中c和b均為常數(shù)。對(duì)于特定退化過(guò)程b一般為定值,文獻(xiàn)[15] 給出了不同退化類型下b的取值,閘門銹蝕與鋼筋銹蝕屬于同類,取b=1。

當(dāng)確定參數(shù)b之后,模型中的參數(shù)c和β可以結(jié)合既有結(jié)構(gòu)的銹蝕觀測(cè)數(shù)據(jù)采用極大似然估計(jì)法或矩估計(jì)法來(lái)進(jìn)行估計(jì)。然而,由于缺乏閘門長(zhǎng)期銹蝕監(jiān)測(cè)數(shù)據(jù),同時(shí)考慮到t時(shí)刻的銹蝕深度與平均銹蝕速率有關(guān),因此利用已有銹蝕率的統(tǒng)計(jì)數(shù)據(jù)對(duì)參數(shù)c和β進(jìn)行反向標(biāo)定[16],即首先假定參數(shù)β為某個(gè)固定值,根據(jù)式(2)并結(jié)合銹蝕速率統(tǒng)計(jì)數(shù)據(jù)即可確定參數(shù)c。

當(dāng)確定參數(shù)之后,閘門銹蝕深度發(fā)展過(guò)程便可以采用Gamma順序采樣法進(jìn)行模擬:

(1)將總時(shí)間T等分為n個(gè)子區(qū)間。

(2)當(dāng)t=ti時(shí),生成服從Ga(α(ti)-α(ti-1),β)的隨機(jī)數(shù)Δdi,Δdi便為時(shí)段ti-ti-1的銹蝕增量。則ti時(shí)刻的銹蝕累積深度為di=di-1+Δdi。

(3)重復(fù)以上過(guò)程,直至達(dá)到總時(shí)間T。

2 健康因子的構(gòu)建

閘門的退化過(guò)程可由多種特征參數(shù)變化來(lái)表示,不同的特征對(duì)于閘門退化過(guò)程的敏感性不同,并且各個(gè)特征之間存在一定的相關(guān)性。若將多個(gè)特征指標(biāo)進(jìn)行單獨(dú)建模,不僅不符合多變量相互耦合實(shí)際情況,也會(huì)增加建模的復(fù)雜性。因此,對(duì)于多特征的情況,常常采用信息融合的方法,利用降維等手段將多維數(shù)據(jù)轉(zhuǎn)化為一維復(fù)合健康指標(biāo),這種方法不僅充分考慮了多變量相互關(guān)聯(lián)的情況,也更有利于對(duì)閘門退化過(guò)程進(jìn)行建模。因此,將多種特征指標(biāo)進(jìn)行融合得到能夠表征閘門退化過(guò)程的健康因子是進(jìn)行壽命預(yù)測(cè)的關(guān)鍵。

2.1 數(shù)據(jù)獲取

閘門銹蝕會(huì)使構(gòu)件截面面積變小,空間結(jié)構(gòu)整體效應(yīng)發(fā)生改變,從而影響閘門的強(qiáng)度和剛度。同時(shí),閘門銹蝕會(huì)引起自振頻率和振型等模態(tài)參數(shù)的改變。為此,本文綜合考慮了閘門靜力學(xué)和動(dòng)力學(xué)特性,選取了應(yīng)力、變形、主梁撓度,“干模態(tài)”下前10階自振頻率、位移模態(tài)、應(yīng)變模態(tài),以及考慮流固耦合作用下的“濕模態(tài)”前10階自振頻率、位移模態(tài)、應(yīng)變模態(tài),共63種特征指標(biāo)來(lái)共同反映閘門退化過(guò)程。

2.2 特征指標(biāo)篩選

特征篩選是指從特征集中選擇出最優(yōu)的子集來(lái)反映閘門退化過(guò)程?？紤]到單調(diào)性強(qiáng)的特征能夠很好地表征閘門的退化趨勢(shì),因此首先對(duì)特征進(jìn)行單調(diào)性篩選。Spearman相關(guān)系數(shù)可用于評(píng)估2個(gè)變量之間相關(guān)關(guān)系的程度,考慮到時(shí)間序列是嚴(yán)格單調(diào)的,因此計(jì)算每一種特征時(shí)序與對(duì)應(yīng)時(shí)間序列的Spearman作為該種特征的單調(diào)性得分。同時(shí),考慮到變化程度大的特征能更好地反映退化過(guò)程。因此,綜合考慮單調(diào)性和離散性對(duì)特征指標(biāo)進(jìn)行篩選。

對(duì)于某一特征參數(shù)序列X=(x1,x2,…,xn),相應(yīng)的時(shí)間序列T1=(t1,t2,…,tn),則:

(1)單調(diào)性指標(biāo)為

(4)

式中:n為數(shù)據(jù)的個(gè)數(shù);分別將X,T1中的元素按照從小到大的順序編秩;Ri表示xi在X中的秩次;Qi表示ti在T1中的秩次。

(2)離散性指標(biāo)為

為發(fā)現(xiàn)最優(yōu)特征指標(biāo),首先選擇出單調(diào)性較強(qiáng)的特征指標(biāo),然后再對(duì)離散性高的特征指標(biāo)進(jìn)行選擇。

2.3 特征融合及健康因子的構(gòu)建

多特征參數(shù)可以很好地表征設(shè)備退化狀態(tài),也會(huì)使設(shè)備建模和數(shù)據(jù)處理復(fù)雜化。因此,要將多維數(shù)據(jù)融合成一維數(shù)據(jù),即通過(guò)采用健康因子來(lái)表征閘門的退化過(guò)程。

主成分分析(Principal Component Analysis,PCA)是一種數(shù)據(jù)融合和特征提取的重要方法,主要思想是通過(guò)尋找一組正交基,將高維特征映射到低維空間,從而達(dá)到降維的目的。其主要步驟為:

(1)將性能退化指標(biāo)構(gòu)成特征數(shù)據(jù)矩陣Y為

(7)

(3)計(jì)算P的特征值m1,m2,…,mp并從大到小進(jìn)行排序,計(jì)算對(duì)應(yīng)的特征向量e1,e2,…,ep。

(4)根據(jù)特征值mq,計(jì)算第q個(gè)主成分對(duì)應(yīng)的貢獻(xiàn)率

(8)

一般取貢獻(xiàn)率ξ>80%的主成分可以較全面地反映出特征的信息[17],將其當(dāng)做特征融合后的健康因子來(lái)反映設(shè)備退化過(guò)程。

3 閘門壽命預(yù)測(cè)

3.1 粒子濾波原理

粒子濾波是一種基于蒙特卡洛(Monte Carlo,MC)理論的貝葉斯估計(jì)算法,可用于處理非線性、非高斯系統(tǒng)的狀態(tài)估計(jì)問(wèn)題。貝葉斯算法是粒子濾波的基礎(chǔ),通常采用狀態(tài)轉(zhuǎn)移方程和觀測(cè)方程來(lái)描述系統(tǒng)的狀態(tài)空間

(9)

式中:xk-1、xk分別表示k-1和k時(shí)刻系統(tǒng)的狀態(tài)值;fk表示系統(tǒng)的狀態(tài)傳遞函數(shù);εk-1表示系統(tǒng)的過(guò)程噪聲;yk表示k時(shí)刻系統(tǒng)的觀測(cè)值;hk表示觀測(cè)函數(shù);vk表示系統(tǒng)的測(cè)量噪聲。

貝葉斯濾波是對(duì)系統(tǒng)的后驗(yàn)分布進(jìn)行求解,首先根據(jù)狀態(tài)轉(zhuǎn)移方程,根據(jù)k-1時(shí)刻的觀測(cè)值,可預(yù)測(cè)得到k時(shí)刻系統(tǒng)狀態(tài)的先驗(yàn)分布

p(xk|y1:k-1)=

(10)

式中:p(xk|xk-1)由狀態(tài)轉(zhuǎn)移方程定義;y1:k-1為1:(k-1)時(shí)刻的系統(tǒng)觀測(cè)值。

在此基礎(chǔ)上,根據(jù)k時(shí)刻觀測(cè)值,利用貝葉斯理論得到k時(shí)刻的后驗(yàn)概率密度分布

式中p(yk|xk)由狀態(tài)測(cè)量方程定義。

由式(10)、式(11)可知,上述過(guò)程存在著復(fù)雜的積分運(yùn)算,往往難以得到精確解。為了解決上述問(wèn)題,粒子濾波算法借鑒蒙特卡羅方法將k時(shí)刻的后驗(yàn)分布離散加權(quán)為

(12)

由式(12)可知,得到粒子對(duì)應(yīng)的權(quán)值分布,就可以實(shí)現(xiàn)對(duì)狀態(tài)概率分布的近似估計(jì),假設(shè)q(x0:k|y1:k)表示過(guò)去時(shí)刻狀態(tài)的后驗(yàn)估計(jì),同時(shí)狀態(tài)和觀測(cè)信息遵循一階馬爾科夫過(guò)程,則序貫重要性采樣粒子權(quán)重可轉(zhuǎn)換為

(13)

然而,隨著迭代次數(shù)的增加,權(quán)重有可能只存在少數(shù)粒子中,大部分粒子的權(quán)重可能小到忽略不計(jì),造成估計(jì)性能的下降。因此,常常通過(guò)采用重采樣方法降低粒子退化問(wèn)題。

3.2 退化模型的建立

閘門在運(yùn)行過(guò)程中由于工作環(huán)境惡劣,其退化過(guò)程呈現(xiàn)出明顯的不確定性。Wiener過(guò)程可以很好地描述復(fù)雜設(shè)備退化過(guò)程中的不確定性,在RUL預(yù)測(cè)中得到了廣泛運(yùn)用?？紤]到以上因素,利用非線性維納過(guò)程來(lái)刻畫閘門隨機(jī)退化過(guò)程,即

X(t)=X(0)+λeθt+σBB(t) 。

(14)

式中:X(t)表示t時(shí)刻的系統(tǒng)狀態(tài);X(0)為系統(tǒng)初始狀態(tài);λ為漂移系數(shù);eθt為時(shí)間t和參數(shù)θ的非線性函數(shù),用于表示退化過(guò)程中的非線性;σB為擴(kuò)散系數(shù);B(t)為標(biāo)準(zhǔn)維納過(guò)程。

因此,閘門的退化過(guò)程可以描述為

xk=xk-1+λk-1(eθk-eθ(k-1))+ηk。

(15)

式中ηk=σB(B(tk)-B(tk-1))服從于正態(tài)分布N(0,σ2Δt),Δt=tk-tk-1。

事實(shí)上,由于測(cè)量結(jié)果與真實(shí)狀態(tài)之間不可避免地存在誤差,因此,測(cè)量狀態(tài)和真實(shí)狀態(tài)之間的關(guān)系可以表示為

yk=xk+v。

(16)

式中v為測(cè)量噪聲,服從于正態(tài)分布N(0,γ2)。

同時(shí),由于過(guò)程噪聲可以由其他參數(shù)的不確定性實(shí)現(xiàn)[18]。因此,根據(jù)退化模型可以得出其狀態(tài)空間模型,即

由式(17)可得,利用Yk=[xk,λk,θk,σk,γk] 可表示出系統(tǒng)的狀態(tài),并結(jié)合第3.1節(jié)所敘述的粒子濾波原理,利用在線數(shù)據(jù)實(shí)時(shí)更新未知參數(shù),便可以實(shí)現(xiàn)RUL預(yù)測(cè)。

3.3 剩余壽命預(yù)測(cè)

在得到狀態(tài)空間方程之后,根據(jù)第3.2節(jié)所述的粒子濾波算法對(duì)模型參數(shù)進(jìn)行更新,可以得到t時(shí)刻的系統(tǒng)狀態(tài),根據(jù)退化模型進(jìn)行遞推,可以對(duì)t+φ時(shí)刻的系統(tǒng)狀態(tài)進(jìn)行預(yù)測(cè)[19]。退化模型遞推公式為

1≤ζ≤φ。

(18)

圖1 RUL預(yù)測(cè)示意圖Fig.1 Schematic diagram of RUL prediction

4 案例驗(yàn)證

4.1 有限元模型建立

某溢洪道露頂式弧形鋼閘門尺寸為12 m×10.5 m(寬×高),主要材料為Q345B,彈性模量取2.06×105N/mm2,密度取7 850 kg/m3,泊松比取0.3,作用水頭取9.4 m,主梁計(jì)算跨度為7.6 m,現(xiàn)場(chǎng)原型如圖2所示。閘門為空間薄壁結(jié)構(gòu),采用Shell181殼單元模擬。閘門濕模態(tài)仿真時(shí),上游水體長(zhǎng)度取超過(guò)閘門高度10倍的長(zhǎng)度,取150 m,采用Fliud30單元模擬,水流和閘門的相互作用表現(xiàn)在閘門迎水面上,將其設(shè)為流固耦合面。

圖2 閘門現(xiàn)場(chǎng)Fig.2 Site photo of the gate

4.2 閘門銹蝕深度演化過(guò)程模擬

根據(jù)文獻(xiàn)[20] 提出的銹蝕線性模型,閘門上的銹蝕深度變化可根據(jù)式(19)計(jì)算,即

(19)

式中:dL(t)為閘門t時(shí)刻的銹蝕量;t為閘門運(yùn)行年限;v0為銹蝕速率;tI為閘門開始銹蝕的年限,一般取3～7 a,本文取5 a。

閘門上的銹蝕分為全面銹蝕和局部銹蝕,根據(jù)工程實(shí)際,全面銹蝕常常發(fā)生在面板、支臂、縱梁、肋板、次梁等結(jié)構(gòu)上,局部銹蝕常常發(fā)生在下主梁積水處、面板邊側(cè)和邊梁漏水處。根據(jù)統(tǒng)計(jì)我國(guó)閘門銹蝕速率的均值在(0.015～0.070)mm/a范圍內(nèi),考慮面板和局部銹蝕區(qū)域與水長(zhǎng)期接觸,因此對(duì)面板和局部銹蝕區(qū)域銹蝕速率取v1=0.07 mm/a,對(duì)于支臂、縱梁等其他區(qū)域取v2=0.06 mm/a。

假定尺度參數(shù)β=0.1,結(jié)合式(2)、式(19)便可以得到形狀參數(shù)α(t)為

α(t)=c(t-5),t>5 。

(20)

式中:v1=0.06 mm/a時(shí),c=0.6;v2=0.07 mm/a,c=0.7。

確定形狀參數(shù)α(t)和尺度參數(shù)β后,便可以利用Gamma順序采樣法模擬閘門銹蝕發(fā)展過(guò)程。值得說(shuō)明的是由于銹蝕過(guò)程存在著不確定性,利用Gamma順序采樣法會(huì)產(chǎn)生多條銹蝕曲線,通過(guò)實(shí)際檢測(cè)的銹蝕數(shù)據(jù)的增加,利用參數(shù)估計(jì)方法可減少這種不確定性。為了數(shù)據(jù)的可重復(fù)性和合理性,本文采用了設(shè)置隨機(jī)種子的方法,找到與線性模型相近的曲線,既能保證其銹蝕過(guò)程的不確定性,又保證其合理性,得到的銹蝕發(fā)展過(guò)程如圖3所示。圖3中線性增長(zhǎng)模型為文獻(xiàn)[20] 提出的銹蝕模型,其無(wú)法體現(xiàn)銹蝕過(guò)程中的不確定性;考慮銹蝕深度發(fā)展過(guò)程的不確定性,利用Gamma過(guò)程模擬銹蝕深度發(fā)展過(guò)程時(shí),其是由一系列不同長(zhǎng)度的微小增長(zhǎng)段構(gòu)成。

圖3 銹蝕深度演化過(guò)程Fig.3 Evolution process of corrosion depth

根據(jù)文獻(xiàn)[21] ,當(dāng)最大應(yīng)力發(fā)生處厚度≤16 mm時(shí),取[σ] =225 MPa。通過(guò)對(duì)不同運(yùn)行年限下的閘門進(jìn)行仿真,結(jié)果表明在29 a時(shí)閘門上最大應(yīng)力σmax>0.95[σ] =213.75 MPa達(dá)到失效的條件,圖4為閘門失效時(shí)刻應(yīng)力云圖。

圖4 閘門應(yīng)力云圖Fig.4 Stress contours of gate

4.3 閘門特征指標(biāo)篩選

利用第2.2節(jié)所敘述的單調(diào)性和離散性指標(biāo)對(duì)特征進(jìn)行篩選,特征指標(biāo)得分如圖5所示。圖5中橫坐標(biāo)1～3分別為應(yīng)力、變形、主梁撓度,4～33分別為干模態(tài)下前10階自振頻率、位移模態(tài)、應(yīng)變模態(tài),34～63分別為濕模態(tài)下前10階自振頻率、位移模態(tài)、應(yīng)變模態(tài)。

圖5 特征指標(biāo)得分Fig.5 Scores of characteristic indices

綜合考慮,先取單調(diào)性得分>0.9的特征指標(biāo),然后再?gòu)膯握{(diào)性篩選之后的指標(biāo)中進(jìn)行離散性篩選?？紤]到10階之后的離散性得分趨于平穩(wěn)且數(shù)值較小,對(duì)健康因子的構(gòu)建影響不大,因此取離散性得分前10特征用于融合健康因子,共篩選出應(yīng)力、變形干模態(tài)下3、4、5、8階位移模態(tài),濕模態(tài)下4、5、8階位移模態(tài),5階應(yīng)變模態(tài)在內(nèi)的10種特征,10種特征的退化曲線如圖6所示。閘門所篩選出的變形、位移模態(tài)、應(yīng)變模態(tài)等特征的初始狀態(tài)云圖如圖7所示。(由于應(yīng)力圖4已給出,本處不再列出)。

圖6 特征指標(biāo)退化曲線Fig.6 Degradation curves of feature indices

圖7 特征指標(biāo)云圖Fig.7 Contours of feature indices

4.4 閘門健康因子的構(gòu)建

根據(jù)2.3節(jié)所論述的PCA方法對(duì)篩選后的特征進(jìn)行融合以構(gòu)建健康因子,利用PCA方法得到的各主成分分量得分如圖8所示。

圖8 各主成分貢獻(xiàn)率Fig.8 Contribution rate of each principal component

其中,第一主成分貢獻(xiàn)率ξ=90.960 1%,滿足ξ>80%的要求,因此,取第一主成分構(gòu)建健康因子,健康因子曲線,如圖9所示。

圖9 健康因子曲線Fig.9 Health index curve

4.5 閘門剩余壽命預(yù)測(cè)

在得出健康因子退化曲線之后,便可以利用粒子濾波方法實(shí)時(shí)跟蹤退化曲線,并進(jìn)行RUL預(yù)測(cè)。在利用粒子濾波方法進(jìn)行壽命預(yù)測(cè)前,需先確定粒子的初始分布。然而事實(shí)上,對(duì)于一個(gè)特定的系統(tǒng)很難進(jìn)行足夠多的壽命測(cè)試以獲取可靠的參數(shù)先驗(yàn)分布[22]。因此,參數(shù)的先驗(yàn)分布一般是根據(jù)經(jīng)驗(yàn)得出,通過(guò)對(duì)閘門腐蝕退化過(guò)程進(jìn)行反復(fù)試驗(yàn)并結(jié)合先驗(yàn)知識(shí),發(fā)現(xiàn)初始分布為均勻分布更加符合閘門退化過(guò)程,具體參數(shù)如表1所示。

表1 初始參數(shù)分布Table 1 Initial parameter distribution

利用第3.1節(jié)所敘述的粒子濾波原理對(duì)閘門進(jìn)行RUL預(yù)測(cè),由于粒子濾波中的粒子是從初始分布中隨機(jī)產(chǎn)生的,導(dǎo)致每次循環(huán)預(yù)測(cè)結(jié)果在微小范圍內(nèi)波動(dòng)。為了更好地說(shuō)明粒子濾波進(jìn)行RUL預(yù)測(cè)的原理和步驟,圖10—圖12是粒子濾波單次循環(huán)預(yù)測(cè)的結(jié)果。粒子濾波實(shí)時(shí)跟蹤退化曲線如圖10所示,由圖10可以看出所提方法能很好地跟蹤退化曲線,并能有效地預(yù)測(cè)閘門RUL。

圖10 實(shí)時(shí)跟蹤退化曲線Fig.10 Tracking of degradation process

狀態(tài)跟蹤集(即不同的預(yù)測(cè)時(shí)刻)和粒子數(shù)是衡量粒子濾波算法有效性的關(guān)鍵參數(shù)。表2為算法在不同的組合方案下的仿真,仿真次數(shù)為100次,從表2可看出在不同組合方案下,該算法的預(yù)測(cè)精度均較高,表明該算法具有較強(qiáng)的穩(wěn)定魯棒性、泛化適用性。

表2 不同組合下RUL預(yù)測(cè)結(jié)果對(duì)比Table 2 Comparison of RUL prediction results under different combinations

綜合考慮預(yù)測(cè)精度和時(shí)效性,將粒子數(shù)設(shè)為1 000,則在每一個(gè)預(yù)測(cè)時(shí)刻可以得到1 000個(gè)RUL預(yù)測(cè)值,從而得出RUL分布直方圖,部分預(yù)測(cè)時(shí)刻的RUL分布頻率直方圖如圖11所示。將粒子濾波預(yù)測(cè)中值當(dāng)做預(yù)測(cè)值,由圖11可知,在95%的置信區(qū)間內(nèi),預(yù)測(cè)值與真實(shí)值吻合程度較高。隨著閘門運(yùn)行年限的增加,在不同預(yù)測(cè)時(shí)刻可以得到預(yù)測(cè)的剩余壽命概率密度函數(shù),如圖12所示。從圖12可以看出,隨著閘門運(yùn)行年限的增加,剩余壽命的概率分布區(qū)間范圍在逐漸減小,說(shuō)明預(yù)測(cè)的不確定性在逐漸降低。

圖11 剩余壽命分布直方圖Fig.11 Histogram of RUL distribution

圖12 剩余壽命概率密度函數(shù)Fig.12 Probability density function of RUL

4.6 RUL預(yù)測(cè)結(jié)果及分析

為了驗(yàn)證所提算法的穩(wěn)定性以及預(yù)測(cè)性能的優(yōu)劣水平,對(duì)每一預(yù)測(cè)時(shí)刻進(jìn)行100次重復(fù)試驗(yàn)。同時(shí)為了更好地指導(dǎo)工程實(shí)際,給出了各個(gè)預(yù)測(cè)剩余壽命區(qū)間出現(xiàn)的概率,部分預(yù)測(cè)時(shí)刻結(jié)果如表3所示,表3中預(yù)測(cè)值是指100次循環(huán)下預(yù)測(cè)RUL的均值,Var是指100次循環(huán)下預(yù)測(cè)RUL的方差,AE代表預(yù)測(cè)值與實(shí)際RUL的絕對(duì)誤差。

表3 100次重復(fù)試驗(yàn)結(jié)果Table 3 Results of 100 replicate tests

從閘門運(yùn)行第10年開始預(yù)測(cè),按照上述方法可得到閘門不同運(yùn)行年限下的RUL預(yù)測(cè)結(jié)果,如圖13所示。由圖13可以看出,RUL預(yù)測(cè)結(jié)果與實(shí)際RUL結(jié)果相近,且真實(shí)RUL均在95%的預(yù)測(cè)RUL置信區(qū)間內(nèi)。

圖13 剩余壽命預(yù)測(cè)結(jié)果Fig.13 RUL prediction results

為了客觀衡量預(yù)測(cè)的精度,引入均方根誤差(Root Mean Square Error,RMSE)、平均絕對(duì)誤差(Mean Absolute Error,MAE)和方差絕對(duì)誤差(Variance Absolute Error,VAE)來(lái)進(jìn)行評(píng)價(jià)(式(21)),評(píng)價(jià)指標(biāo)得分越小說(shuō)明預(yù)測(cè)的精度越高,結(jié)果如表4所示。從表4可以看出3個(gè)評(píng)價(jià)指標(biāo)得分均比較小,說(shuō)明所提方法預(yù)測(cè)閘門RUL具有較高的精度。

表4 評(píng)價(jià)指標(biāo)得分Table 4 Scores of evaluation indices

5 結(jié) 論

(1)利用應(yīng)力、自振頻率、位移模態(tài)、應(yīng)變模態(tài)等多特征參數(shù)相較于單一指標(biāo)更能全面、準(zhǔn)確地反映閘門退化過(guò)程。

(2)利用非線性維納過(guò)程對(duì)閘門退化過(guò)程進(jìn)行建模,并基于粒子濾波對(duì)其進(jìn)行剩余壽命預(yù)測(cè),可以得到剩余壽命的概率分布以描述預(yù)測(cè)結(jié)果的不確定性,從而更好地指導(dǎo)閘門日常管理。同時(shí),由評(píng)價(jià)指標(biāo)RMSE為1.395 5,說(shuō)明所提方法具有較高的預(yù)測(cè)精度。

(3)未來(lái)根據(jù)剩余壽命預(yù)測(cè)結(jié)果,綜合考慮經(jīng)濟(jì)性、安全性等因素以確定合理的維修策略是下一步研究的關(guān)鍵。