基于動態故障樹的煤礦突水風險概率評價

王華鎣，李仲學，趙怡晴，雷　煒，舒　楊

(1.北京科技大學土木與環境工程學院，北京 100083；2.北京科技大學金屬礦山高效開采與安全教育部重點實驗室，北京 100083)

?

基于動態故障樹的煤礦突水風險概率評價

王華鎣1，2，李仲學1，2，趙怡晴1，2，雷煒1，2，舒楊1，2

(1.北京科技大學土木與環境工程學院，北京 100083；2.北京科技大學金屬礦山高效開采與安全教育部重點實驗室，北京 100083)

摘要：故障樹分析是煤礦突水的常用風險評估方法，考慮隨著時間的推移基本事件和布爾邏輯關系變化的動態故障樹是煤礦安全事故研究的熱點領域，而其中一個難題是在動態故障樹仿真時，難以獲得基本事件在不同時間點的失效或修復概率。本文邀請專家評價我國煤礦突水動態故障樹分析中基本事件的失效或修復可能性，采用三角模糊數對專家的評價語言量化，通過反模糊化和矩量法，計算基本事件失效和修復時間的Weibull分布參數，依據計算和查閱的數據，對煤礦突水動態故障樹模型進行Monte Carlo仿真，結果表明我國煤礦在動工2.5年后突水事故發生概率在0.001上下波動。

關鍵詞：煤礦突水；動態故障樹；三角模糊數

煤礦突水主要是指圍巖內的地下水以超常流速或體積涌入到井巷空間及井下作業場所，可能對人員、設備、作業、環境等造成損害、影響或風險的現象。突水一旦發生或處置不當，就會導致較為嚴重的事故后果。據2013年資料，我國煤礦有礦井12985個，其中極復雜和復雜水文地質類型的礦井有905個[1]。隨著煤礦安全管理的加強，我國煤礦突水事故次數從2002年的162起下降到2013年的21起，死亡人數從2005年的605人下降到2013年的89人，但水害一直是煤礦地下開采的重要危險源之一，在各種事故類型中，突水事故死亡人數和發生起數仍列第3位(瓦斯和頂板事故居前兩位)。

故障樹分析(Fault Tree Analysis，FTA)是煤礦突水風險評估的常見方法之一[2-4]，由于煤礦突水基本事件(Basic Event，BE)難以確切界定或相關失效數據難以獲取，突水FTA研究往往停留在最小割集、最小徑集、結構重要度等方面，涉及煤礦突水風險概率的研究較為少見。一般的FTA以BE失效概率在風險評估對象中固定不變、且布爾邏輯關系中未考慮時間影響因素為假設條件，然而在實際的系統風險評估中，BE失效概率和邏輯門關系都會受到系統時間的影響[5-6]。而動態故障樹(Dynamic Fault Tree，DFT)[5]在原有的靜態布爾邏輯門的基礎上，考慮多個事件之間發生故障的先后順序，運用動態門予以表征，包括功能相關門(Functional Dependency，FDEP)、優先與門(Priority AND，PAND)、順序相關門(Sequence Dependency，SEQ)、備件門(SPARE)等。運用DFT進行風險概率評估時，有多種計算方法，如Markov過程模擬、貝葉斯網絡模擬[7]、Monte Carlo模擬以及數值積分等[8]。

1煤礦突水DFT模型

影響煤礦突水的因素較多，且經常相互影響，形成隱患，若處置不當，便可引發突水事故，并導致嚴重后果或危害。為了分析和表征突水因素、隱患、事故及后果的相互作用關系，評估相關風險，可以采用FTA模型。譬如，王長申等針對一般煤礦安全評價不足以分析突水危險性的狀況，基于煤礦突水事故致因理論建立了突水專項評價的兩級FTA模型[9]；再如，陳建民等為了對煤礦突水事故進行一般性評價，總結了煤礦的突水水源、導水通道等，建立了具有一般性的四級FTA模型[10]。至目前為止，煤礦突水FTA模型分析多為靜態的，未考慮時間因素對模型的影響。

但是，在現實系統中，系統狀態不僅與相關事件的組合有關，而且會與事件發生的順序(譬如時間順序)相關。以煤礦井下突水為例，礦井涌水量急劇增加事件、排水能力不足事件與應急措施失效事件結合起來，就會導致突水，可以認為是靜態故障樹概念；而以上三種事件不僅存在結合而且存在事件發生順序的關系，即礦井涌水量急劇增加事件先于排水能力不足事件發生，且排水能力不足事件先于應急措施失效事件發生，才會導致突水，類似情況就需要運用動態故障樹建模予以表征。本文在現有煤礦突水風險評價的靜態FTA模型基礎上，通過增設FTA動態門，把傳統的FTA擴展到DFT，建立了煤礦突水風險概率評價的DFT模型，既能夠體現突水故障樹的靜態特征，也可以反映動態特征，如圖1所示。

