基于非接觸磁特性參數(shù)的埋地管道可靠度MCMC定量評價模型

邢海燕，王松弘澤，孫曉軍，段成凱，劉 傳，弋 鳴，劉偉男

(1.東北石油大學 機械科學與工程學院，黑龍江大慶 163318；2.哈爾濱工程大學 動力與能源工程學院，哈爾濱 150001)

0 引言

隨著我國油氣資源的開發(fā)建設，根據(jù)“十三五”規(guī)劃，2020年底我國油氣長輸管道總里程已達到16萬km，但目前其中大約一半已服役超過20年，老化問題十分突出，且埋地油氣長輸管道運行環(huán)境復雜，一旦發(fā)生泄漏、爆炸等事故，往往造成災難性的后果[1-2]。因此，可靠度評價對于埋地管道安全運行至關重要。但是，對于在役運行的埋地管道，進行全面開挖檢測不經(jīng)濟、也不現(xiàn)實，而常用的無損檢測方法在不開挖條件下對埋地管道進行可靠度評價存在困難，渦流、漏磁、超聲等方法在非接觸檢測領域，分別存在趨膚效應、周向磁化困難、阻抗不匹配等檢測障礙[3]；漏磁、超聲等管道內(nèi)置式檢測，由于空間和角度的限制經(jīng)常造成卡堵問題、檢測成本非常高。

在20世紀末興起的金屬磁記憶檢測技術(MMM)，既可發(fā)現(xiàn)宏觀缺陷；又能發(fā)現(xiàn)因異常應力集中引起的早期損傷，在無損檢測領域顯示其獨特優(yōu)勢。XU等[4]研究了不同疲勞階段磁信號間的變化規(guī)律，得出磁性信號隨位錯釘扎因子的增加而增強的結論；KOLOKOLNIKOV等[5]研究了弱磁下的焊接過程中，應力變化與自磁場強度之間關系；劉斌等[6-7]根據(jù)電子自旋理論和洪德法則，基于全勢線性綴加平面波法(FLAPW)建立了磁力學模型，得到屈服極限前后的磁記憶信號特征，并定量解釋應力集中與自漏磁場之間的關系；HU等[8]研究了對接焊縫宏觀缺陷影響因素，獲得不同靈敏度方向和剝離值的磁記憶信號的特性；胥永剛等[9]提出了基于固有時間尺度分解(ITD)的磁記憶特征提取方法，研究低速重載齒輪故障的早期檢測；蘭清生等[10]通過研究李薩如圖當量面積變化，探討試件不同的應力集中狀態(tài)；易方等[11]提出一種基于模糊核支持向量機的管道缺陷識別方法；DUBOV[12]針對冶金生產(chǎn)中缺陷的不確定性難題，提出一種基于MMM檢測的缺陷診斷方法；楊理踐等[13]基于密度泛函理論，建立鐵磁晶體磁記憶模型,采用正交化平面波贗勢法研究塑性變形后磁記憶信號變化特征；邢海燕等[14]利用MMM方法，建立了不同載荷下角焊縫缺陷處應力集中程度的定量表征模型。

針對不開挖條件下、在役運行的埋地管道可靠度定量評價難題，提出基于金屬磁記憶檢測技術的可靠度定量評價方法。為此，進行埋地管道非接觸磁記憶檢測試驗，首先引入Kendall相關系數(shù)，分析缺陷尺寸與不同磁特性參數(shù)之間的相關性，并獲取不同壓力條件下ΔH(y)與d、k(x,y)與L之間的非線性回歸擬合修正式；之后結合Modified B31G評價標準，利用埋地管道的多種非接觸磁特性參數(shù)，在蒙特卡洛法基礎上，引入能夠反映管道缺陷動態(tài)時變性的馬爾科夫鏈，建立基于非接觸磁特性參數(shù)的埋地管道可靠度馬爾科夫鏈蒙特卡洛(Markov Chain Monte Carlo,MCMC)定量評價模型。通過算法比較與現(xiàn)場試驗，驗證模型的有效性，為實際工程中的埋地管道可靠度定量評價提供一種新的方法。

1 試驗方法

