數(shù)據(jù)驅(qū)動的設(shè)備時變可靠度評價及故障預(yù)測

林子昕，王少琦，安維中，別海燕

（中國海洋大學(xué)化學(xué)化工學(xué)院，山東青島266100）

化工行業(yè)復(fù)雜的生產(chǎn)工藝條件對設(shè)備的要求很高，設(shè)備是工業(yè)生產(chǎn)運行的基礎(chǔ)，設(shè)備運行的可靠性直接影響裝置的連續(xù)穩(wěn)定運行。由于設(shè)備的缺陷直接造成事故或者間接造成的事故會給企業(yè)造成巨大的損失，同時由于設(shè)備缺陷造成裝置停工，影響企業(yè)正常生產(chǎn)，帶來的經(jīng)濟損失也是巨大的。因此，保證化工生產(chǎn)裝置的長周期穩(wěn)定運行越來越重要，設(shè)備的可靠性必須得到重視。

可靠性評估是指根據(jù)設(shè)備在整個壽命周期內(nèi)的可靠性結(jié)構(gòu)、分布模型以及有關(guān)的可靠性信息，利用數(shù)理統(tǒng)計的方法對所評定產(chǎn)品的可靠性指標(biāo)給出估計的過程[1]。可靠性評估分為單元可靠性評估和系統(tǒng)可靠性評估，一般來說，“單元”是指作為單獨研究或單獨試驗對象的任何元器件、零件，甚至一臺完整的設(shè)備；而“系統(tǒng)”則包括了組成系統(tǒng)的各單元以及各單元間的可靠性結(jié)構(gòu)[2]。近年來，可靠性評估技術(shù)迅速發(fā)展，現(xiàn)有的單元可靠性評價方法可大致分為基于物理失效模型、數(shù)據(jù)驅(qū)動的失效統(tǒng)計模型和基于經(jīng)驗?zāi)Ｐ偷目煽啃栽u估方法。在實際應(yīng)用中，由于難以獲取復(fù)雜、高可靠性設(shè)備失效機理的物理模型，數(shù)據(jù)驅(qū)動的可靠性評估與故障預(yù)測方法成為近年來研究的主流，通過分析設(shè)備失效數(shù)據(jù)發(fā)現(xiàn)其失效分布類型，從而進行統(tǒng)計推斷。

相比于電子和微電子產(chǎn)品的可靠性分布模型近似服從于指數(shù)分布[3]，一些學(xué)者[4-5]通過研究指出，機械設(shè)備產(chǎn)品的可靠性分布模型大多數(shù)服從威布爾分布模型。文獻[6]采用威布爾分布與極大似然估計相結(jié)合的方法建立數(shù)學(xué)模型來間接評估尾軸承的可靠壽命；文獻[7]通過改進經(jīng)驗分布函數(shù)的建立途徑，優(yōu)化威布爾分布參數(shù)的求解過程，進而開發(fā)一種液壓支架可靠性評估的改進威布爾分布法。威布爾分布能夠很好地擬合不同類型失效率的產(chǎn)品壽命數(shù)據(jù)，在設(shè)備可靠性建模中應(yīng)用最廣泛。

在工程領(lǐng)域中，各種主客觀因素的影響，結(jié)構(gòu)的抗力性能在服役期間逐漸減弱，常把這種考慮結(jié)構(gòu)抗力隨時間變化的可靠度稱為時變可靠度[8]。近年來，時變可靠度的計算日益受到重視。文獻[9]考慮了一般大氣環(huán)境下銹蝕因素對鋼筋混凝土梁可靠度的影響，根據(jù)極限狀態(tài)方程建立時變可靠度模型。文獻[10]研究了飛機鎖緊機構(gòu)轉(zhuǎn)動關(guān)節(jié)的磨損問題，結(jié)果表明關(guān)節(jié)磨損隨時間的變化關(guān)系對鎖緊機構(gòu)的可靠性有重要影響。

由于化工設(shè)備的運行工況、載荷、應(yīng)力、強度等參數(shù)隨時間或載荷作用次數(shù)等壽命指標(biāo)在不斷變化，化工設(shè)備的可靠性是一個典型的時變過程，傳統(tǒng)的可靠性分析理論與方法不能較好地體現(xiàn)設(shè)備的動態(tài)特性，難以準(zhǔn)確地反映設(shè)備的可靠性隨壽命指標(biāo)的變化規(guī)律。本文在傳統(tǒng)威布爾（Weibull）分布可靠性評估模型基礎(chǔ)上引入隨時間變化的設(shè)備可靠性影響因素，建立多變量故障速率模型，對設(shè)備進行時變可靠度評價以及故障預(yù)測。

