張益炬 付 秋

1.西南石油大學(xué)石油工程學(xué)院，四川 成都 610500；

2.浙江偉星新型建材股份有限公司，浙江 臨海 317000；

3.中國(guó)石化華北分公司，河南 鄭州 450006

0 前言

油氣長(zhǎng)輸管線的安全問(wèn)題一直是國(guó)內(nèi)外石油天然氣行業(yè)關(guān)注的問(wèn)題， 而腐蝕引起的管線失效是失效模式中最重要的原因。 這種失效模式需要對(duì)在役管線進(jìn)行及時(shí)的安全可靠性評(píng)估［1］。李斌等人［2］用JC 算法手工計(jì)算評(píng)價(jià)某長(zhǎng)輸管線腐蝕剩余強(qiáng)度可靠性并提出用計(jì)算機(jī)計(jì)算的觀點(diǎn)，Wang Changlong 等人［3］用LabVIEW 軟件實(shí)現(xiàn)JC算法對(duì)面板堆石壩泄洪道可靠性計(jì)算； 章慧健［4］、 馮曉波等人［5］運(yùn)用MATLAB 軟件和隨機(jī)數(shù)發(fā)生器來(lái)實(shí)現(xiàn)蒙特卡洛算法計(jì)算工程結(jié)構(gòu)可靠度； 李遠(yuǎn)瑛等人［6］通過(guò)MATLAB 軟件來(lái)實(shí)現(xiàn)JC 算法、響應(yīng)面算法和蒙特卡洛算法的結(jié)構(gòu)可靠度計(jì)算；法國(guó)學(xué)者Fernández M S B 等人［7］用MATLAB 軟件實(shí)現(xiàn)蒙特卡洛算法對(duì)測(cè)量不確定度的評(píng)估；陳健等人［8］用MATLAB 軟件實(shí)現(xiàn)蒙特卡洛算法應(yīng)用于箱梁施工可靠度計(jì)算；Yu Deng 等人［9］利用蒙特卡洛算法抽樣不確定性來(lái)降低HRAS （Health Risk Assessment System）的模糊性，進(jìn)行農(nóng)村飲用水水質(zhì)健康風(fēng)險(xiǎn)評(píng)估；南伊利諾伊大學(xué)的Sheng Yanyan［10］運(yùn)用MATLAB 軟件實(shí)現(xiàn)蒙特卡洛算法對(duì)IRT（Item Response Theory）模型估計(jì)；另外還有學(xué)者運(yùn)用MATLAB 軟件實(shí)現(xiàn)蒙特卡洛算法對(duì)均土坡質(zhì)可靠度、混凝土結(jié)構(gòu)可靠度、煤巷錨桿支護(hù)結(jié)構(gòu)可靠性、公路邊坡穩(wěn)定性可靠度的模擬計(jì)算［11-16］；桂勁松等人［16］用MATLAB 軟件實(shí)現(xiàn)基于最優(yōu)化原理的蒙特卡洛算法；孫尚新等人［17］用MATLAB 軟件實(shí)現(xiàn)蒙特卡洛算法對(duì)大型復(fù)雜系統(tǒng)平均壽命的評(píng)定。由于MATLAB 軟件對(duì)長(zhǎng)輸管線腐蝕缺陷可靠性計(jì)算與軟件計(jì)算相結(jié)合的研究較少，本文在前人研究的基礎(chǔ)上，利用MATLAB 軟件實(shí)現(xiàn)了JC 算法和蒙特卡洛算法計(jì)算長(zhǎng)輸管線腐蝕缺陷的可靠性指標(biāo)β 及其驗(yàn)算點(diǎn)。

1 實(shí)例分析

1.1 基礎(chǔ)資料

