屈持,王海清,姜巍巍,孫浩,張景康
(1 中國石油大學(華東) 機電工程學院,安全科學與工程系,山東青島266580; 2 中國石化青島安全工程研究院,山東青島266071; 3 中海油安全技術服務有限公司,天津300450)
在對可修系統的故障數據進行可靠性分析時,隨機點過程中的“修復如新”和“修復如舊”模型被視為可修系統評估模型的兩種特殊情況[1-5],而對于大多數安全關鍵設備來說[6-12],不完全維修更符合工程實際。不完全維修是可修系統在經歷維修活動后,系統改善效果介于“修復如新”和“修復如舊”之間的一種維修模型,通過量化維修活動對設備可靠性的影響,實現了對安全關鍵設備的可靠性分析,這有助于保障生產裝置的平穩運行。
目前常用的不完全維修模型主要有兩種:廣義更新過程和比例強度模型,其中廣義更新過程由Kijima 等[13]提出,該模型在“更新過程”和“非齊次泊松過程”的基礎上,引入役齡回退因子的概念,通過系統的役齡減小程度來刻畫維修效果;Cox[14]提出了比例強度模型,其通過故障強度減小系數來衡量系統的維修效果,但由于該模型的非連續性,很難獲得諸如故障強度、瞬時故障間隔時間等可靠性指標封閉形式的解,導致其難以在實際問題中運用。Syamsundar 等[15]考慮將最小維修模型和比例強度模型結合起來對工業系統進行可靠性評估,但存在參數求解復雜的問題;張根保等[16]和Wang 等[17-18]在不完全維修條件下分別利用廣義比例強度過程和對數線性過程對多臺數控機床的故障數據進行可靠性評估,然而并未考慮多個故障數據間差異對可靠性評估的影響。
鑒于此,本文針對不完全維修活動中的可修系統,提出了一種混合Kijima Ⅰ虛擬役齡模型,量化了維修活動對設備有效役齡的影響,考慮了不同故障原因下設備故障數據差異性導致可靠性評估結果偏差大的影響,并應用非線性約束優化的方法解決了模型參數難以求解的問題,最后將文中模型應用到某LNG 接收站卸船系統安全聯鎖回路中的壓力變送器進行驗證,將該模型下的可靠性指標與混合樣本下的可靠性指標進行對比,以表明混合Kijima Ⅰ模型的優勢與數值性能。
Kijima Ⅰ虛擬役齡模型是不完全維修的重要類型,該模型假設維修只作用于當前故障間隔時間,不能減少設備維修之前累積的磨損[19]。若單個可修設備的第k 個故障間隔時間為xk= tk- tk-1,k=1,2,…,n 且t0= 0,在經歷t1時刻的維修活動后,其壽命狀態與e1時刻的壽命狀態等效,e1<t1,這里e1為經過t1維修時刻的等效役齡。此外,每次維修活動后的等效役齡ek不僅受到之前等效役齡的影響,還需考慮當前故障間隔時間的等效役齡,即將第k 個故障間隔時間xk減少到qkxk,當設備在tk時刻發生不完全維修后,其等效役齡表示為

式中,qk表示第k 個維修活動的役齡回退因子,該值量化了設備的維修程度,本文假設每次維修活動的役齡回退因子均為定值,即qk= q,因此有

此時,設備的運行時間與等效役齡的關系如圖1所示。
在設備第k - 1 次維修后,設備的第k 個故障間隔時間xk的累積條件分布函數為

圖1 設備運行時間與虛擬役齡關系Fig.1 Relationship between equipment operating time and virtual service age

由于在各種故障過程模型中的靈活性和有效性,Weibull 分布被廣泛應用于可靠性與壽命研究中,因此本文選取Weibull 分布作為潛在故障分布[20]。單臺設備的故障間隔時間服從兩參數Weibull分布,其累積分布函數為

其相應的概率密度函數為

