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

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

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

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

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

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

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

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

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

1.1 威布爾分布數學模型

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

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

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

1.2 參數估計

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

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

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

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

2 時變可靠度評價模型

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

2.1 灰色關聯分析

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

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

2.2 響應面法

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

響應面法分為多項式響應面法和神經網絡響應面法，神經網絡響應面法計算過程十分復雜，工程上的計算多采用多項式響應面法。多項式響應模型見式(8)。

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

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

通過可決系數R2對響應關系式進行檢驗，見式(9)。

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

2.3 時變可靠度模型

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

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

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

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

3 應用實例

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

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

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

圖1 為形狀參數β與各影響因素的關聯系數圖。圖2 為尺度參數η與各影響因素的關聯系數圖。

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

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

圖3 響應關系式擬合的形狀參數與原始數據對比圖

圖4 響應關系式擬合的尺度參數與原始數據對比圖

由圖可知，擬合數據曲線與原始數據曲線基本一致，分別計算曲線的可決系數R2，分別為R2（β響應式）=0.9862、R2（η響應式）=0.9037?？蓻Q系數均大于0.90，響應關系式擬合效果較好。

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

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

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

表3 機泵1運行條件

表4 機泵的可靠性計算

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

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

4 結論

本文在傳統威布爾分布的基礎上，改進模型參數的計算方式，建立了數據驅動的設備時變可靠度評價模型，并對設備的故障進行預測?？捎糜谟嬎阈略O備可靠度和MTBF的計算以及運行工況隨時間發(fā)生變化的設備時變可靠度的計算和故障預測。在運行工況發(fā)生改變的情況下，設備的時變可靠度計算模型可對設備故障時間進行預測。

符號說明

H—— 揚程，m

p1—— 泵的進壓，MPa

p2—— 泵的出壓，MPa

Δp—— 泵的壓差，MPa

R—— 可靠度

T—— 流體溫度，°C

t—— 設備故障時間，d

V—— 體積流量，L/h

W—— 泵的功率，W

λ—— 失效率，%

β—— 威布爾分布形狀參數

η—— 威布爾分布尺度參數，d