根據(jù)某現(xiàn)場(chǎng)受腐蝕長(zhǎng)輸管線測(cè)試數(shù)據(jù)資料，用彈性極限法［4］進(jìn)行數(shù)據(jù)概率分布擬合得服從對(duì)數(shù)正態(tài)分布R、服從正態(tài)分布S，R 與S 為相互獨(dú)立變量，該受腐蝕長(zhǎng)輸管線的功能函數(shù)為Z=g R，（ ）

S =0，μR=100， μS=50，δR=0.12，δS=0.15 為其均值和變異系數(shù)，求可靠性指標(biāo)β 以及驗(yàn)算點(diǎn)［3］。

1.2 改進(jìn)蒙特卡洛算法與JC 算法

1.2.1 改進(jìn)蒙特卡洛算法可靠性計(jì)算

蒙特卡洛（Monte Carlo）算法是一種直接隨機(jī)抽樣方法，對(duì)不服從正態(tài)分布函數(shù)的變量，可以不進(jìn)行正態(tài)當(dāng)量化而直接進(jìn)行求解，如某隨機(jī)變量X 服從對(duì)數(shù)正態(tài)分布，可以直接借助MATLAB 軟件隨機(jī)數(shù)發(fā)生器產(chǎn)生一組所需服從對(duì)數(shù)正態(tài)分布隨機(jī)向量?；贛ATLAB 軟件的蒙特卡洛算法編程簡(jiǎn)潔，求解效率高，精確度高。 蒙特卡洛算法實(shí)際是一種隨機(jī)抽樣（樣本抽取數(shù)量越大結(jié)果越精確但耗時(shí)增加）， 以是否符合編程中設(shè)定條件為 “篩網(wǎng)”來(lái)“過(guò)濾”出總樣本中符合條件樣本數(shù)目以達(dá)到求解可靠度的方法。 蒙特卡洛算法借助MATLAB軟件實(shí)現(xiàn)步驟如下：

a）確定計(jì)算可靠度公式R（t）=P（R-S=0）；

b） 確定R 和S 的概率密度函數(shù)（pdf）：f（R）和f（S），及累積概率分布函數(shù)（cdf）：F（R）和F（S）；

c）生成R 與S 在0～1 之間的均勻分布隨機(jī)數(shù)，樣本數(shù)量N 可以取10 000 以上；

d） 檢驗(yàn)生成的每對(duì)R 隨機(jī)數(shù)與S 隨機(jī)數(shù)是否滿足R≥S，如果滿足則取1，不滿足則取0；

e）統(tǒng)計(jì)并累加共取1（滿足R≥S）的次數(shù)N1；

f） 可靠度R（t）即為滿足R≥S 的次數(shù)N1在總樣本數(shù)量N 中出現(xiàn)的頻率N1/N，失效概率Pf=1-（N1/N）。

上述蒙特卡洛算法只能求出失效概率或者可靠度，不能求得驗(yàn)算點(diǎn)。 由于可靠性指標(biāo)β 是原點(diǎn)坐標(biāo)距極限狀態(tài)曲面所有點(diǎn)中最短直線長(zhǎng)度，驗(yàn)算點(diǎn)是極限狀態(tài)曲面上此時(shí)所對(duì)應(yīng)的點(diǎn)，因此這是個(gè)最短距離的最優(yōu)化問(wèn)題。 可以采用蒙特卡洛算法在極限狀態(tài)曲面進(jìn)行隨機(jī)抽樣，再用MATLAB 軟件中“find”命令找出最短距離和在曲面上的對(duì)應(yīng)點(diǎn)，即為可靠性指標(biāo)β 和驗(yàn)算點(diǎn)［18-19］。采用基于上述最優(yōu)化原理的改進(jìn)蒙特卡洛算法進(jìn)行編程求解，改進(jìn)蒙特卡洛算法程序如下：

clc，clear

n=input（‘input how many example you want：’）；% 設(shè)定抽樣次數(shù)

f=exp（normrnd（4.598，0.119 6，1，n））；%產(chǎn)生相應(yīng)抽樣次數(shù)對(duì)數(shù)正態(tài)分布隨機(jī)數(shù)