1 基于威布爾分布可靠性的評估模型

用于評估機械設(shè)備或零件可靠性的數(shù)學(xué)模型包括指數(shù)分布、正態(tài)分布、對數(shù)正態(tài)分布和威布爾分布，對比其他數(shù)學(xué)模型，威布爾分布兼容性好，對各種類型的數(shù)據(jù)擬合能力強，可以較全面地描述設(shè)備不同失效期的失效過程與特征。

1.1 威布爾分布數(shù)學(xué)模型

二參數(shù)威布爾分布模型的失效率函數(shù)按式(1)計算。

式中，t為故障時間，d（或h、s 等）；β為威布爾分布形狀參數(shù)；η為威布爾分布尺度參數(shù)。

在威布爾分布中最重要的參數(shù)是形狀參數(shù)，它可以決定概率密度曲線的基本形狀，尺度參數(shù)起放大或縮小曲線的作用，但不影響分布曲線的形狀。通過改變形狀參數(shù)可以表示不同階段設(shè)備的失效情況，同時，改變形狀參數(shù)的值威布爾分布也可以作為其他分布的近似，如可將形狀參數(shù)設(shè)為合適的值近似正態(tài)、對數(shù)正態(tài)、指數(shù)等分布。

1.2 參數(shù)估計

計算威布爾分布的模型參數(shù)有多種方法，目前常用的方法有極大似然估計法、矩估計法和最小二乘法估計法等。文獻[11]比較了以上各種計算方法的特點，結(jié)果表明極大似然估計法計算得出的數(shù)據(jù)更能滿足實際要求，因此本文用極大似然估計法來求解基于威布爾分布可靠性評估模型的形狀參數(shù)β和位置參數(shù)η。

設(shè)總體X～Weibull(β,η),X1,X2,…,Xn為來自總體的n個獨立的樣本，則其似然函數(shù)見式(3)。

將式(3)取對數(shù)，可得對數(shù)似然函數(shù)，見式(4)。

使用Newton-Raphson 迭代法計算上述非線性方程組，即可求出極大似然法估計的模型形狀參數(shù)β和尺度參數(shù)η，最終得到基于威布爾分布的可靠性評價模型。

2 時變可靠度評價模型

在化工行業(yè)中，即使一批結(jié)構(gòu)形式、尺寸、材料和加工方法上都相同的設(shè)備，它們的使用壽命卻各有長短，不盡相同，其原因是多方面的，工作環(huán)境如溫度和壓力的變化會對設(shè)備的壽命產(chǎn)生很大的影響。因此，在傳統(tǒng)威布爾分布的可靠度評價模型的基礎(chǔ)上進行改進，引入設(shè)備可靠性的影響因素，采用數(shù)據(jù)驅(qū)動的方法，建立一種基于威布爾分布的設(shè)備時變可靠度評價模型。

2.1 灰色關(guān)聯(lián)分析

灰色關(guān)聯(lián)分析是對一個系統(tǒng)發(fā)展變化態(tài)勢的定量描述和比較的方法，研究系統(tǒng)中特征因素與影響因素之間的關(guān)系，計算各影響因素序列與特征序列的相關(guān)程度，依照關(guān)聯(lián)度大小得出結(jié)論。化工設(shè)備生產(chǎn)運行中影響因素較多，本文通過灰色關(guān)聯(lián)分析計算關(guān)聯(lián)度來確定影響設(shè)備可靠性的主要因素。關(guān)聯(lián)系數(shù)計算公式見式(6)。

式(7)為數(shù)列xi對特征數(shù)據(jù)序列x0的關(guān)聯(lián)度。即把各個關(guān)聯(lián)系數(shù)集中為一個平均值，以此來反映整體關(guān)聯(lián)程度[12]。

2.2 響應(yīng)面法