圖1　煤礦突水DFT模型

本文DFT模型構建的基本思路是：依據我國煤礦突水的相關文獻、規范、事故調查等基礎信息及數據，建立一般FTA模型，然后對FTA模型進行動態布爾邏輯關系分析，形成DFT模型。關于FTA模型的建立，主要參照《煤礦防治水手冊》、《煤礦安全規程》、《煤礦防治水規定》等，結合事故調查、文獻查閱、專家咨詢等獲取的信息，以突水水源充分積聚作為突水事故發生的先決條件，構建了涉及地表水、孔隙水、裂隙水、巖溶水和老空水等水源的四級FTA模型。關于突水事件動態布爾邏輯關系分析，主要是引入了順序相關門SEQ和優先與門PAND，其中SEQ門用來表征僅當排水能力不足事件先于應急措施失效事件發生，且礦井涌水量急劇增加事件先于排水能力不足事件發生，煤礦突水事故才會發生的情況，擴展了一級靜態與門為動態門；PAND用來表征僅當導水通道活化事件先于充水誘因事件發生時，突水事件才會發生的情況，擴展了兩個二級靜態與門中的礦井涌水量急劇增加為動態門。

如圖1所示，通過DFT分析，識別得到煤礦突水的人為、技術、環境和管理等方面的22種相關隱患因素(或者基本事件，即圖中的BE)，其中人為隱患因素9種，分別以BE6、BE7、BE8、BE12、BE13、BE15、BE17、BE19、BE21表示；技術隱患因素2種，分別以BE18、BE20表示；環境隱患因素7種，分別以BE1、BE2、BE3、BE4、BE5、BE9、BE10表示；以及管理隱患因素4種，分別以BE11、BE14、BE19、BE22表示。

2基于三角模糊數的突水BE概率評價

鑒于煤礦突水事故及其種類的偶發性，通常難于獲得同類BE的歷史數據或失效/修復概率分布，因此，本文采用三角模糊數表征突水BE的概率。

2.1突水BE概率的三角模糊數表征

由于煤礦突水BE的發生具有較大偶然性，其概率分布通常難以確切獲知，即使是領域專家也難于對突水BE的失效概率或修復概率給出客觀的數值表征，因此，主觀的語言描述往往是處理類似問題的唯一途徑。就此而言，模糊集的理論及方法有著獨到的優勢，本文即使采用三角模糊數方法將專家對突水BE發生可能性或者概率的主觀評價語言轉化為定量表征。

(1)

將評語集V的評價等級元v劃分為很低、低、中等、高和很高等五級，Purda等在采用三角模糊數度量核能工廠事故發生可能性時，將專家的評價語言分為7級，其中“中等”三角模糊數為(0.35，0.5，0.65)[11]，核能工廠事故與煤礦突水事故具有相似的特點——事故發生概率小，后果嚴重，本文評語為5級，將“中等”三角模糊數隸屬區間定為(0.3，0.5，0.7)，類似地，得到如下的煤礦突水BE可能性評價等級的三角模糊數：突水DFT中BE失效可能性和修復可能性為很高和很低幾率較小，高、低次之，大多趨于中等左右，經過多次試算，最終確定的三角模糊數，如式(2)～(6)所示。

(2)

(3)

(4)

(5)

(6)

假設領結分析模型中有m個事件需要n個專家對其進行評價，可得到m×n階矩陣，見式(7)。

(7)

式中：μmen為第n位專家對第m個BE的可能性評價語言的三角模糊數，如這位專家的評語為“低”，μmen=μ低=μ2=(0.08，0.23，0.38)。

2.2突水BE概率計算

(8)

式中：rmen為第n位專家對第m個BE的模糊評價值。

(9)

(10)

式(9)中，f為BE概率實際值；r為模糊評價值。

利用式(9)和模糊評價值計算BE概率的實際值，得到矩陣P(式(11))。

(11)

式中，pmen為第n位專家評價的第m個BE概率實際值。

2.3BE概率的時間分布

事件在時刻t的失效率fp是指已經運營到t時刻，在t時刻后的單位時間內發生事件所描述的事故的概率fp[14]，換而言之在1單位時間內發生fp次事件所描述的事故，通常fp<=1，即在1/f單位時間內發生1次BE所描述的事故。同理，BE所描述的事故發生后，在1/rp(rp為修復率)單位時間內成功修復1次事故。根據以上推論，計算系統BE失效或修復時間，如式(12)所示。