H=f；

fx_1=normcdf（log（f），4.598，0.119 6）；% 轉(zhuǎn)化為正態(tài)累積概率分布

y1=norminv（fx_1，0，1）；%求正態(tài)分布下分位數(shù)［7］，即隨機(jī)數(shù)概率

fx_2=normcdf（H，50，7.5）；

y2=norminv（fx_2，0，1）；

b1=sqrt（y1.^2+y2.^2）；%所有隨機(jī)數(shù)可靠性指標(biāo)求解；beta=min（b1）%可靠性指標(biāo)最小值β；c=find（b1==beta）；%確定對(duì)應(yīng)位置；f（c），H（c）%輸出驗(yàn)算點(diǎn)。改進(jìn)蒙特卡洛算法借助MATLAB 軟件計(jì)算結(jié)果見(jiàn)表1。

表1 不同抽樣次數(shù)下蒙特卡洛算法計(jì)算結(jié)果（程序計(jì)算10 次均值）

1.2.2 JC 算法可靠性計(jì)算

JC 算法是由國(guó)際結(jié)構(gòu)安全度聯(lián)合委員會(huì)（JCSS）推薦的方法，簡(jiǎn)稱JC 算法，它屬于當(dāng)量正態(tài)化法。 該方法主要對(duì)不服從正態(tài)分布的隨機(jī)變量先進(jìn)行正態(tài)當(dāng)量化后，采用改進(jìn)的一次二階矩法求驗(yàn)算點(diǎn)和結(jié)構(gòu)可靠性指標(biāo)β。 正態(tài)當(dāng)量化后的隨機(jī)變量在設(shè)計(jì)驗(yàn)算點(diǎn)處分布的密度函數(shù)尾部面積保證與正態(tài)當(dāng)量化前相等以使正態(tài)當(dāng)量化前后失效概率不變。 JC 算法手工計(jì)算可靠度時(shí)，不僅復(fù)雜繁瑣，而且容易出錯(cuò)；而利用MATLAB 軟件計(jì)算能解決此問(wèn)題，且計(jì)算效率高，滿足任意可靠性指標(biāo)的前后誤差精度，避免手工計(jì)算的單一性。

JC 算法借助MATLAB 軟件實(shí)現(xiàn)［17］步驟如下：

a）存在列向量數(shù)組（i=1，2，…，n-1，n）服從某個(gè)非正態(tài)分布類型（如對(duì)數(shù)正態(tài)分布，極值Ⅰ型分布等）統(tǒng)計(jì)參數(shù)為μXi與σXi的極限狀態(tài)方程f（X1，X2，…，XN-1，XN）=0，假定取μXi=x*i。

clear，clc；

n=input（‘input level of error accuracy you want：’）；%

設(shè)定誤差精度；

mu_R_0=100；mu_p1_0=50；

sigma_R_0=1.2；sigma_p1_0=7.5；

beta_0=0；beta=1；%根據(jù)問(wèn)題要求設(shè)定初始值；

while abs（beta-beta_0）＞1e-n%根據(jù)步驟5 設(shè)置程序誤差精度作為while 循環(huán)條件；

beta_0=0；

beta_0=beta；

sigma_lnR=0.118 6；

mu_lnR=4.693；

sigma_R=mu_R_0*sigma_lnR；

mu_R=mu_R_0*（1-log（mu_R_0）+mu_lnR）；

cos_theta_R=-sigma_R/sqrt（sigma_R^2+sigma_p1_0^2）；

cos_theta_p1=sigma_p1_0/sqrt（sigma_R^2+sigma_p1_0^2）；

beta =［（mu_R_0 -mu_p1_0）+（mu_R -mu_R_0） -（mu_p1_0-mu_p1_0）... ］/［sqrt（sigma_R^2+sigma_p1_0^2）］；

new_mu_R_0=mu_R+beta*sigma_R*cos_theta_R；