近年來，元模型被廣泛用于通過使用代理模型預(yù)測系統(tǒng)響應(yīng)來降低可靠性分析的計算成本[13]，響應(yīng)面法等基于回歸的元模型代理模型得到了廣泛應(yīng)用。響應(yīng)面法（Response Surface Method）是利用數(shù)學(xué)和統(tǒng)計學(xué)技術(shù)解決復(fù)雜系統(tǒng)的輸入（隨機變量）與輸出（系統(tǒng)響應(yīng)）關(guān)系的方法。

響應(yīng)面法分為多項式響應(yīng)面法和神經(jīng)網(wǎng)絡(luò)響應(yīng)面法，神經(jīng)網(wǎng)絡(luò)響應(yīng)面法計算過程十分復(fù)雜，工程上的計算多采用多項式響應(yīng)面法。多項式響應(yīng)模型見式(8)。

式中，ai,k，bi,j,k,n(i=1,2,...,n)為表達式的待定系數(shù)。

從響應(yīng)面函數(shù)表達式可看出，若隨機變量個數(shù)取N，共有(N+1)(N+2)/2 個待定系數(shù)。綜合響應(yīng)式函數(shù)的計算量、計算精度和計算穩(wěn)定性等方面因素，響應(yīng)面函數(shù)最常用的是二次多項式形式，式(8)中K= 2即為二次多項式響應(yīng)模型。

通過可決系數(shù)R2對響應(yīng)關(guān)系式進行檢驗，見式(9)。

一些時變可靠性分析方法采用替代模型來逼近復(fù)雜系統(tǒng)的隱式狀態(tài)函數(shù)[14]，本文采用響應(yīng)面法建立失效率模型中形狀參數(shù)β和尺度參數(shù)η與設(shè)備運行影響因素的關(guān)系，將影響因素引入設(shè)備失效率計算模型，建立多變量失效速率模型。

2.3 時變可靠度模型

在傳統(tǒng)威布爾分布的可靠性評價模型中，用極大似然法可計算得到設(shè)備的形狀參數(shù)β和尺度參數(shù)η，根據(jù)響應(yīng)面法可計算得到模型參數(shù)與影響因素的響應(yīng)關(guān)系式，分別見式(10)和式(11)。

將式(10)和式(11)代入威布爾分布的失效率模型可得出故障速率模型，見式(12)。

計算時變可靠度需考慮設(shè)備在不同時間段內(nèi)的運行工況，運行條件不同，故障速率模型也不同，設(shè)備時變可靠度計算模型如下。

若t∈[0,t1],可直接計算可靠度，見式(13)。

3 應(yīng)用實例

在化工生產(chǎn)過程中，機泵設(shè)備始終都發(fā)揮著至關(guān)重要的作用，以某公司催化裂化裝置中的12 類機泵運行數(shù)據(jù)為例，對設(shè)備時變可靠度進行評價。根據(jù)機泵的歷史故障時間數(shù)據(jù)，采用極大似然法計算失效率函數(shù)中的形狀參數(shù)β和尺度參數(shù)η以及平均故障間隔時間，如表1所示。

機泵運行過程中影響因素較多，如泵的揚程，泵中流體的流量，流體的溫度、密度，泵的進壓、出壓以及泵的功率等，機泵運行狀況見表2。采用灰色關(guān)聯(lián)分析法比較影響因素與模型參數(shù)β和η之間關(guān)聯(lián)度的大小，選擇關(guān)聯(lián)度較大的因素為模型的主要影響因素。

以機泵的形狀參數(shù)β和尺度參數(shù)η為特征參考數(shù)列，分別計算其與泵的進壓、泵的出壓、泵的功率、泵的揚程、進出泵的流體壓差、流體溫度和流體流量等影響因素數(shù)列的關(guān)聯(lián)度。

圖1 為形狀參數(shù)β與各影響因素的關(guān)聯(lián)系數(shù)圖。圖2 為尺度參數(shù)η與各影響因素的關(guān)聯(lián)系數(shù)圖。

從圖中可以看出，除了泵的進壓與模型參數(shù)的關(guān)聯(lián)度較低，其他影響因素與模型參數(shù)的關(guān)聯(lián)度均較高，通過計算平均關(guān)聯(lián)度，得到影響因素與模型參數(shù)的關(guān)聯(lián)度。模型參數(shù)與各影響因素關(guān)聯(lián)度大小順序為密度ρ＞出壓p2＞揚程H＞功率W＞流量V＞溫度T＞進壓p1。