(12)

式中，tmen為第n位專家評價的第m個BE風險失效或修復時間。

近年來二參數Weibull分布廣泛應用在設備壽命等研究中[15]，Maiti等在對地下煤礦風險進行定量評估時，發現井下機械事故、運輸事故等失效間隔時間(time between failures，TBF)通常服從weibull分布[16，17]，假設事故第一次失效時間為BE所描述事故的狀態改變起點，TBF就是BE的失效時間，因此本文BE的失效時間按照Weibull分布擬合，為使修復時間分布特點一致，BE修復時間亦采用Weibull分布擬合，獲得BE在不同時間點失效或修復的概率。

二參數Weibull分布概率密度公式見式(13)，分布函數見式

(13)

(14)

式中：η為Weibull分布尺度參數；β為Weibull分布形狀參數。

二參數Weibull分布參數估算方法較多，有圖解法、極大似然法、矩量法(Moment method)、最小二乘法、經驗法等[18-19]，Paulo等在對Brazil地區風速利用Weibull分布擬合進行預測時，將7種Weibull分布參數估計值與當地歷史風速統計數據進行比較，就擬合后的平均值和標準差，矩量法估算結果較其他方法更接近歷史真實數據[20]。本文采用矩量法對BE失效或修復時間進行Weibull分布參數估計，具體方法見式(15)～(18)。

(15)

(16)

(17)

(18)

對式(12)中矩陣T的每行進行Weibull參數估計，將得到每個BE失效時間或修復時間的分布參數，其中形狀參數行列式B=(β1，…，βm)T，尺度參數行列式N=(η1，…，ηm)T。

3煤礦突水DFT的Monte Carlo仿真

依據BE概率、Weibull分布參數及統計數據，對我國煤礦突水DFT進行Monte Carlo仿真，最后比較煤礦突水DFT仿真結果、FTA計算結果和近年來統計數據。

3.1煤礦突水DFT的仿真方法

在煤礦突水DFT的Monte Carlo仿真中，抽取的隨機數有兩種：①不同時間失效/修復概率不同的動態BE描述的事故處于發生/修復狀態下的最小時間，該隨機數按Weibull分布抽取；②失效概率不隨時間變化的靜態BE描述的事故狀態隨機數(0，發生；1，不發生)，每次仿真時按固定的概率抽取。Monte Carlo仿真時間步長分為固定時間步長和可變時間步長，本文中選擇固定時間步長隨仿真過程而推進，編寫程序管理模擬時間T仿真系統的實際時間。

3.2BE概率數據及處理

環境因素隱患BE1、BE2、BE3、BE4、BE5、BE9、BE10可通過統計資料取得，但其他隱患數據無法通過歷史統計數據獲得，根據SEQ和PAND的邏輯關系，利用前文第2章的處理方法，得到BE在不同時間點失效或修復概率。

依據模糊可靠性評價方法，首先邀請10位煤礦防治水專家對煤礦突水DFT模型中部分BE在單位時間(小時，h)內失效及修復可能性進行評價，得到專家模糊評價表(表1、表2)。

根據式(2)～(10)，依據表1、表2求得BE的失效或修復概率矩陣P。以隱患BE11為例，獲得專家對“鉆孔封閉不良”的失效概率分布，見式(19)。

p(0.77)，p(0.77)，p(0.95)，p(0.95)，

p(0.95)，p(0.95)]

=[0.1373，0.1373，0.029，0.029，0.029，0.029，

0.1373，0.1373，0.1373，0.1373]

(19)

根據矩陣P計算BE的平均失效概率，結果見表3。

在DFT模型分析中，需要獲得BE的失效或修復時間分布。根據式(12)～(18)，利用矩陣，計算專家對BE的失效和修復時間分布評估參數。對式(19)進行上述處理，得到“鉆孔封閉不良”事件失效時間Weibull分布參數為β=1.3211，η=19.7467。

表1　失效可能性專家模糊評價語言

表2　修復可能性專家模糊評價語言

表3　靜態BE的失效概率統計及計算值

DFT模型一級出現了SEQ，因此對煤礦突水全部隱患所描述的BE均應考慮不同時間點的失效概率。根據SEQ及PAND的邏輯關系，在對煤礦突水DFT進行Monte Carlo仿真時，只需要掌握BE15、BE16、BE17、BE18在不同時間點的修復率。按照以上分析，得到煤礦突水BE失效時間和修復時間Weibull分布參數，見表4。

