基于拓撲網絡算法的煤礦瓦斯爆炸風險度量

成連華，郭阿娟，劉 黎，郭慧敏

(1.西安科技大學 安全科學與工程學院，陜西 西安 710054；2.西安科技大學 西部礦井瓦斯智能抽采工程研究中心，陜西 西安 710054)

0 引 言

近年來，中國煤炭行業安全生產狀況持續好轉，但隨著煤炭開采向深部轉移，煤層瓦斯含量急劇升高，瓦斯爆炸風險依然存在。國家礦山安全監察局統計數據顯示，2010—2020年煤礦事故死亡10 203人，因瓦斯爆炸事故死亡2 762人，約占煤礦事故死亡人數的27.10%，瓦斯爆炸風險嚴重威脅工人的生命安全。因此，對瓦斯爆炸風險演化過程定量分析，對制定預防措施具有重要現實意義。

目前國內外學者圍繞煤礦瓦斯爆炸風險評價展開了大量研究。LI等基于改進層次分析法、貝葉斯網絡等方法，評估了煤礦瓦斯爆炸風險[1-2]；TONG等運用概率風險評估模型、蒙特卡羅法分析煤礦瓦斯爆炸事故中礦工不安全行為所造成的事故風險[3-4]；田水承等運用文獻分析法和熵權法確定礦工不安全狀態評價指標綜合權重，并進行實例分析[5]；王婉青和李紅霞分別基于IAHP方法、模糊綜合評價法分析了瓦斯爆炸事故風險中的關鍵風險，同時對事故后果提出預防措施[6-7]。但大多數風險評價方法集中于對事故結果進行分析，風險演化過程的定量分析有待進一步研究。拓撲網絡是一種用來分析復雜系統相互作用及演化的工具，由于其契合復雜系統的特點，在城市交通[8-9]、電力系統[10]、房建事故[11]等諸多領域得到應用。如ZHOU等構建地鐵施工事故網絡(SCAN)探究26種事故之間的聯系[12]；丁明等運用動態故障樹理論，提出停電事故鏈網絡模型，并評估支路的關鍵風險程度[13]；周新宇等構建房建事故行為致險網絡模型，解析了人的不安全行為與事故間的作用機理[14]。然而當前拓撲網絡應用于煤礦瓦斯爆炸風險分析研究主要集中于無向無權網絡，對于有向加權網絡研究還處于探索階段。

鑒于此，運用拓撲網絡，將影響瓦斯爆炸事故的風險因素作為節點，因素之間的演化關系作為連接邊，通過計算風險演化過程中直接和間接傳播概率，構建煤礦瓦斯爆炸風險拓撲網絡度量模型。同時對節點出入度、聚類系數及最短路徑等拓撲參數進行分析，確定事故演化最短路徑和關鍵風險因素，以期為煤礦瓦斯爆炸風險量化分析提供參考。

1 瓦斯爆炸風險演化分析

1.1 瓦斯爆炸事故案例收集

以國家礦山安全監察局官網中統計的煤礦瓦斯爆炸事故調查報告作為數據來源，通過搜集2010—2020年的瓦斯爆炸事故調查報告，梳理了152起瓦斯爆炸事故案例，主要包含事故基本情況、事故經過、事故原因等基本信息。

1.2 瓦斯爆炸風險因素分析

基于層次全息建模(hierarchical holographic modeling，HHM)方法[15]從多維度、多層次進行事故過程識別。結合安全系統工程理論從人、物、管理、環境、技術方法5個維度進行分類，歸納分析出67個風險因素、2個風險事件。

統計瓦斯爆炸風險出現頻次，按公式(1)計算發生概率Pi。瓦斯爆炸事故Ck中如果出現風險因素類別i，則Ck=1，否則Ck=0。其中k為某個瓦斯爆炸事故案例；k={1,2,…,n}，n=152。

(1)

由于在一次瓦斯爆炸事故中可能存在多個風險，因此概率和大于1。統計分析瓦斯爆炸風險結果見表1。

表1 風險因素分析與編碼Table 1 Analysis and coding of risk factor

1.3 瓦斯爆炸風險演化路徑構建

任何一起事故過程都可用一個或幾個事故鏈描述。在瓦斯爆炸風險演化過程中，當風險節點狀態偏離時，風險累積達到風險閥值，觸發整個系統發生連鎖反應，進而導致事故發生。