式中,α為形狀參數,θ-1/α為Weibull模型的尺度參數。
所以設備在經歷不完全維修活動后,Kijima Ⅰ虛擬役齡模型下的Weibull累積條件分布函數為

其相應的條件概率密度函數為

所以定時截尾下的故障間隔時間似然函數為

其對數似然函數為

一般情況下,分別對式(9)中的α、θ求導并令其為零,可得對應的方程組,通過解方程的方法求得參數的極大似然估計值。然而上述方程組沒有封閉形式的解,而且難以通過數值解法求得參數估計值,因此可利用非線性約束規劃的方法進行求解[21],以對數似然函數值最大為目標函數,結合參數的約束條件,求得該模型的參數估計值。
假設第r 臺設備在第j 個統計區間[0,Tj]內發生nj個故障,其故障發生時刻分別為tk,j(k =1,2,…,nj; j = 1,2,…,r),令xk,j表示故障時刻tk-1,j與tk,j之間的故障間隔時間,即xk,j= tk,j- tk-1,j,若在周期性監測結束時,設備并沒有發生故障,此時最后一次故障間隔時間為Tj- tk,j。
在對可修系統的故障數據進行可靠性分析時,需先通過故障數據趨勢檢驗確定出分析數據所用的模型。由于累積故障強度函數對時間圖檢驗法的適用性[22],在故障數據的可靠性分析中得到廣泛使用,其表達式如式(10)所示[23]。

式中,dk,j表示在tk,j時仍在工作的設備數;若tk,j為時間截尾數據,則dk= dk-1- 1,反之則dk= dk-1。
對該方法的散點圖即(tk,j,S(tk,j))而言,當散點近似于一條直線時,表示設備故障數據沒有趨勢,上凸表示設備故障強度減小,下凹則表示設備故障強度增大。
赤池信息準則(AIC)和貝葉斯信息準則(BIC)作為可靠性模型選擇的常用方法之一[24],其數值越小表明模型選擇越精確,表達式為[25]:

式中,p 為所選模型參數的個數,v 為故障數據的個數,maxlnL為故障數據的極大似然函數。
安全關鍵設備在運行過程中會由于各種因素發生不同故障模式的情況,針對相似工況下多樣本的設備壽命數據分析[26-27],經典的統計方法并未考慮不同故障模式下壽命數據存在差異的影響,而混合分布模型主要用于非同質產品按照一定比例混合成總體的建模過程,適用于設備存在一種或多種故障模式的情況。因此,在多樣本的安全關鍵設備壽命數據下,本文提出一種混合Kijima Ⅰ模型,其累積條件分布函數為

對于來自相似工況下的安全關鍵設備壽命時間認為屬于同一分布,但由于不同故障原因導致分布參數不同,因而可用1~9 標度法從設備故障頻率和故障后果兩個方面判定出各故障模式的重要度,進而求得各混合參數,表1 給出了確定混合參數的評分準則。

表1 1~9標度評分準則Table 1 1—9 scale scoring criteria
通常情況下可靠性指標可以給出系統變化的定量描述,其中平均故障間隔時間(MTBF)和指定時間下的可靠度(R)是安全關鍵設備可靠性的重要評估指標,針對壓力變送器的故障數據,在獲得混合Kijima Ⅰ模型的參數估計值后,就可以對可靠性指標進行計算。
可靠壽命的表達式為:

設備某一時刻的可靠度為:

基于混合Kijima Ⅰ模型的可靠性評估流程如圖2所示。
液化天然氣(liquid natural gas,LNG)接收站的功能是將產地輸出的LNG 經專用運輸船運輸到接收站并在此進行氣化外輸[28],其工藝流程主要包括卸船、儲存、低高壓運輸、再冷凝、火炬放空、氣化外輸工藝和輔助設施(排污、后臺管理工藝等),圖3為某LNG接收站的工藝流程。