在表4中，不同BE在不同時間點的失效概率和修復概率相差很大，而以往的靜態FTA以事件發生概率恒定為前提假設，是不合系統實際情況的，DFT模型考慮了時間因素對系統的動態影響，更和服實際情況。礦山系統開始基建后，上述BE的發生概率開始變化，因此基建開始時間點是BE所描述隱患的潛伏起始點，也是Weibull分布的原點。

表4　礦井突水BE中失效及修復時間的Weibull分布參數

3.3煤礦突水風險概率評價

依據國家煤礦安全監督局2013年的《全國煤礦水文地質類型劃分結果資料匯編》及前文第2章介紹的評價方法，統計得到地表水積聚等事件發生的概率，如表3所示。

BE所描述的事故發生概率自其存在礦山開采系統時開始增加，如開始鉆孔工作后，鉆孔封閉不良事故發生概率增加，礦山開始基建工作后，與突水DFT模型中BE相關的采礦工作都將開始，煤礦突水隱患狀態將改變，BE失效概率開始增加，進而導致BE修復概率的變化，由此將Monte Carlo仿真的時間起點選取為煤礦基建工作開始點。

Monte Carlo仿真結果見圖2，從圖2中可以看到，在煤礦動工的22000個小時(2.5年)之前，煤礦突水事故發生概率不斷上升，2.5年后事故發生概率在0.001上下波動。我國煤礦在理想的環境中基建計劃時間一般為2～3年，礦山建設期間，由于地質資料的不詳細、安全管理的不完善等，突水事故發生概率不斷增加，礦山運行階段，各個生產線、管理等方面都趨于穩定狀態，但容易出現管理疲勞等現象，突水事故發生概率在0.001左右波動。

根據FTA頂事件發生概率計算公式，計算得到煤礦突水的發生概率為0.000027。利用國家煤礦安全監督局公布的突水事故年發生次數和我國煤礦礦井總個數，計算得到近兩年煤礦突水事故的年發生率為0.0017325，見圖3，可以發現DFT仿真結果明顯較FTA計算結果符合實際情況。

圖2　煤礦突水事故發生概率分布

圖3　煤礦突水發生概率DFT與FTA結果對比

由圖3可知，FTA計算結果遠遠小于DFT仿真結果，這是由于FTA靜態邏輯門是DFT動態邏輯門的特殊狀態，反而言之，FTA的與門要求底事件在同一時間點發生，頂事件才會發生，而在實際情況中，底事件在某一段時間內處于失效狀態，頂事件就會發生，如DFT的PAND門，因此導致FTA計算結果較DFT仿真結果遠遠小于實際情況。

4結論

本文建立了煤礦突水DFT分析模型，對突水事故發生概率進行Monte Carlo仿真。

1)對無法獲得失效率及修復率數據的BE，可采用三角模糊數評價及反模糊化方法，該方法運用模糊集理論，盡可能消除了專家主觀模糊性的誤差因素。

2)根據DFT分析模型，對我國煤礦突水事故發生概率進行Monte Carlo仿真，發現在現階段煤礦防治水工作計劃處于良好實施的情況下，煤礦突水發生概率將穩定于0.001，而FTA分析計算煤礦突水發生概率為0.000027，據統計計算，2012年我國煤礦突水發生率為0.001848，2013年為0.001617，由此顯示DFT較FTA更加符合我國煤礦突水風險評估的實際情況。

參考文獻

[1]武強，趙蘇啟，孫文潔，等.中國煤礦水文地質類型劃分與特征分析[J].煤炭學報，2013，38(6)：901-905.

[2]高岳.用事故樹法確定煤礦防治水安全檢查表[J].煤田地質與勘探，2010，38(6)：58-61.

[3]許江濤，鄧寅生，文廣超，等.事故樹分析法在礦井水害防治中的應用[J].西安科技大學學報，2009，29(4)：405-408.

[4]王長申，孫亞軍，杭遠.基于事故樹分析的煤礦潛在突水危險評價研究[J].巖石力學與工程學報，2009，28(2)：398-305.

[5]N.Khakzad，F.Khan，P.Amyotte.Dynamic risk analysis using bow-tie approach[J].Reliability Engineering and System Safety，2012(104)：36-44.

[6]K.Durga Rao，V.Gopika，V.V.S.Sanyasi Rao，et al.Dynamic fault tree analysis using Monte Carlo simulation in probabilistic safety assessment[J].Reliability Engineering and System Safety，2009(94)：872-883.