new_mu_p1_0=mu_p1_0+beta*sigma_p1_0*cos_theta_p1；

mu_R_0=mu_R；

sigma_R_0=sigma_R； %為步驟2～4 的MATLAB 軟件實(shí)現(xiàn)

end

new_mu_R_0，new_mu_p1_0%輸出驗(yàn)算點(diǎn)

JC 算法借助MATLAB 軟件計(jì)算的結(jié)果見(jiàn)表2。

表2 不同循環(huán)條件精度下JC 算法計(jì)算結(jié)果（程序計(jì)算10 次均值）

2 結(jié)果對(duì)比分析

由表1～2 可看出， 改進(jìn)蒙特卡洛算法在抽樣100 萬(wàn)次的情況下已達(dá)到滿意的計(jì)算精度，JC 算法在循環(huán)條件精度為10-5時(shí)也能滿足計(jì)算要求。 對(duì)應(yīng)條件下的改進(jìn)蒙特卡洛算法和JC 算法計(jì)算結(jié)果與李斌等人［2］手工計(jì)算結(jié)果進(jìn)行對(duì)比，結(jié)果見(jiàn)表3。

表3 不同計(jì)算方法結(jié)果對(duì)比（程序運(yùn)行10 次取均值）

由表3 計(jì)算結(jié)果可知，JC 算法與改進(jìn)蒙特卡洛算法的計(jì)算結(jié)果與李斌等人手工計(jì)算相對(duì)精確解吻合良好，表明基于MATLAB 軟件計(jì)算的兩種方法能高效準(zhǔn)確地應(yīng)用于長(zhǎng)輸管線腐蝕缺陷的可靠性評(píng)估。 長(zhǎng)輸管線在受到不同簡(jiǎn)單外力（如冰載荷等）作用下，對(duì)程序做適當(dāng)?shù)恼{(diào)試就可得到所需計(jì)算模型。

改進(jìn)蒙特卡洛算法較JC 算法編程語(yǔ)言更簡(jiǎn)潔。 由表1 可知如果抽樣次數(shù)過(guò)大（某些情況下為滿足求解精度需要較大的抽樣次數(shù)）或計(jì)算模型較復(fù)雜時(shí)，對(duì)于抽樣次數(shù)呈10 倍增長(zhǎng)的改進(jìn)蒙特卡洛算法運(yùn)行時(shí)間會(huì)很長(zhǎng)，計(jì)算效率會(huì)降低；而JC 算法在滿足高精度求解情況下計(jì)算時(shí)間明顯較少， 對(duì)比結(jié)果見(jiàn)圖1。 在這種情況下JC 算法計(jì)算性能比改進(jìn)蒙特卡洛算法更穩(wěn)定。兩種算法特點(diǎn)不同，改進(jìn)蒙特卡洛算法較JC 算法程序語(yǔ)言簡(jiǎn)潔明了，MATLAB 軟件強(qiáng)大的計(jì)算功能對(duì)不同情況下長(zhǎng)輸管線腐蝕缺陷的可靠性計(jì)算能起到積極推動(dòng)作用。

圖1 改進(jìn)蒙特卡洛算法與JC 算法計(jì)算性能對(duì)比圖

在100 萬(wàn)次的抽樣精度下， 對(duì)改進(jìn)蒙特卡洛算法程序稍微修改，添加循環(huán)程序來(lái)計(jì)算在不同承載情況下該長(zhǎng)輸管線腐蝕缺陷的失效概率與可靠性，并繪制兩者關(guān)系圖，見(jiàn)圖2。

圖2 100 萬(wàn)次抽樣精度時(shí)不同受力下改進(jìn)蒙特卡洛算法可靠性計(jì)算趨勢(shì)圖

3 結(jié)語(yǔ)

a）JC 算法和改進(jìn)蒙特卡洛算法都能很好滿足工程計(jì)算精度，與李斌等人手工計(jì)算相對(duì)精確解吻合良好。兩種算法都能很好地運(yùn)用到工程計(jì)算中。

b）改進(jìn)蒙特卡洛算法較JC 算法編程語(yǔ)言更簡(jiǎn)潔，JC 算法較蒙特卡洛算法更通俗易懂。

c）在相同計(jì)算精度要求下改進(jìn)蒙特卡洛算法程序運(yùn)行時(shí)間更長(zhǎng)，JC 算法運(yùn)行更穩(wěn)定。

［1］ 方華燦. 油氣長(zhǎng)輸管線的安全可靠性分析［M］. 北京：石油工業(yè)出版社，2002.5.Fang Huacan. Reliability Analysis of Safety of Oil Gas Longdistance Transport Pipes［M］. Beijing：Petroleum Industry Press，2002.5.

［2］ 李 斌，陳 曦，王 帥. 基于JC 算法的腐蝕管道剩余強(qiáng)度可靠性評(píng)價(jià)［J］. 石油化工腐蝕與防護(hù)，2009，26（1）：62-64.Li Bin，Chen Xi， Wang Shuai.Residual Strength Reliability Evaluation of Corroded Pipeline Based on JC Method ［J］.Corrosion & Protection in Petrochemical Industry， 2009，26（1）：62-64.

［3］ Wang Changlong. Structural Reliability Calculation Based on JC Method ［C］. IITA International Conference On Services Science，Management and Engineering，2009.

［4］ 章慧健. 工程結(jié)構(gòu)可靠度計(jì)算的Matlab 實(shí)現(xiàn)［J］. 四川建筑，2007，27（1）：154-156.Zhang Huijian. Engineering Reliability Calculation Using Matlab ［J］. Sichuan Architecture，2007，27 （1）：154-156.

［5］ 馮曉波， 楊 樺. 用Matlab 實(shí)現(xiàn)蒙特卡羅法計(jì)算結(jié)構(gòu)可靠度［J］. 中國(guó)農(nóng)村水利水電，2002，（8）：50-51.Feng Xiaobo，Yang Hua.Monte Carlo Method with Matlab for Calculating Structural Reliability ［J］. China Rural Water and Hydropower，2002， （8）：50-51.

［6］ 李遠(yuǎn)瑛， 張德生. 基于Matlab 的結(jié)構(gòu)可靠度分析方法研究［J］. 嘉應(yīng)學(xué)院學(xué)報(bào)（自然科學(xué)版），2011，29（2）：45-48.Li Yuanying， Zhang Desheng. The Analysis of Plane Frame Structure Displacement Reliability Based on MATLAB ［J］.Journal of Jiaying University （Natural Science Edition）， 2011，29 （2）：45-48.

［7］Fernández M S B，Calderon J M A，Díez P M B.Imple-mentation in Matlab of the Adaptive Monte Carlo Method for the Evaluation of Measurement Uncertainties ［J］. Accred. Qual.Assur.，2009，（14）：95-106

［8］ 陳 健，盤欽卿，王衛(wèi)星.Matlab 實(shí)現(xiàn)蒙特卡洛法在箱梁施工可靠度中的應(yīng)用［J］. 重慶工商大學(xué)學(xué)報(bào)（自然科學(xué)版），2011，28（1）：90-93.Chen Jian，Pan Qinqing，Wang Weixing.Application of Monte Carlo Method Realized by Matlab to the Reliability of Boxgirder Construction ［J］.Journal of Chongqing Technology and Business University （Natural Science Edition）， 2011，28 （1）：90-93.

［9］ Yu Deng， Ni Fuquan， Yao Zhenguang. The Monte Carlo-Based Uncertainty Health Risk Assessment Associated With Rural Drinking Water Quality ［J］. Journal of Water Resource and Protection，2012，4（9），772-778.

［10］ Sheng Yanyan. Markov Chain Monte Carlo Estimation of Normal Ogive IRT Models In Matlab［J］.Journal of Statistical Software，2008，25（8）：1-15.

［11］ 楊 明，張可能，劉宇飛，等. 基于Matlab 的均質(zhì)土坡穩(wěn)定可靠度的蒙特卡羅模擬［J］. 地質(zhì)與勘探，2009，45（1）：89-94.Yang Ming， Zhang Keneng， Liu Yufei， et al. Monte Carlo Simulation of Homogeneous Soil-Slope Reliability Based on Matlab［J］.Geology and Exploration，2009，45 （1）：89-94.

［12］ 王偉宏. 基于Matlab 實(shí)現(xiàn)混凝土結(jié)構(gòu)可靠度分析［J］. 路基工程，2008，9（4）：161-162.Wang Weihong. Realizing by Matlab Concrete Structure Reliability Analysis ［J］. Subgrade Engineering， 2008，9 （4） ：161-162.

［13］ 朱川曲，張道兵，朱海燕. 基于蒙特卡羅法的煤巷錨桿支護(hù)結(jié) 構(gòu) 可 靠 性 分 析［J］. 中 國(guó) 安 全 科 學(xué) 學(xué) 報(bào)，2008，18（4）：146-150.Zhu Chuanqu， Zhang Daobing， Zhu Haiyan. Reliability Analysis on Bolt Support Structure of Coal Roadway Based on Monte-Carlo ［J］. China Safety Science Journal， 2008，18（4）：146-150.

［14］ 梁利端，李玲玲. 淺談構(gòu)件可靠度基于Matlab 的蒙特卡羅法［J］. 蘭州工業(yè)高等?？茖W(xué)校學(xué)報(bào)，2008，15（3）：7-9.Liang Liduan， Li Lingling. Component Reliability Based on the Monte Carlo of MATLAB ［J］. Journal of Lanzhou Polytechnic College，2008，15 （3）：7-9.

［15］ 周泉宇，楊仕教，羅 輝. 談公路邊坡穩(wěn)定性可靠度的蒙特卡羅計(jì)算方法［J］. 山西建筑，2009，35（3）：14-15.Zhou Quanyu，Yang Shijiao，Luo Hui.Discussion on Highway Slope Stability Reliability Based on Monte Carlo Calculation Method ［J］.Shanxi Architecture，2009，35 （3）：14-15.

［16］ 桂勁松， 康海貴. 結(jié)構(gòu)可靠度計(jì)算的最優(yōu)化方法及其matlab 實(shí)現(xiàn)［J］. 四川建筑科學(xué)研究，2004，30（2）：18-20.Gui Jingsong， Kang Haigui. The Structural Reliability Calculation Base on Matlab［J］. The Building Science Research of Sichuan，2004，30 （2）：18-20.

［17］ 孫尚新，高精先，張 林. 用Matlab 實(shí)現(xiàn)蒙特卡羅法評(píng)定大型復(fù)雜系統(tǒng)平均壽命［J］. 中國(guó)安全生產(chǎn)科學(xué)技術(shù)，2006，2（6）：24-26.Sun Shangxin， Gao Jingxian，Zhang Lin. The MTTF'S Evaluation of Large-scale Complex System with MATLAB and Monte Carol［J］.Journal of Safety Science and Technology，2006，2 （6）：24-26.

［18］ 周 品，趙新芬.Matlab 數(shù)學(xué)建模與仿真［M］. 北京：國(guó)防工業(yè)出版社，2009.4.Zhou Pin， Zhao Xinfen. Matlab Mathematical Modeling and Simulation ［M］. Beijing： National Defense Industry Press，2009.4.

［19］ 任虹宇， 楊 洋， 張 旭. 輸氣站站內(nèi)埋地管道防腐層PCM 檢測(cè)［J］. 天然氣與石油，2012，30（6）：58-59.Ren Hongyu， Yang Yang， Zhang Xu. PCM Inspection on Corrosion Control Coating of Buried Pipeline in Gas Transportation Station［J］.Natural Gas and Oil，2012，30 （6）：58-59.