考慮到12 類機泵數(shù)據(jù)有限，選擇泵中流體的密度、泵的出壓p2和揚程H，作為機泵運行中的主要影響因素，采用響應(yīng)面法建立選擇的影響因素數(shù)據(jù)與設(shè)備形狀參數(shù)和尺度參數(shù)之間的響應(yīng)關(guān)系式。響應(yīng)關(guān)系式中變量個數(shù)為3，需要求解的二次響應(yīng)式的響應(yīng)系數(shù)個數(shù)是10，在響應(yīng)關(guān)系式中，aij為變量一次項和平方項的系數(shù)；bij為變量交叉項的系數(shù)，c為常數(shù)項。在Matlab中，計算得到響應(yīng)系數(shù)，通過Matlab擬合的響應(yīng)關(guān)系式計算的數(shù)據(jù)與設(shè)備形狀參數(shù)和尺度參數(shù)數(shù)據(jù)的對比如圖3和圖4所示。

圖3 響應(yīng)關(guān)系式擬合的形狀參數(shù)與原始數(shù)據(jù)對比圖

圖4 響應(yīng)關(guān)系式擬合的尺度參數(shù)與原始數(shù)據(jù)對比圖

由圖可知，擬合數(shù)據(jù)曲線與原始數(shù)據(jù)曲線基本一致，分別計算曲線的可決系數(shù)R2，分別為R2（β響應(yīng)式）=0.9862、R2（η響應(yīng)式）=0.9037。可決系數(shù)均大于0.90，響應(yīng)關(guān)系式擬合效果較好。

根據(jù)響應(yīng)系數(shù)，得出威布爾分布形狀參數(shù)β和尺度參數(shù)η與影響因素（出壓、揚程、壓差）的關(guān)系式，見式(19)和式(20)。

在失效率函數(shù)中，引入設(shè)備運行的影響因素，將所求得的關(guān)系式代入設(shè)備失效率模型中，可得出故障速率模型和時變可靠度模型，分別見式(21)和式(22)。

若機泵的運行條件隨時間發(fā)生變化，可根據(jù)上述兩個模型計算設(shè)備時變可靠度。表3為機泵1的運行條件，根據(jù)所給條件計算機泵1的時變可靠度和MTBF。假設(shè)機泵在0～100d運行參數(shù)不變，其中密度為1000g/cm3，出壓為6.2MPa，揚程為600m，在100～200d 和200～300d 時，由于物性參數(shù)和操作參數(shù)發(fā)生變化，3 種影響因素均有一定程度的增大，計算其模型參數(shù)，時變可靠度和平均剩余壽命如表4所示。

表3 機泵1運行條件

表4 機泵的可靠性計算

從表4 中可以看出，機泵運行200～300d 后，時變可靠度為0.3373，此時機泵的可靠度較低，在200～300d 時剩余壽命為負(fù)值，說明此時機泵可能已達到故障狀態(tài)。

若不考慮機泵運行參數(shù)的變化，在泵中流體密度1000g/cm3、泵的出壓6.2MPa、揚程600m操作條件下，運行300d后，計算機泵的可靠度為0.4757，在運行工況近似的情況下，時變可靠度計算結(jié)果為0.3373，證明了時變可靠度計算方法的可行性。

4 結(jié)論

本文在傳統(tǒng)威布爾分布的基礎(chǔ)上，改進模型參數(shù)的計算方式，建立了數(shù)據(jù)驅(qū)動的設(shè)備時變可靠度評價模型，并對設(shè)備的故障進行預(yù)測。可用于計算新設(shè)備可靠度和MTBF的計算以及運行工況隨時間發(fā)生變化的設(shè)備時變可靠度的計算和故障預(yù)測。在運行工況發(fā)生改變的情況下，設(shè)備的時變可靠度計算模型可對設(shè)備故障時間進行預(yù)測。

符號說明

H—— 揚程，m

p1—— 泵的進壓，MPa

p2—— 泵的出壓，MPa

Δp—— 泵的壓差，MPa

R—— 可靠度

T—— 流體溫度，°C

t—— 設(shè)備故障時間，d

V—— 體積流量，L/h

W—— 泵的功率，W

λ—— 失效率，%

β—— 威布爾分布形狀參數(shù)

η—— 威布爾分布尺度參數(shù)，d