[7]N.Khakzad，F.Khan，P.Amyotte.Dynamic safety analysis of process systems by mapping bow-tie into Bayesian network[J].Process Safety and Environmental Protection，2013(91)：46-53.

[8]A.Suprasad，D.Glenn，H.Eileen.A new approach to solve dynamic fault trees[C].In：Annual IEEE reliability and maintainability symposium，2003：374-379.

[9]王長申，孫亞軍，杭遠.安全檢查表法評價中小煤礦潛在突水危險性[J].采礦與安全工程學報，2009(9)：297- 303.

[10]J.Chen，R.Yang.Analysis of mine water inrush accident based on FTA[J].Procedia Environmental Sciences，2011(11)：1150-1554.

[11]J.H.Purda，J.Lu，G.Zhang，et al.A fuzzy reliability assessment of basic events of fault trees through qualitative data processing[J].Fuzzy Sets and Systems，2014(243)：50-69.

[12]D.Wang，P.Zhang，L.Chen.Fuzzy fault tree analysis for fire and explosion of crude oil tanks[J].Journal of Loss Prevention in the Process Industries，2013(26)：190-1398.

[13]孔峰.模糊多屬性決策理論方法及應用[M].北京：中國農業科學技術出版社，2008：79-85.

[14]許金華，孫德強，范英，等.基于FTA仿真的三高氣田事故風險概率[J].系統工程理論與實踐，2012，3(4)：877-884.

[15]R.Jiang，D.N.P.Murthy.A study of Weibull shape parameter：Properties and significance[J].Reliability Engineering and System Safety，2011，96(12)：1619-1626.

[16]J.Maiti，V.Khanzode，P.Ray.Severity analysis of Indian coal mine accidents-a retrospective study for 100 years[J].Safety Science，2009，47(7)：1033-1042.

[17]V.V.Khanzode，J.Maiti，P.K.Ray.A methodology for evaluation and monitoring of recurring hazards in underground coal mining[J].Safety Science，2011(49)：1173-1179.

[18]S.A.Akdag，A.Dinler.A new method to estimate Weibull parameters for wind energy applications[J].Energy Conversion and Management，2009(50)：1761-1766.

[19]H.Saleh，A.Abou El-Azm Aly，S.Abdel-Hady.Assessment of different methods used to estimate Weibull distribution parameters for wind speed in Zafarana wind farm，Suez Gulf，Egypt[J].Energy，2012(44)：710-719.

[20]P.A.C.Rocha，R.C.de Sousa，C.F.de Andrade，et al.Comparison of seven numerical methods for determining Weibull parameters for wind energy generation in the northeast region of Brazil[J].Applied Energy，2012(89)：395-400.

收稿日期：2015-11-18

作者簡介：王華鎣(1990-)，女，碩士，2013年畢業于武漢理工大學礦物資源專業，現就讀于北京科技大學采礦系。E-mail：why_ustb@163.com。 通訊作者：趙怡晴(1980-)，女，講師，主要從事資源經濟與管理、礦山安全預警、礦山安全標準方面的研究工作。E-mail：yiqingzhao@foxmail.com。

中圖分類號：TD745

文獻標識碼：A

文章編號：1004-4051(2016)07-0102-07

DFT based evaluation of risk probabilities for underground coal mine water inrush

WANG Hua-ying1，2，LI Zhong-xue1，2，ZHAO Yi-qing1，2，LEI Wei1，2，SHU Yang1，2

(1.Schod of Civil & Environmental Engineering，University of Science and Technology Beijing，Beijing 100083，China；2.Key Laboratory of High Efficiency Mining and Safety for Metal Mines(Ministry of Education)，University of Science and Technology Beijing，Beijing 100083，China)

Abstract:Fault tree analysis(FTA) has been widely utilized as a tool for water inrush risk assessment in coal mines，and dynamic fault tree analysis(DFT) is still being studying，which takes the fact into consideration that the probability of basic events(BE) and Boolean logic will be changing.While one problem is that it’s difficult to obtain failure and repair probability data at different time when simulates DFT.This study models DFT of water inrush in Chinese coal mine and invites experts to estimate the failure and repair probability of BE，then triangular fuzzy numbers are adopted to quantify experts’ estimate.Using defuzzification and moment method，the shape and scale parameters of Weibull distribution about failure and repair time are achieved.The result of Monte Carlo simulation for DFT illustrates that the probability of underground coal mine water inrush is about 0.001 in 2.5 years since the mine begin its capital construction.

Key words：water inrush in coal mine；dynamic fault tree analysis；triangular fuzzy number