試驗材料采用輸油管道常用的L245N鋼管，尺寸為?325 mm×15 000 mm，壁厚10 mm，一端焊有盲板，另一端留有接口與液壓泵連接，對管道試件施加不同等級內(nèi)壓：0,3,5,8,10,12,15 MPa(其中0,5,10,15 MPa為分析數(shù)據(jù)，其余為檢測數(shù)據(jù))。為研究基本條件下，不同尺寸缺陷的非接觸磁特性參數(shù)表征規(guī)律，如圖1所示，在第1組管件沿管長方向預制4個圓形缺陷，以模擬實際工況中常見的點蝕腐蝕，長度L均為15 mm，深度d分別為2,4,6,8 mm；在第2組管件以同樣間隔預制4個缺陷，深度d均為4 mm，長度L分別為5,10,15,20 mm，每組管件數(shù)量均為6個。試驗中將2組管件埋地1 m深，在地面上做好軸線標記，利用TSC-5M-32型磁記憶儀和11-6W型非接觸探頭組成檢測系統(tǒng)，對不同內(nèi)壓下的管件沿軸線的3條測試路徑進行非接觸掃描檢測，獲取不同尺寸缺陷的非接觸磁記憶信號，試驗布置如圖2所示。

(a)第1組管件

(b)第2組管件圖1 試驗用管件示意Fig.1 Schematic diagram of tubular specimen

(a)測試路徑布置

(b)試驗布置圖2 試驗布置示意Fig.2 Schematic diagram of experimental layout

2 試驗結果分析

2.1 磁特性參數(shù)提取

通過試驗得到不同內(nèi)壓P下，所對應不同缺陷長度L和缺陷深度d的非接觸磁記憶數(shù)據(jù)，提取如下磁特性參數(shù)(如圖3所示)。

切向磁場強度峰峰值ΔH(x):

ΔH(x)=H(x)max-H(x)min

(1)

式中,H(x)max,H(x)min分別為切向磁記憶信號在缺陷附近區(qū)域的最大值與最小值。

法向磁場強度峰峰值ΔH(y):

ΔH(y)=H(y)max-H(y)min

(2)

式中,H(y)max,H(y)min分別為法向磁記憶信號在缺陷附近區(qū)域的最大值與最小值。

合成磁場強度峰峰值ΔH(x,y):

(3)

切向磁場強度變化率k(x):

(4)

式中,x為坐標位置；Δx1為切向磁記憶信號在缺陷附近區(qū)域取最大值與最小值時的位置坐標之差。

法向磁場強度變化率k(y):

k(y)=[H(y)max-H(y)min]/Δx2

(5)

式中,Δx2為法向磁記憶信號在缺陷附近區(qū)域取最大值與最小值時的位置坐標之差。

合成磁場強度變化率k(x,y):

(6)

(a)原始信號

(b)特性參數(shù)提取圖3 磁記憶原始信號與特性參數(shù)提取Fig.3 Extraction of magnetic original signals andcharacteristic parameters

2.2 缺陷深度與磁特性參數(shù)的關系

經(jīng)試驗數(shù)據(jù)分析可知，第2.1節(jié)中各磁特性參數(shù)與缺陷深度d、缺陷長度L對應關系的相關程度存在差異。為此，在數(shù)據(jù)總體分布未知情況下，引入Kendall相關系數(shù)[15]，對不同缺陷深度和長度磁特性參數(shù)的相關性進行研究。

首先，由試驗獲取0 MPa內(nèi)壓、缺陷長度L=15 mm、缺陷深度分別為d=2,4,6,8 mm對應的ΔH(x),ΔH(y),ΔH(x,y),k(x),k(y)，k(x,y)數(shù)值，通過統(tǒng)計產(chǎn)品與服務解決方案軟件(Statistical product and service solutions，SPSS)計算得到不同缺陷深度d的Kendall相關系數(shù)值；其次，按上述方法繼續(xù)施加不同等級內(nèi)壓P=5,10,15 MPa，完成不同內(nèi)壓下各磁特性參數(shù)Kendall系數(shù)值的計算。以壓力P為橫坐標，Kendall值為縱坐標，得到不同內(nèi)壓下各磁特性參數(shù)的Kendall相關系數(shù)，如圖4所示。由圖4可知，缺陷長度L一定，隨壓力P和缺陷深度d的改變，相對于其他磁特性參數(shù)，唯有法向磁場強度峰峰值ΔH(y)在各壓力等級下與d相關系數(shù)均接近1，證明ΔH(y)可敏感反映并表征d的變化規(guī)律。在此基礎上，給出ΔH(y)與d之間的對應關系，如圖5所示。

圖4 不同內(nèi)壓下各磁特性參數(shù)Kendall相關系數(shù)Fig.4 Kendall correlation coefficient of each magneticcharacteristic parameter at different internal pressures

圖5 ΔH(y)與d的試驗數(shù)據(jù)關系曲線Fig.5 Curves of experimental data for ΔH(y)-d correlation

對圖5的ΔH(y)與d關系曲線進行非線性回歸擬合，可得：

+40.89

(7)

