張 杰,張佩穎,虞維超,刁 逢,王圣潔,宮 敬
(1.中國石油大學(北京) 油氣管道輸送安全國家工程實驗室,北京 102249;2.中國石油天然氣股份有限公司北京油氣調(diào)控中心,北京 100007)
管道輸送是油氣資源輸送最有效的方式,管道失效數(shù)據(jù)統(tǒng)計結果表明,管道腐蝕失效占管道失效總量的近30%,是造成管道失效的主要原因之一[1]。為了提高管道的安全運行能力,國內(nèi)外學者對腐蝕管道的結構可靠性進行了大量研究[2-4]。腐蝕管道結構可靠性評價是評估腐蝕對管道完整性破壞程度的一種有效方法,生產(chǎn)和運行過程中所涉及的不確定性都通過概率方法納入其中[5-7]。腐蝕管道結構可靠性分析通?;诤g缺陷管道的極限狀態(tài)函數(shù),通過建立管道腐蝕缺陷生長模型,采用蒙特卡羅等方法計算管道的失效概率,進而評價管道在全生命周期內(nèi)的可靠性[8-10]。
建立適當?shù)墓艿栏g缺陷生長模型對腐蝕管道結構可靠性評價的準確性和有效性至關重要:一方面,過高估計管道腐蝕速率會導致對管道進行不必要的檢測和維護,造成管道管理成本過高;另一方面,過低估計管道腐蝕速率則會由于未能及時對管道進行預防性維護而導致管道失效,造成嚴重的后果[11]。管道腐蝕缺陷的生長過程本質(zhì)上是隨機的,包括時間和空間的變異性。時間變異性意味著管道腐蝕缺陷的增長速率會隨時間變化,空間變異性意味著不同管道腐蝕缺陷的增長速率不同。文獻[12]采用隨機過程的方法建立了管道腐蝕缺陷的概率增長模型。
由于設計、制造和運行工況的原因,管道的結構參數(shù)、缺陷的尺寸參數(shù)之間相互關聯(lián),多個腐蝕缺陷的增長相互依賴,在統(tǒng)計意義上均表現(xiàn)出一定程度的相關性[13-16]。張鵬等[14]研究了單個缺陷的尺寸參數(shù)之間的相關性對腐蝕管道結構可靠性的影響,但未考慮多個腐蝕缺陷生長過程的相關性;Li等[15]將缺陷之間的相關性視為概率交集進行研究,提出了一種管道腐蝕缺陷的概率增長計算方法,但沒有模擬多個腐蝕缺陷相互依賴生長的物理過程。此外,管道失效數(shù)據(jù)是管道結構退化到臨界值的直接體現(xiàn),蘊含了寶貴的管道完整性信息。為了準確地反映腐蝕管道的安全水平,腐蝕管道結構可靠性評價中涉及的參數(shù)不確定性、隨機變量之間的相關性以及失效數(shù)據(jù)均需要加以考慮。然而,在現(xiàn)有的腐蝕天然氣管道結構可靠性評價中,鮮有研究全面考慮隨機變量之間相關性的影響,更沒有涉及到管道失效數(shù)據(jù)對腐蝕管道結構可靠性評價結果的更新和校正?;诖?,本文建立了考慮失效數(shù)據(jù)和隨機變量之間相關性的腐蝕管道可靠性評價方法,以完善目前的研究,為提高腐蝕管道的安全管理提供理論基礎。
腐蝕管道結構的可靠性評價是一個典型的預測問題,其利用管道在線檢測得到的有限數(shù)據(jù),評估未來一段時期內(nèi)管道結構的可靠性,并預測管道系統(tǒng)達到極限狀態(tài)的時間,從而制定相應的預防性維修策略。
管道腐蝕缺陷存在兩種極限狀態(tài):小孔泄漏失效和爆裂失效。對于一個給定的管道腐蝕缺陷,給出管道在時間t時對應小孔泄漏失效和爆裂失效的極限狀態(tài)方程g1(t)和g2(t)分別如下[17-18]:
g1(t)=0.8wt-dmax(t)
(1)
g2(t)=pb(t)-pop(t)
(2)
(3)
式中:0.8為安全系數(shù);wt為管道壁厚(mm);t為上一次在線檢測后經(jīng)過的時間(a);dmax(t)為時間t時管道腐蝕缺陷的最大深度(mm);pop為管道運行壓力(MPa);pb(t)為腐蝕缺陷處的管道爆破壓力(MPa);χ為模型誤差;σu為抗拉強度(MPa);D為管道直徑(mm);l(t)為時間t時管道腐蝕缺陷的軸向長度(mm)。
g1(t)≤0和g2(t)≤0分別代表含腐蝕缺陷管道小孔泄漏失效和爆裂失效的發(fā)生。含腐蝕缺陷管道小孔泄漏失效和爆裂失效對應的數(shù)學表達式如下[19]:
Psl(t)=Prob[(g1(t)≤0)∩(g2(t)>0)]
(4)
Pbu(t)=Prob[(g1(t)>0)∩(g2(t)≤0)]
(5)
式中:Psl和Pbu分別為含腐蝕缺陷管道的小孔泄漏失效概率和爆裂失效概率。
在給定的管道腐蝕缺陷處發(fā)生小孔泄漏和爆裂失效被認為是互斥事件[20],因此含腐蝕缺陷管道的總失效概率Pf為發(fā)生小孔泄漏和爆裂失效的概率之和:
Pf(t)=Psl(t)+Pbu(t)
(6)
腐蝕缺陷深度是影響腐蝕管道結構可靠性評估的一個關鍵參數(shù),而腐蝕缺陷長度的增長對腐蝕管道結構可靠性評價結果的影響很小[20],故假設管道腐蝕缺陷長度不隨時間增長。本文采用線性隨機變量模型來描述管道腐蝕缺陷深度的生長,將多個腐蝕缺陷的年深度增量量化為具有相關性且服從特定分布的多維隨機變量:
di,max(t)=di,max(0)+Xi(t)
(7)
其中:
(8)
式中:di,max(t)為管道腐蝕缺陷i在時間t時的最大深度(mm);di,max(0)為管道腐蝕缺陷i的初始最大深度(mm),通過內(nèi)檢測數(shù)據(jù)來確定;Xi(t)為管道腐蝕缺陷i的總深度增量(mm);ΔXi,j(Δt)表示管道腐蝕缺陷i在第j年的年深度增量,即深度腐蝕速率(mm/a);Δt為模擬時間步長,取1 a;T為腐蝕管道結構可靠性的評價周期(a)。
腐蝕管道結構可靠性評價涉及的隨機變量之間的相關性可分為兩個層次:第一層為管道結構參數(shù)、缺陷尺寸參數(shù)等隨機變量之間的相關性;第二層為多個腐蝕缺陷生長過程的相關性。前者是一個二維隨機變量之間的相關性問題,而后者則是與腐蝕缺陷數(shù)量有關的多維隨機變量之間的相關性問題。具體來說,本文考慮了管徑與壁厚、腐蝕缺陷初始長度與初始深度以及多個腐蝕缺陷年深度增量之間的相關性,它們都可以采用具有特定相關系數(shù)矩陣的Copula函數(shù)來描述和解決。
Copula函數(shù)是一個包含n(n≥2)個標準均勻分布變量Ui(i=1,2,…,n)的聯(lián)合分布函數(shù)[21]。本文采用高斯Copula函數(shù)來描述本研究所涉及的隨機變量之間的相關性:
C(u1,u2,…,un)=Φn(Φ-1(u1),Φ-1(u2),…,Φ-1(un);R)
(9)
式中:C(u1,u2,…,un)表示Copula函數(shù),其中ui是Ui的值(i=1,2,…,n);Φn(·;R)表示具有相關系數(shù)矩陣R(n×n)的n維標準正態(tài)分布函數(shù),其中Φ-1(·)是一維標準正態(tài)分布函數(shù)的逆函數(shù);R的非對角元素rij(i,j=1,2,…,n;i≠j)代表Φ-1(ui)與Φ-1(uj)之間的線性相關系數(shù)。
利用高斯Copula函數(shù)可以非常便捷地生成具有特定相關系數(shù)的n(n≥2)維隨機變量樣本,以本研究為例,首先通過具有相關系數(shù)矩陣R的高斯Copula函數(shù)Φn(·);R)生成相關樣本ui(i=1,2,…,n);然后利用等概率密度邊緣變換方法[見公式(10)],可生成具有特定分布的相關隨機變量樣本xi(i=1,2,…,n)。
Fi(xi)=ui(i=1,2,…,n)
(10)
式中:Fi(·)表示累積分布函數(shù)。
圖1 含n個腐蝕缺陷的管段失效概率模擬計算流程圖Fig.1 Flow chart for simulating calculation of failure probability of the pipeline segment with n corrosion defects注:pi,b(t)為腐蝕缺陷i處的管道爆破壓力(MPa);Psl,pipe、Pbu,pipe和Pf,pipe分別表示含腐蝕缺陷管道的小孔泄漏失效概率、爆裂失效概率和總失效概率。
貝葉斯更新的對象為概率分布的參數(shù)值,故先基于秩回歸方法確定失效概率的最佳分布,再通過貝葉斯方法更新腐蝕管段失效概率分布函數(shù)的參數(shù)。本文選取5種候選分布來描述腐蝕管段的失效概率,并從中得到最佳分布。具體步驟如下[1]:
(1) 選擇指數(shù)分布、正態(tài)分布、對數(shù)正態(tài)分布、威布爾分布和耿貝爾分布作為腐蝕管道失效概率的候選分布,每種候選分布的累積分布函數(shù)和處理后的線性方程列于表1。
表1 描述腐蝕管道失效概率候選分布的累積分布函數(shù)和線性方程Table 1 Cumulative distribution functions and linear equations of candidate distributions describing the failure probability of corroded pipelines
(2) 將腐蝕管段逐年的失效概率{t,P(t)}代入各累積分布函數(shù)的線性方程中,采用秩回歸方法獲得線性方程系數(shù)A和B的值,進一步計算得到各累積分布函數(shù)的相應參數(shù)值。
(3) 采用下式計算各線性方程的相關系數(shù),相關系數(shù)絕對值最大的分布則為描述腐蝕管道失效概率的最佳分布。采用其累積分布函數(shù)F(t|θ)表示腐蝕管段的失效概率,其中θ為上述步驟(2)得到的該累積分布函數(shù)的相應參數(shù)值。線性方程相關系數(shù)ρ的計算公式為
(11)
在貝葉斯方法中,失效概率分布函數(shù)的參數(shù)θ被視為隨機變量,并具有先驗分布函數(shù)p(θ)。更新后的參數(shù)θ的分布稱為后驗分布,其概率密度函數(shù)可根據(jù)Bayes定理[23]計算如下:
L(tf1,tf2,…,tfn|θ)p(θ)
(12)
其中:
(13)
(14)
式中:∝表示正比關系;tfi(i=1,2,…,n)為管道失效數(shù)據(jù),即管道由于腐蝕缺陷而失效的時間;L(tf1,tf2,…,tfi|θ)為管道失效數(shù)據(jù)的抽樣密度函數(shù);m(tf1,tf2,…,tfi)為管道失效數(shù)據(jù)的邊緣密度函數(shù),其不依賴于參數(shù)θ;p(θ|tf1,tf2,…,tfi)為參數(shù)θ的后驗分布函數(shù)。
為了便于描述,將參數(shù)θ設為q維實向量,同時為了得到每個未知參數(shù)θi的概率特性,推導得出θi的邊緣密度函數(shù)如下:
p(θi|tf1,tf2,…,tfi)
(15)
式中:θ(-i)表示在θ中除去參數(shù)θi的剩余參數(shù)向量。
馬爾科夫蒙特卡洛模擬技術是一種解決上述復雜高維函數(shù)積分問題的有效方法。該方法通過迭代得到參數(shù)向量θ的一個仿真序列θ(j)(j=1,2,…,),如果j足夠大,θ(j)將收斂于一個由后驗分布中生成的隨機序列。因此,執(zhí)行大量迭代后,算法產(chǎn)生的仿真序列θ(m),θ(m+1),…,θ(m+k),可認為是來自后驗分布的樣本;相反,θ(1),θ(2),…,θ(m-1)被認為是算法的老煉期,不能作為代表后驗分布的樣本。圖2展示了馬爾科夫蒙特卡洛模擬技術的原理。
圖2 馬爾科夫蒙特卡洛模擬技術的原理圖Fig.2 Principle diagram of Markov Monte Carlo simulation technology
本次研究涉及的模型和方法較多,本研究整體方法流程圖以及各個模型和方法在不同階段的應用,見圖3。
圖3 本研究整體方法流程圖Fig.3 Overall method flow chart of this article
本文采用我國某天然氣管道的腐蝕失效數(shù)據(jù),應用上述方法對一長2 km且包含5個腐蝕缺陷的管段進行了可靠性預測與更新,以驗證本文方法的可靠性和可行性。根據(jù)2006年的在線檢測報告結果,得到該管道和腐蝕缺陷基本參數(shù)的概率特征,見表2。腐蝕管道結構可靠性評價中所涉及的隨機變量之間相關性的相關系數(shù)取值,見表3。該管道于2016年(上一次在線檢測后的第10年)發(fā)生一起腐蝕失效事件,標記時間tf1=10(指上一次在線檢測的時間是2006年,對應的標記時間為第0年)。
表2 管道和腐蝕缺陷基本參數(shù)的概率特征Table 2 Probability characteristics of basic parameters of the pipeline and the corrosion defects
表3 隨機變量之間相關性的相關系數(shù)取值Table 3 Values of correlation coefficients between random variables
圖4、圖5和圖6展示了基于Copula函數(shù)生成的隨機變量相關樣本,結果表明產(chǎn)生的隨機變量相關樣本確實具有特定的相關性,從而驗證了Copula函數(shù)生成相關樣本的有效性。
圖4 管徑與壁厚的相關樣本Fig.4 Correlated samples of diameters and wall thickness
圖5 腐蝕缺陷初始長度與初始深度的相關樣本Fig.5 Correlated samples of defect initial lengths and initial depths
圖6 不同腐蝕缺陷深度腐蝕速率之間的相關樣本Fig.6 Correlated samples of defects depth growth rate
通過100萬次的蒙特卡洛模擬試驗,在50 a預測周期內(nèi)分別對考慮隨機變量之間相關性和不考慮隨機變量之間相關性的腐蝕管段發(fā)生的小孔泄漏失效概率、爆裂失效概率和總失效概率進行了計算,其計算結果見圖7。
圖7 隨機變量之間的相關性對腐蝕管道失效概率的影響Fig.7 Effect of correlations between random variables on pipeline failure probability
由圖7可見,考慮隨機變量之間的相關性會引起腐蝕管段小孔泄漏失效概率和總失效概率的降低,而對腐蝕管段爆裂失效概率的影響則顯現(xiàn)出雙向性,即在大約27 a之前,會引起腐蝕管段爆裂失效概率的降低,在27 a之后又會引起腐蝕管段爆裂失效概率的增加。這種現(xiàn)象的出現(xiàn)與失效模式對應的管道極限狀態(tài)方程和實際管道失效數(shù)據(jù)有關,這是由于工程實際中,與管道結構可靠性評價有關的一些參數(shù)之間確實存在相關性,所以考慮參數(shù)相關性的腐蝕管道失效概率的預測結果更為準確。模擬參數(shù)相關性的關鍵在于相關系數(shù)的確定,但目前尚未有明確的數(shù)關系數(shù)取值方法,通常根據(jù)管道內(nèi)檢測數(shù)據(jù)進行統(tǒng)計或經(jīng)驗判定。
利用管道腐蝕失效數(shù)據(jù)對腐蝕管道結構可靠性進行更新的前提是將管道失效概率擬合為具有一定參數(shù)的概率分布。為了選擇最優(yōu)的概率分布函數(shù)來描述腐蝕管段的總失效概率,本文采用秩回歸方法計算出了描述管道失效概率的5種候選分布函數(shù)的相關系數(shù),見表4。
表4 描述管道失效概率的5種候選分布函數(shù)的相關系數(shù)計算結果Table 4 Calculation results of correlation coefficients of the five candidate distribution functions describing pipeline failure probability
由表4可以得到描述腐蝕管道失效概率的最佳分布及其分布函數(shù)。由于后文的主要研究目標為通過失效數(shù)據(jù)來更新和校正腐蝕管道結構的可靠性評價結果,而不關注管道具體的失效模式,故本文基于腐蝕管道總失效概率來說明貝葉斯推斷在本研究中的具體應用。表4顯示正態(tài)分布是描述腐蝕管段總失效概率的最佳分布,故采用表1中的線性方程計算得到了相應的分布函數(shù)的參數(shù)值為μ=17.003 3、σ=3.285 8,則腐蝕管段的總失效概率分布函數(shù)F(t|θ)的表達式為
(16)
采用貝葉斯方法,通過引入失效數(shù)據(jù)對公式(16)所描述的正態(tài)分布函數(shù)中的均值μ和標準差σ參數(shù)進行更新。由于μ和σ都具有正的增量,故采用伽馬分布作為它們的先驗分布。μ和σ在失效分布函數(shù)F(t|θ)中的先驗信息可以用來定義它們的先驗分布函數(shù),采取均值和均值的10%分別作為它們的先驗分布的均值和標準差。因此,p(μ)和p(σ)分別可表示如下:
(μ>0,α1>0,λ1>0)
(17)
(σ>0,α2>0,λ2>0)
(18)
式中:α1=170.033;λ1=10;α2=32.858;λ2=10;Γ(·)表示伽馬函數(shù)。
根據(jù)管道腐蝕失效數(shù)據(jù)(失效時間為2016年,標記時間tf1=10),結合貝葉斯定理,參數(shù)μ和σ的后驗密度函數(shù)可表示為
p(μ,σ|tf1)∝fsl(tf1|μ,σ)×
(19)
式中:fsl(tf1|μ,σ)是關于參數(shù)μ和σ的似然函數(shù);sl表示小孔泄漏事件。
采用OpenBUGS軟件對本算例進行貝葉斯更新的模型示意圖,見圖8。由于這是一個多參數(shù)模型,因此為每個參數(shù)設置了兩條馬爾科夫鏈,每條鏈迭代50萬次產(chǎn)生后驗樣本,認為前5 000次迭代為老煉期。表5列出了參數(shù)μ和σ的后驗均值和95%置信區(qū)間。
圖8 基于OpenBUGS的貝葉斯更新模型示意圖Fig.8 Diagram of OpenBUGS Bayesian updating model
表5 參數(shù)μ和σ的后驗均值和95%置信區(qū)間Table 5 Posterior mean and 95% confidence interval of parameters μ and σ
圖9 失效數(shù)據(jù)對管道失效概率的更新和校正Fig.9 Update and correction of failure probability of pipeline with the failure data
由圖9可見,貝葉斯更新前,管道失效時間的均值點估計為17.003 3 a,管道于第10年發(fā)生了一起小孔泄漏失效,由此可知管道的腐蝕情況比預期嚴重;基于失效數(shù)據(jù),采用貝葉斯方法更新當前的管道失效概率評估結果,結果表明同一時間更新后的管段失效概率比更新前更高,更新后的管道失效時間均值點估計為16.14 a,這與前述分析相一致,即管道腐蝕情況比預期嚴重??梢?,更新后的管道失效概率引入了失效數(shù)據(jù)蘊含的有效信息,能更準確地反映管道的腐蝕狀態(tài)。需要說明的是,失效數(shù)據(jù)使后驗均值低于先驗均值,但后驗均值不可能等于10 a,因為貝葉斯更新的意義是在先驗分布的基礎上,引入最新數(shù)據(jù),得到后驗分布,先驗信息仍然是主導后驗分布的主要原因。就本研究而言,管道是一個包含多個腐蝕缺陷的系統(tǒng),先驗分布是對管道系統(tǒng)整體的認知,失效數(shù)據(jù)則是關于管道上某個腐蝕缺陷的信息,因此失效數(shù)據(jù)可對管道結構可靠性評估結果起到更新和校正的作用。
為了完善目前腐蝕管道結構可靠性評價模型和方法,本文建立了考慮失效數(shù)據(jù)和隨機變量之間相關性的腐蝕管道結構可靠性評價方法,能夠更準確地反映腐蝕管道的腐蝕狀態(tài)和安全水平,并得到如下結論:
(1) 腐蝕管道參數(shù)之間存在多種相關性,若不加以考慮,會造成腐蝕管道結構可靠性評價結果的不準確。本文研究表明:多種相關性的綜合作用會引起腐蝕管道小孔泄漏失效概率和總失效概率的降低,對腐蝕管道爆裂失效概率的影響顯現(xiàn)出雙向性。
(2) 貝葉斯方法可以結合腐蝕管道運行期間的失效數(shù)據(jù),實時地對腐蝕管道結構可靠性預測結果進行更新和校正。算例中,管道失效時間的先驗均值估計為第17 a,引入失效數(shù)據(jù)后,其后驗均值估計為第16 a,反映了管道的腐蝕情況較預期更為嚴重。
本文對腐蝕管道結構可靠性評價模型和方法進行了改進和完善,但仍有許多其他問題需要考慮和研究,比如目前國內(nèi)管道失效相關數(shù)據(jù)的積累還不夠充分,因此需要管道運營者注重收集和整理管道失效相關數(shù)據(jù),并建立管道腐蝕缺陷數(shù)據(jù)庫,助力管道的安全管理。