依照事故發展流程，確定風險節點的時序關系和邏輯順序，將風險節點依序相連形成事故鏈[16-17]，實現對瓦斯爆炸風險演化過程的描述。以山西平遙二畝溝煤業瓦斯爆炸事故為例進行分析，對事故直接原因①和②進行分析，可先后提取出D17與B22的2個風險因素。事故經過表明：炮采組工人張某(無爆破工特種作業人員證件)未執行“一炮三檢”和“三人連鎖爆破”制度違章爆破，爆破明火引爆9103工作面采空區涌入瓦斯，發生瓦斯爆炸，可先后提取出A16，A17和D21的3個風險因素和F11這1個風險事件。最終從該案例中提取出8個風險因素、1個風險事件和4條事故鏈，見表2。

表2 瓦斯爆炸事故信息描述及事故鏈構建Table 2 Information description of gas explosion accident and construction of accident chain

分析152起煤礦瓦斯爆炸事故案例，將含有共同風險因素的事故鏈融合到一個全局網絡中[18]，運用系統動力學軟件構建煤礦瓦斯爆炸風險演化路徑，如圖1所示。

圖1 煤礦瓦斯爆炸風險演化路徑Fig.1 Risk evolution path of gas explosion in coal mine

2 瓦斯爆炸風險度量模型構建

2.1 風險網絡連邊計算

引入交互作用矩陣(factor interaction matrix，FIM)理論計算風險連邊權重，以分析系統影響因素的交互作用關系，如圖2所示。

圖2 因子交互作用矩陣

因子交互矩陣共有n個風險因素，主對角上各風險因素表示為Pi(i=1,2,…,n)，非對角上方區域表示對角線上其他風險因素對相應Pi的影響xij，非對角下方區域表示對角線上Pi對其他各風險因素的影響xji。每一行對應的數值相加表示為Ci，每一列的對應數值相加表示為Ei，則風險網絡連邊權重wi的計算表達式為

(2)

由于文章篇幅所限，以安全檢查不到位→線頭外露→電路老化→電氣火花→引爆瓦斯(D19→B14→B11→B22→F11)風險演化路徑為例進行計算，由FIM理論計算可得風險網絡連邊權重wi=(0.55，0.48，0.35，1.08)。

2.2 風險網絡迭代加權

設計結構矩陣(design structure matrix，DSM)是一種用來表示系統中元素及相互作用的結構化建模方法，反映風險傳播網絡中一個元素通過所有可能路徑影響另一個元素的集成風險傳播概率，可直觀對復雜系統進行可視化分析[19-20]。筆者將迭代加權思想與DSM相結合，構建集成風險傳播DSM矩陣。

實際上，風險傳播不僅存在直接傳播，也存在多條路徑的間接傳播關系。R_DSM(Ri,Rj)為Rj對Ri的直接傳播概率(傳播路徑是1階)。

SP(1)(Ri,Rj)=R_DSM(Ri,Rj)

(3)

Rj對Ri的2階路徑風險傳播(經過中間風險因素Rp)概率為

(4)

Rj對Ri的3階路徑風險傳播(經過中間風險因素Rp和Rq)概率為

(5)

風險演化過程中可能存在多條1階、2階或高階傳播路徑[21]。綜合考慮所有風險傳播路徑，得到Rj對Ri的集成風險傳播概率CP為

(6)

以上述風險網絡連邊計算的風險演化路徑(D19→B14→B11→B22→F11)為例進行分析，RD19直接影響下一階RB14，通過RB14間接影響RB11(路徑長度為2)，如圖3(a)所示。由于RD19導致RB14發生概率為0.55，而RB14導致下一階RB11發生概率為0.48，由式(4)計算得RD19與RB11風險間接(2階)傳播概率為：R_DSM(RB11,RD19)=0.48×0.55=0.264。演化路徑R1→R2→…→R5集成傳播概率為：CP(R5,R1)=1-(1-0.26×0.35)(1-0.55×0.48×1.08)=0.35。如圖3(b)所示。

圖3 風險集成迭代加和概率Fig.3 Risk-integrated iterative summation probability

2.3 風險拓撲網絡形成

將風險迭代加和權重轉化成CSV文件，并運用Gephi 0.9.2軟件構建煤礦瓦斯爆炸風險加權拓撲網絡，如圖4所示。

圖4 風險加權拓撲網絡模型Fig.4 Risk-weighted topological network model

從圖4可以看出，煤礦瓦斯爆炸風險拓撲網絡模型是一個完整的連通網絡。圓圈代表67個風險節點，連線代表214條網絡連邊，連線的粗細程度代表節點間的集成傳播概率。如B25線條最粗，可直觀看出與該風險鏈接路徑最多，風險權值變化最大。通過定量描述瓦斯爆炸風險加權網絡中風險因素之間的復雜關系，為煤礦瓦斯爆炸風險度量提供了一個新的視角。

3 拓撲網絡風險度量分析