考慮ΔH(y)公式值與實際值間存在誤差，引入系數(shù)λ,μ分別對P和d進行修正，如下式：

ΔH(y,λ,μ)=2.95exp(0.13λP)-[2.95

×exp(0.13λP)+23.91]/[1

(8)

為保證獲取最優(yōu)的修正系數(shù)λ,μ，將圖5中數(shù)據(jù)點所對應的d,P，ΔH(y)值與式(8)相結合，獲取一系列修正系數(shù)λ,μ。根據(jù)每組λ,μ所對應式(8)中ΔH(y)計算值與試驗值間的相對誤差，繪制不同修正系數(shù)λ,μ與相對誤差的三維圖，見圖6?？梢钥闯?，相對誤差先減小、后增大，即存在相對誤差極小值，優(yōu)選最小相對誤差對應的修正系數(shù)：λ=0.97,μ=0.95。

圖6 不同修正系數(shù)對應的相對誤差三維圖Fig.6 3D diagram of relative errors for differentcorrection coefficients

將試驗中3,8,12 MPa壓力下的ΔH(y)值用于修正公式的驗證，結果如圖7所示。分析得出，同一壓力和缺陷深度條件下，ΔH(y)實測值與公式值的誤差最大為6.07%，說明修正公式有效。

圖7 ΔH(y)實測值與預測值對比Fig.7 Comparison between measured and predictedvalues of ΔH(y)

2.3 缺陷長度與磁特性參數(shù)的關系

引入Kendall相關系數(shù)，按上述方法完成對試驗所獲各磁特性參數(shù)與缺陷長度L之間相關性的定量研究，獲取不同內(nèi)壓下各磁特性參數(shù)的Kendall相關系數(shù)，結果見圖8?？梢钥闯觯毕萆疃萪一定，隨壓力P和缺陷長度L的改變，相對于其他磁特性參數(shù)，合成磁場強度變化率k(x,y)在各壓力等級下與L相關系數(shù)均接近1，證明k(x,y)可敏感性反映并表征L的變化規(guī)律。進一步給出了k(x,y)與L之間的對應關系，如圖9所示。

圖8 缺陷長度Kendall相關系數(shù)Fig.8 Kendall correlation coefficient of defect length

圖9 k(x,y)與L的試驗數(shù)據(jù)關系曲線Fig.9 Curves of experimental data for k(x,y)-L correlation

對圖9所示的k(x,y)與L關系曲線進行非線性回歸擬合，同樣考慮k(x,y)公式值與實際值間存在誤差，引入系數(shù)α,β分別對P和L進行修正，如下式：

k[(x,y),α,β]=[0.30exp(0.16αP)+2.21]

(9)

為了保證最優(yōu)修正系數(shù)α,β，采用第2.2節(jié)的相同方法，繪制不同組α,β下對應的相對誤差三維圖，如圖10所示，可以看出，相對誤差先減小、后增大，存在誤差的極小值，優(yōu)選修正系數(shù)：α= 1.04,β= 0.98。

將試驗中3,8,12 MPa壓力下的k(x,y)值用于修正公式的驗證，結果如圖11所示。分析得出，同一壓力和缺陷長度條件下k(x,y)的實測值與公式值的誤差最大為8.79%，修正公式有效。

圖10 不同修正系數(shù)對應的相對誤差三維圖Fig.10 3D diagram of relative errors for differentcorrection coefficients

圖11 k(x,y)實測值與預測值對比Fig.11 Comparison between measured andpredicted values of k(x,y)

3 基于非接觸磁特性參數(shù)的埋地管道可靠度MCMC模型建立

3.1 磁特性參數(shù)分布

由上述得知，缺陷深度d和缺陷長度L為影響管道失效的關鍵性因素。為建立基于磁特性參數(shù)的可靠度模型，將兩組管件240組ΔH(y)實測數(shù)據(jù)進行統(tǒng)計分析，得到統(tǒng)計結果如圖12所示。

圖12 ΔH(y)統(tǒng)計分布圖Fig.12 Statistical distribution of ΔH(y)

同理，將兩組管件240組k(x,y)實測數(shù)據(jù)進行統(tǒng)計分析，得到統(tǒng)計結果如圖13所示?？梢钥闯靓(y)符合正態(tài)分布、k(x,y)符合對數(shù)正態(tài)分布，可以基于上述磁特性參數(shù)建立埋地管道的可靠性模型。

圖13 k(x,y)統(tǒng)計分布圖Fig.13 Statistical distribution of k(x,y)

3.2 可靠度定量模型建立

對于油氣長輸管線來說，其破損通常與管道的失效應力r、實際應力s直接相關，將r和s看成隨機變量，則描述管道失效的極限狀態(tài)函數(shù)Z為：

Z=r-s

(10)

若Z>0，則管道失效應力r大于實際應力s，管道可靠；反之，則管道將損壞。

在油氣管道評價中，國際上廣泛采用的Modified B31G準則有著高可靠度與低誤差的特性[16]。為使所建可靠性模型具有實際應用價值，引入該準則：

(11)

式中，r為環(huán)向失效應力,MPa;σf為流變應力,MPa,σf=1.1σs(σs為屈服應力);t為管道公稱壁厚,mm;D為管道外徑,mm；d為缺陷深度,mm;M為Folias系數(shù);L為缺陷軸向長度,mm。

當L2/Dt≤50時：

(12)

當L2/Dt>50時：

(13)

綜合式(11)～(13)可知，缺陷深度d和缺陷長度L是影響管道失效的重要因素，且對于在役運行埋地管道而言，該缺陷變量難以獲取。結合非接觸磁記憶檢測優(yōu)勢與第2.2,2.3節(jié)內(nèi)容，建立埋地管道可靠性評價模型。

為此，分別求解ΔH(y)關于d、k(x,y)關于L的偏導：

(14)

(15)

從式(14)(15)可知，偏導存在且大于零，即反函數(shù)存在：

d(ΔH(y),P)=(1.06-0.0408P)ln{[16.89-ΔH(y)]

÷[ΔH(y)-2.95exp(0.126P)

-40.89]}+4.00

(16)

L(k(x,y),P)=(7.88-0.1485P)ln{[k(x,y)

-30.93]/[0.30exp(0.1664P)

+2.21]}+2.46

(17)

用式(16)(17)分別替換式(11)～(13)中d和L，聯(lián)立式(10)(11)，即可得到基于非接觸磁特性參數(shù)的可靠性模型：

(18)

(19)

其中：L2(k(x,y),P)/Dt≤50

(20)

其中：L2(k(x,y),P)/Dt>50

3.3 基于MCMC的可靠度計算

建立可靠度模型后，在計算可靠度時，一般利用蒙特卡洛方法。但由于常規(guī)的蒙特卡洛方法在計算過程中，隨機產(chǎn)生的樣本之間存在獨立性，不考慮當前狀態(tài)對未來狀態(tài)的影響，往往需要抽取大量樣本，收斂速度較慢。

由于馬爾科夫鏈能夠考慮到樣本之間的關聯(lián)性與時變性，且收斂速度較快，故在傳統(tǒng)蒙特卡洛法的基礎上引入馬爾科夫鏈，求解埋地管道的可靠度，利用當前管道運行狀態(tài)對后續(xù)時刻的可靠度進行評價，解決動態(tài)性問題，克服工程實際中歷史樣本難以獲取的問題。計算可靠度模型式(18)概率值，即可靠度R(t)：

R(t)=P(Z=r-s>0)

=Pr(r)Ps(s)

(21)

式(21)中Pr(r),Ps(s)分別為失效應力r及實際應力s的累積概率，其計算公式分別為：

(22)

(23)

基于Metropolis-Hasting (M-H)采樣的可靠度模型參數(shù)后驗分布計算流程(如圖14所示)如下：(1)設定狀態(tài)轉移次數(shù)閾值m，需要的樣本個數(shù)N；(2)從前述先驗概率分布采樣得到初始狀態(tài)值x0，設定t=0；(3)從條件概率分布Q(x∣xt)中采樣得到樣本值x*，從均勻分布中采樣u～[0,1]。如果u<α(xt,x*)，則接受轉移(即xt+1=x*)，否則不接受轉移；(4)重復步驟(3)直至達到狀態(tài)轉移次數(shù)閾值，生成的樣本集即為待求參數(shù)Pr(r),Ps(s)的后驗分布π(x)。

于是有N個r及N個s值，構成了r和s的隨機抽樣數(shù)組。若ri>si，則可靠；最后，確定ri>si的總次數(shù)Ni，則有：

R(t)=Ni/N

(24)

圖14 Metropolis-Hasting采樣流程圖Fig.14 Metropolis-Hasting sampling flow chart

3.4 可靠度計算結果及分析

利用試驗管道基本參數(shù)與試驗數(shù)據(jù)，對不同壓力下的R(t)進行計算，獲得埋深1 m，不同ΔH(y)與k(x,y)的可靠度R(t)云圖。

圖15 在不同ΔH(y)與k(x,y)條件下的R(t)云圖

由圖15可看出，在一定壓力下，埋地管道可靠度R(t)隨著磁特性參數(shù)數(shù)值的增加而減??；磁特性參數(shù)一定時，R(t)隨壓力的增加而顯著減小，結果符合工程實際中R(t)的變化規(guī)律。

3.5 算法比較

相同條件下，采用傳統(tǒng)MC算法計算的可靠度與MCMC方法計算的結果進行對比校驗。圖16示出相同條件下兩種方法計算的可靠度收斂曲線,MCMC方法在收斂速度上明顯優(yōu)于傳統(tǒng)MC方法，后者在抽樣進行到26 000次才顯現(xiàn)收斂趨勢；而前者在抽樣進行到5 500次左右時即成功收斂，兩種方法抽樣次數(shù)相差近5倍。

同時，MCMC方法在收斂穩(wěn)定性上也顯著優(yōu)于傳統(tǒng)MC方法，分別對MC與MCMC方法進行10次重復計算，由圖17對比看出，MCMC方法計算結果的離散程度明顯小于傳統(tǒng)MC方法，穩(wěn)定性更好。對10次重復試驗獲得的可靠度R(t)進行量化計算，MCMC方法的方差僅為0.041 7，而MC方法的方差為0.102 6。

圖16 MC與MCMC之間的可靠度R(t)收斂曲線對比Fig.16 Comparison of R(t) convergence curves betweenMC and MCMC

(a)MC方法

(b)MCMC方法圖17 MC與MCMC收斂穩(wěn)定性對比Fig.17 Convergence stability comparison betweenMC and MCMC

MCMC方法在收斂速度和穩(wěn)定性上的優(yōu)勢，使埋地管道可靠度計算所需時間大大縮短、精度大大提高，提升了該方法的工程實用性與有效性。

4 可靠度模型驗證

為了對所建可靠度模型進行現(xiàn)場驗證，對大慶采油五廠杏南中轉站某埋地輸油管道進行現(xiàn)場試驗。該埋地管道平均埋深為1 m，材質(zhì)為L245N管線鋼，其相關參數(shù)如表1所示。

表1 管道參數(shù)Tab.1 Parameters of the pipeline

通過現(xiàn)場試驗測得的數(shù)據(jù)，在利用工程圖紙排除管線中彎頭、法蘭等干擾因素對信號的影響后，利用第3.3節(jié)中構建的可靠度模型，對該埋地管道進行可靠度預測。其中可靠度異常區(qū)段如表2第2列所示，進一步對預測的可靠度異常區(qū)段進行開挖驗證，同時進行超聲波測厚復檢，結果如表2第4列所示，其中區(qū)段610～620 m出現(xiàn)泄漏，開挖后現(xiàn)場照片見圖18，利用管道實際缺陷參數(shù)，采用Modified B31G標準對異常區(qū)段進行可靠度驗證計算，其結果見表2第3列所示。

表2 埋地管道可靠度異常區(qū)段Tab.2 Reliability abnormal section of buried pipeline

由表2第2，3列可以看出，基于磁特性參數(shù)、利用所建的可靠度模型所得的可靠度，與基于實際缺陷尺寸、利用Modified B31G標準計算的驗證可靠度之間的平均誤差約為8.24%，因此本文所提出的埋地管道可靠度定量評價模型具有有效性，為在不開挖情況下，計算在役埋地管道可靠度，提供了一種新的方法。

圖18 泄漏區(qū)段管道開挖結果Fig.18 Excavation verification photo of the pipelinein the leaking section

5 結論

(1)引入Kendall相關系數(shù)，對試驗所獲非接觸磁特性參數(shù)的相關性進行定量研究，發(fā)現(xiàn)ΔH(y)和k(x,y)可分別敏感地反映并表征d,L變化規(guī)律。

(2)分別非線性擬合ΔH(y)與d、k(x,y)與L之間的回歸函數(shù)，并優(yōu)選修正系數(shù)，經(jīng)檢測樣本驗證得知，公式修正后最大誤差分別為6.07%和8.79%。

(3)基于非接觸磁特性參數(shù)，結合Modified B31G評價標準，在蒙特卡洛法基礎上，引入馬爾科夫鏈，建立埋地管道可靠度定量評價模型。與MC方法相比，MCMC方法計算可靠度所需抽樣次數(shù)大幅減少，穩(wěn)定性更好。

(4)現(xiàn)場試驗結果表明，開挖后利用管道實際缺陷尺寸，利用Modified B31G評價標準直接計算的驗證可靠度，與本文所建的可靠度定量評價模型計算的可靠度之間平均誤差約為8.24%，這為實際工程埋地管道非接觸檢測與可靠度定量評價提供了一種新的思路。