圖2 可靠性評估流程Fig.2 Flow chart of reliability assessment
作為LNG 接收站首要工藝的卸船流程,在整個接收站流程中占有較大比重[29]。卸船時,先連接卸料臂,然后進行氮氣置換,LNG 經船上的輸送泵,經卸料臂及其支管匯集到總管,并通過總管輸送到LNG 儲罐中;儲罐中的蒸發氣(BOG)通過蒸發氣回流臂返回到LNG 船艙,從而維持船艙內的氣相壓力平衡;卸料完成后,利用N2分液罐內的氮氣吹掃殘留在卸料臂、卸料管線上的LNG。在上述操作過程中,要注意觀察各工藝流程中壓力變送器的變化情況,一旦出現儀表故障或壓力波動異常的工況,應立即啟動泄壓系統進行降壓,保障工藝流程的正常運行[30]。在卸料過程中儲罐產生的蒸發氣會引起一定的壓力變化,此外隨著卸料工藝的發展,運輸船、卸料管線、BOG 回流管線等都會發生壓力波動,因此系統地分析卸船工藝中的壓力變送器對保障整個卸船裝置的平穩運行具有現實意義[31],圖4 為LNG卸船流程中的部分壓力變送器分布。
隨著自動化領域數字化、網絡化、智能化的快速發展,智能儀表面臨著功能安全失效因素增多和信息安全攻擊等威脅。壓力變送器是卸船過程中檢測、控制、管理功能實現的重要設備,通過分析卸料管線安全聯鎖回路中的壓力變送器來驗證該算法的有效性和準確性。由于LNG 輸送泵和外輸管線的操作壓力都比較高,當壓力變送器發生故障不能實時傳送壓力信號時,很容易導致泄漏事故發生。通過查找相關標準[32-33],設置形狀參數和尺度參數的初值分別為1.2 和125000,假設役齡回退因子q 為0.15,利用Monte Carlo 模擬獲得三臺壓力變送器同一工況下的故障樣本數據,其數據預處理的結果如表2 所示,其中140000+表示該設備運行至140000 h 時并未發生故障,累積強度函數對時間圖如圖5所示,根據散點圖的凹凸情況,表明當前壓力變送器的故障強度在減小,可應用Kijima Ⅰ虛擬役齡模型進行分析。

圖3 LNG接收站工藝流程Fig.3 Process flow chart of LNG receiving station

圖4 LNG卸船工藝部分流程Fig.4 Partial flowchart of LNG unloading process
通過對比完全維修模型、最小維修模型[34]和Kijima Ⅰ虛擬役齡模型下的參數估計結果,見表3,表明在系統可修復的條件下,Kijima Ⅰ虛擬役齡模型優于Weibull分布和NHPP過程,且更符合實際。
為進一步驗證Kijima Ⅰ模型的有效性,根據式(6)和表3中的參數估計初值,利用Monte Carlo 仿真方法[35-36]得到多個故障間隔時間的抽樣公式,并利用Kijima Ⅰ模型進行驗證計算。不同樣本下的評估結果如表4 所示,可以看出,隨著樣本量的增大,Monte Carlo 仿真數據的參數估計值越來越接近真值,其相對誤差也隨之減小,顯示了該模型的有效性。

式中,Uk為(0,1)均勻分布上的隨機數,所以故障時間的抽樣公式為


表2 模擬故障數據及預處理結果Table 2 Simulated fault data and preprocessing

表3 模型對比結果Table 3 Model comparison results