3.1 風險節點度分析

節點度分析是網絡拓撲結構最重要的指標之一，反映了節點在網絡中的直接影響力和重要性[22]。節點的度值越大，表明節點越重要。節點的度可表示為

(7)

通過Gephi 0.9.2軟件分析得到節點度分布如圖5所示，橫坐標為度貢獻值，縱坐標為度分布數量。平均度為2.315，表明每個風險因素至少與2個不安全因素存在鏈接關系。

圖5 節點度分布Fig.5 Distribution of node degree

選取風險節點中出、入度數排前的節點共20個，度分布如圖6所示。從圖6可以看出，B20，B34的入度和出度最高，表明誘發該風險的途徑較多，因此對該風險進行有效控制，抑制瓦斯爆炸發生的效果最佳。同理，A12，B17，B24等度分布較高的風險節點也是風險防控中需重點關注對象。

圖6 風險節點出、入度統計Fig.6 Outbound and inbound statistics of risk nodes

3.2 風險聚類分析

聚類系數是網絡節點中鄰接節點之間相互連接的比例，反映了網絡節點集團化程度。對于與k個節點相連的節點i而言，它最多可與這些鄰居節點構成k(k-1)/2條連邊，則實際存在連邊Ei與最多可能連邊個數之比是節點i的聚類系數Ci。

Ci=2Ei/k(k-1)

(8)

通過Gephi 0.9.2軟件分析得到風險網絡模型的聚類系數為0.212，表明在整個網絡演化中，導致瓦斯爆炸的成因復雜，但大部分風險因素之間的相互影響并不明顯。

篩選節點加權聚類系數排前的節點共20個，如圖7所示。從圖7可以看出，B20，B33聚類系數最高，表明該節點集團化嚴重，若對該節點及鄰接節點構成的子系統采取斷鏈控制措施，可有效減緩事故發生。

圖7 節點加權聚類系數Fig.7 Node-weighted clustering coefficients

3.3 最短路徑分析

節點度分布及聚類系數分別分析了風險節點的重要程度及聚類程度，但未考慮上一節點引發后續節點的風險傳遞。通過Gephi 0.9.2軟件分析得到平均最短路徑為2.506，表明瓦斯爆炸風險拓撲網絡中一個節點狀態發生變化平均僅通過2.506步就能夠引起鄰接其它節點發生變化。

Dijkstra算法是荷蘭著名計算機科學家E.W.Dijkstra于1959年提出關于圖論中求最短路徑常用算法[23]。借鑒Dijkstra算法，編制計算程序，求解瓦斯爆炸風險最短路徑，得到風險節點的最短路徑及其最短路徑值(風險權值)，見表3。

表3 瓦斯爆炸風險演化最短路徑分析Table 3 Analysis of the shortest path in risk evolution of gas explosion

通過分析表3可知，風險演化路徑越短，權重值越大，但在風險路徑傳遞過程中風險權值逐漸減少，說明該過程存在風險損失；B31→F12風險路徑最短，風險權重變化顯著，表明B31(沖擊波蔓延)對瓦斯爆炸風險演化影響較大；雖然B11→B16→B25、A20→B20→F11和B17→B16→B25這3條最短演化路徑都是由3個風險因素構成，但初始風險因素不同，風險演化路徑不同，所以風險權值也不同。此外，風險演化最短路徑也從側面反映了瓦斯爆炸事故發生的可能性。在實際瓦斯爆炸風險管控工作中，優先對最短風險路徑上的風險因素加以控制，可以提高解決事故的效率。

4 結 論

1)通過統計分析2010—2020年共計152起瓦斯爆炸事故案例，采用HHM方法從人、物、管理、環境、技術方法5個維度進行風險分析，得到67個風險因素，形成214條瓦斯爆炸風險演化路徑鏈。

2)以拓撲網絡算法為基礎，引入DSM對瓦斯爆炸風險進行迭代加權計算，并運用Gephi 0.9.2軟件構建瓦斯爆炸風險拓撲網絡度量模型，實現瓦斯爆炸風險演化過程的量化分析，以揭示事故發生的潛在特征。

3)通過分析瓦斯爆炸風險有向帶權網絡參數，發現風險節點聚集程度較低，呈現明顯小世界網絡特征，表明在整個風險演化過程中大部分風險因素之間相互影響不明顯，存在顯著傳遞關系。

4)運用Dijkstra算法搜尋風險演化最短路徑，得到7條最短風險演化路徑和3個關鍵風險因素，優先對瓦斯爆炸最短風險路徑上關鍵風險加以管控，對預防瓦斯爆炸事故發生、保障企業安全運行具有現實意義。