圖5 累積故障強度函數對時間圖Fig.5 Cumulative failure intensity function versus time
由于卸船過程中的儲罐壓力波動較大,且一般最大壓力值不能超過25 kPa,因此保證壓力變送器的可靠性就顯得尤為重要。結合國內某液化天然氣公司提供的LNG接收站十年間的事故記錄,收集到一批相似工況下的單晶硅諧振式壓力變送器危險失效的故障數據,具體故障原因及故障時間如表5所示。
從表5 中可知,卸船系統中的壓力變送器故障數據從故障原因上可以分為輸出不穩定故障樣本、輸出偏差大故障樣本以及惡意網絡攻擊故障樣本,這三類樣本都屬于壓力變送器的危險失效類別,因此可用混合Kijima Ⅰ模型對不完全維修活動下的壓力變送器進行分析。隨后采用1~9標度法得出壓力變送器故障原因的判斷矩陣,對其歸一化處理后,得到通過一致性檢驗的權重矩陣,壓力變送器的混合參數以及各分部參數如表6所示。

表4 不同樣本的評估結果對比Table 4 Comparison of evaluation results of different samples

表5 壓力變送器樣本數據Table 5 Pressure transmitter sample data
根據式(14)以及表6 中的數據,可得該卸船系統中的壓力變送器的混合累積條件分布函數為

表6 混合虛擬役齡模型分布參數Table 6 Distribution parameters of mixed virtual service age model

為進一步說明混合分布模型的有效性,現忽略壓力變送器發生危險失效的原因,將表5 中的故障樣本合并為一個樣本進行分析,根據式(13),可得混合分布與混合樣本各自對應的中位壽命:t(0.5) =36527 h,t(0.5)′ = 27854 h。若以壓力變送器的中位壽命作為其進行不完全維修活動的節點,則混合樣本會以壓力變送器信息安全失效為主要原因增加企業的成本支出,不符合實際情況;以混合分布模型得到的中位壽命既考慮了壓力變送器易發生的功能安全失效,又兼顧了信息安全失效的特殊情況,相對可以延長企業的檢維修周期,減小企業的安全成本支出。
圖6進一步給出了混合Kijima Ⅰ模型分別與各類別樣本以及混合樣本分布的可靠度對比曲線,可以看出混合分布很好地擬合了三類樣本的分布,綜合衡量了壓力變送器三種故障原因的嚴重度和發生頻率。對于壓力變送器的校驗,國家計量監督部門規定石化行業壓力變送器的檢驗周期為一年,若以一年為時間節點,分別得到壓力變送器各類別樣本對應的可靠度:R1(8760) = 0.83,R2(8760) =0.8618,R3(8760) = 0.9599。若以樣本3 為主,運行一年時的可靠度很高,主要原因在于網絡安全對工控系統攻擊的頻率不高,若在此時檢維修會造成企業安全經費的過渡投入;考慮到壓力變送器常見的故障原因,若以樣本1 和樣本2 為主,一年時的可靠度較低,若根據企業的大修計劃(3 年)對壓力變送器進行檢修,現場可能會提前發生不安全事故,造成財產損失或者人員傷亡。因此,企業應根據壓力變送器不同的運行環境,合理調整檢維修計劃,保證設備的機械完整性水平。
(1)針對傳統可靠性評估模型未考慮維修活動對設備可靠性影響的問題,本文在不完全維修的條件下,結合設備的故障趨勢,提出了一種符合工程實際的Kijima Ⅰ虛擬役齡模型,量化了維修活動對設備實際役齡的恢復效果,應用非線性約束規劃法解決了模型參數求解困難的問題,進而實現了對石化安全關鍵設備的可靠性評估,有助于保證設備以較高的可靠度運行。
(2)在多樣本故障數據的設備可靠性評估中,考慮到功能故障和網絡攻擊等不同故障原因導致設備故障數據間存在差異性,本文進一步提出了混合Kijima Ⅰ虛擬役齡模型,解決了故障數據間差異性導致評估結果偏差大的問題,對比混合樣本的可靠性評估指標,文中模型的可靠性評估結果合理性更高,并為設備的差異化維修提供了一定的理論指導,避免了企業安全成本的過渡投入。

圖6 混合Kijima Ⅰ模型與各樣本的可靠度曲線Fig.6 Reliability curve of mixed Kijima Ⅰmodel and each sample