沖擊地壓“雙驅動”智能預警架構與工程應用

陳 結,杜俊生,蒲源源 ，姜德義,齊慶新

(1.重慶大學 煤礦災害動力學與控制國家重點實驗室，重慶 400044; 2. 重慶大學 資源與安全學院，重慶 400044; 3. 煤炭科學研究總院 深部開采與沖擊地壓研究院，北京 100013)

煤炭是我國的主要能源，即使在“雙碳”戰略的大背景下，直到2030年煤炭在我國一次能源結構中占比仍將達到46%。淺部煤炭資源日益枯竭使我國煤礦開采逐漸轉入深部。目前超過一半的煤炭資源儲存在1 000 m以深的地層中。為保證持續發展所需的能源供給，勢必加快向深部開采的步伐。深部開采過程中，煤巖體應力環境比淺部更加惡劣，具有沖擊危險性的礦井數量增多，嚴重威脅煤礦安全生產。近年來，我國主要產煤基地均有嚴重的沖擊地壓事故發生，造成了巨大的人員財產損失。因此，沖擊地壓的防災減災問題研究仍是煤炭行業剛性需求及科研熱點。

沖擊地壓預測預警是防災減災的重要環節，可以為沖擊地壓現場治理提供指導，其核心目標是在災害發生之前準確判定災害的時間、區域與危險性等級。目前絕大多數沖擊地壓預測方法是對采動區域的危險等級進行評估。沖擊地壓危險性預測初期，主要是通過采動區域狀態評價指標參數對區域危險等級進行靜態評價。目前主要方法為綜合指數法、可能性指標法、地質動力區劃法、應力集中評價法及多因素耦合法等。通過液壓支架壓力監測、巷道錨桿監測、鉆屑法監測及煤巖體應力監測對沖擊地壓危險性進行動態分析。隨后，微震、地音、電磁輻射、聲發射等地球物理方法廣泛用于沖擊地壓危險性動態評價。相比于靜態評價法，地球物理法可以實時捕捉煤巖體破裂、應力集中等現場信息，具有更好的現場工程應用價值。由于沖擊地壓災害的致災因素具有多元性和不確定性，通常將多參量指標結合起來進行更為全面的監測預警。例如，竇林名等綜合采場應力監測與微震監測，形成了“震動場-應力狀態”耦合一體化監測預警系統，利用獲得的波速(梯度)異常系數與煤體應力增量對沖擊危險性進行實時預警。姜福興等以巷道掘進期間的微震、地音和應力信息以及獲取的真實地質資料構建了回采工作面沖擊危險性動態評估方法，對區域危險性進行評價。劉少虹等建立了一種基于電磁波和地震波CT聯合探測沖擊危險性的評價方法，實現了掘進工作面及回采工作面大面積沖擊危險性等級評價。目前沖擊地壓危險性的動、靜態這種物理模型的評價方法已有相關應用，但通過動-靜態協同的多參量物理評價方法進行沖擊危險性聯合監測預警的工程應用并不常見。

沖擊地壓預測預警更為關鍵的問題是災害時間與區域的量化預測。目前存在的主要問題是缺乏現場監測參量與災害事件、區域之間的映射關系，導致時間與區域的預測很難量化。學者對于沖擊地壓時間與區域的研究主要集中在災害發生以后，對孕災過程中能量、應力、微震等參量進行回溯分析，研究各參量的時空演化規律，其中絕大部分屬于定性或者規律性描述。例如，李康等研究了臨空巷道沖擊地壓耦合致災因素的時空演化過程，提出了巷道靜載應力、動載應力在孕災過程中的時空演化規律。趙志鵬等研究了近直立煤層沖擊地壓發生前微震事件的時空演化規律，發現微震頻次與能量在沖擊地壓發生前的異常變化，定性描述了微震能量偏差值與災害發生的關系。張月征等利用鉆孔應變數據的分形盒維數變化表征了沖擊地壓時空演化過程，發現沖擊地壓發生前，相關區域應變的分形盒維數會出現“升高→降低”變化。HE等研究了某急傾斜工作面發生沖擊地壓前的微震事件空間分布，揭示了孕災過程中災害狀態量的變化情況。然而，上述研究很難直接使用孕災過程中的現場參量前溯預測災害的時間與區域，即無法構建監測參數與災害時間區域的定量映射模型，實現真正意義上的災害預測。

綜上所述，傳統的數學物理模型很難量化監測參數與災害時間區域的定量映射關系，主要原因包括：① 模型泛化性差。沖擊地壓的數學物理預測模型建立在災害機理(假說)上，而目前尚無統一的沖擊地壓災害機理。基于機理的數學物理模型僅適用于特定條件下的沖擊地壓預測；② 模型魯棒性差。數學物理預測模型通常考慮有限幾個影響因素，模型中參數所占權重大，一旦有參數獲取不精確，將直接影響到模型的預測效果；③ 預測動態性差。采掘過程中沖擊地壓主控因素會動態變化，而顯式的數學物理預測模型很難實時反映參數變化，模型的預測動態性無法保證。但是，以深度學習、知識圖譜等為代表的數據驅動方法可以很好地克服數學物理模型的以上弱點，其不關注先驗機理知識，不關注參數間的耦合關系，模型可以動態更新，非常適用于沖擊地壓時間與區域的定量預測。近2年來，已有不少學者將機器學習應用于沖擊地壓的研究，包括沖擊危險性評價、沖擊地壓微震信號處理等方面。然而，目前還沒有直接利用數據驅動的方法進行沖擊地壓時間區域定量預測的報道。

鑒于此，筆者提出一種全新的沖擊地壓“雙驅動”(物理驅動及數據驅動)智能預警架構：在物理驅動框架下，建立基于貝葉斯概率的動靜態協同物理評價模型，對沖擊地壓危險性等級進行實時評價；在數據驅動框架下，利用深度學習模型實現未來微震事件時間置的精準預測，繼而量化沖擊地壓的區域與時間。基于“雙驅動”架構，開發相應的工作面沖擊地壓智能預警平臺，實現預測結果的三維動態可視化，更加直觀地展示沖擊地壓孕災期間危險特性信息，方便工程技術人員準確、快速制訂防災減災措施。

1 沖擊地壓“雙驅動”智能預警技術架構

“雙驅動”智能預警架構基于集成綜合預警思想，利用煤巖體特性、開采參數、實時監測數據及歷史災害資料建立了沖擊地壓危險等級評價的物理驅動模型和災害時間、區域預測的數據驅動模型。其預警的整體流程如圖1所示。

圖1 沖擊地壓“雙驅動”智能預警構架

在物理驅動方面，利用動靜態協同算法，采用采區現場勘測得到煤巖體結構特性及采場空間特征信息及煤巖體力學參數，計算沖擊地壓靜態綜合指數和可能性指數，形成沖擊地壓的靜態評價體系；同時利用煤礦井下安全監測系統獲取的應力、微震參數及反演的波速異常系數，形成動態評價體系。然后利用貝葉斯概率模型對動靜態各指標進行危險等級的后驗概率計算。最后，根據沖擊地壓設置的等級閾值確定風險等級。

在數據驅動方面，利用現場微震在線監測系統獲取的微震事件屬性參數(時間、能量、坐標)作為深度學習模型的訓練集。充分考慮微震事件的長短程時序依賴規律，建立具有線性與非線性多模塊耦合的深度學習模型MSNet。將工作面連續的微震事件時間序列作為模型輸入，預測未來微震事件發生的時間及位置，繼而對沖擊地壓潛在的時間與區域進行定量預測。最終，通過物理和數據“雙驅動”模型實現沖擊地壓綜合智能預警。

最后，基于提出的沖擊地壓“雙驅動”預警架構，開發了相應的智能預警平臺，實現沖擊危險性等級、沖擊危險性區域及潛在災害時間的實時預警和三維動態展示。預警平臺主要分為5個模塊：現場數據獲取與傳輸、數據存儲與前處理、物理與數據驅動模型、實時預警與定量預測、三維顯示及警情發布。

圖2展示了雙驅動模型的沖擊地壓災害預警平臺的總體框架。從以下4個層面進行簡要描述：① 在現場數據獲取層面，獲取井下采場區域煤巖體賦存特征信息；利用物探及鉆探方法獲取采動區域及圍巖體周圍的相關結構特征；通過在線監測系統監測鉆孔應力、微震等實時數據，結合煤礦地質資料及歷史數據，形成完整的多源異構數據體系；② 在數據存儲與前處理層面，利用物聯網實時傳輸技術，將井下在線監測系統獲取的數據、災害區域的靜態性指標參數、地質參數及歷史數據等實時傳輸并匯集至指定的數據倉庫，通過數據去噪、刪除、清洗等形成可用數據；在物理與數據驅動模型層面，分別建立沖擊地壓危險性等級動靜態協同評價及微震事件定量預測模型；③ 在實時預警層面，根據“雙驅動”架構的輸出結果，確定沖擊地壓的等級、區域及時間；④ 在三維顯示及警情發布層面，借助可視化平臺對井下采區三維構造、監測數據及預測預警結果進行平臺顯示。當預警判別結果表示存在警情時，通過預警信息發布系統對結果進行多終端實時發布，并采取對應的應急避險預案和防控措施。

圖2 沖擊地壓綜合預警平臺總體框架

2 沖擊地壓危險性動靜態協同評價模型

2.1 靜態評價

綜合指數法是《防治煤礦沖擊地壓細則》推薦采用的方法。本研究從地質因素和開采技術兩方面進行綜合指數評價計算，并從應力和煤層沖擊傾向性進行可能性指標計算，基于綜合指數評價計算和可能性指標計算獲得區域靜態評價結果。

2.1.1 綜合指數評價計算

從地質因素和開采擾動2個方面的綜合指數來討論沖擊危險性靜態指標，地質因素的危險性指數的計算公式如下：

(1)

式中，為第個地質影響因素的實際指數；max為第個地質影響因素中最大的危險指數；為地質影響因素的數量。

同理，根據細則對沖擊地壓開采擾動因素危險性指數進行分析,計算公式如下：

(2)

式中，為第個開采擾動影響因素的實際指數；max為第個開采擾動影響因素中最大的危險指數；為開采擾動影響因素的數量。

最后，綜合地質因素及開采擾動因素指標可描述為

=max(,)

(3)

其中，為沖擊地壓綜合指數的分數值。沖擊地壓區域危險程度則通過進行分類，具體的等級劃分方法見表1。需要注意的是，本文后續的鉆孔應力、微震參數各指標及綜合評價的風險評級按照表1的劃分標準進行。

表1 靜態指標綜合指數評價標準

2.1.2 可能性指標計算

一般可能性指數分為2部分，第1部分為應力對沖擊地壓的貢獻，第2部分為煤層本身沖擊傾向性對沖擊地壓的貢獻。

應力對沖擊地壓事件的貢獻度為

(4)

式中，=[]，為應力集中系數，為巖石容重，N/m，為開采深度，m，[]為煤層單軸抗壓強度，MPa。

煤層沖擊傾向性對發生沖擊地壓的貢獻度W可表述為

(5)

式中，為沖擊傾向性的彈性能量指數。

沖擊地壓的可能性綜合指數為

(6)

沖擊地壓可能性指標的等級評價可按照表2的標準進行等級劃分。

表2 沖擊地壓可能性指數評價標準

2.2 基于動態應力監測的沖擊危險性評價

2.2.1 鉆孔應力指標

為掌握采場區域的應力狀態，對工作面上下巷鉆孔應力的數據進行分析。基于時間序列的某測點所反應的圍巖應力等級狀態通過式(7)表示。

=++

(7)

式中，為鉆孔應力綜合表征指標；為鉆孔應力值的危險指標；為應力增幅值的危險指標，即當前測點應力值與最初設點應力值之差；為應力增速值的危險指標；，，分別為，，的權重比，根據現場應力特征及專家意見，本文分別取0.5，0.3，0.2。

由式(8)計算得出。

(8)

式中，為測點的應力，MPa。

由式(9)計算得出。

(9)

式中，為測點的當前應力與開始統計的應力之差，MPa。

由式(10)計算得出。

(10)

式中，為測點當前應力與上一刻的應力差。

若式(8)～(10)的，，取值不大于0時，，，取值為0。

2.2.2 波速異常系數

為了探測煤層區域應力場的相對大小，利用地震CT對煤層監測區域的波速場進行反演，并基于反演結果計算波速異常系數()。

(11)

2.3 微震多維信息動態評價

(1)微震值。

根據古登堡和里克特研究得到的地震頻度和地震震級之間的冪率關系，學者們通過變換將地震震級和頻度聯系起來，建立如下方程：

lg(≥)=-

(12)

其中，和均為常數，為統計區域內地震活動水平的表征參數，為該區域內大小震級數量之間的一種比例關系參數，大震級事件相對較多時，則變小，可以一定程度上表征該區域的應力變化和煤巖體的破裂情況。

為研究微震震級與頻度的關系，通過最小二乘法擬合計算可得到，，其中

(13)

式中，為第檔的微震震級，在整個時間窗內，對每個事件能級進行大小排序，以最小的微震事件能級作為起始震級，震級檔間距Δ設為0.2；為第檔震級內的微震實際事件數；為事件窗口內震級的總檔數。

(2)微震()值。

微震的頻次、震級與微震活動具有正相關性。微震頻次和震級越大，沖擊危險性也越大。研究表明，()是一個能夠反映某一區域微震活動程度的定量參數，可以較為直觀反映微震活動的增強或平靜趨勢。()定義如下：

(14)

式中，為第個微震事件的能量；lg為第個微震事件的能級；為時間窗內微震事件總次數。

(3)累計差量指標。

微震累計差量是通過分析某時間段內微震事件率的累計差量變化來判定采區圍巖體的穩定性，設累計差量指標為，大事件累計差量指標為，事件率累計差量指標為，能率累計差量指標為。

根據馬雄忠等提出的聲發射累計差量指標算法來計算微震事件率累計差量指標，設時間窗內的微震事件率數據序列為,,…,-1,(≥4)，其中，為時序為時的事件率。

對序列項前后相減，得

令=max{Δ,0},=min{Δ,0},

將正差量的和設為，負差量的絕對值和設為，即

(15)

(4)等效能級參數。

某段時間內發生微震事件的能級采用如下計算方式：

(16)

式中，為時間窗內微震事件的總能量；為時間窗內第個事件的能量；為微震事件總數。

等效能級參數()是反映新的微震事件相對于正常水平的偏離程度的指標，具體公式如下：

(17)

2.4 基于貝葉斯方法的動靜態協同評價模型

根據前文提出的沖擊地壓靜態評價指標(，)、應力動態評價指標(，)及微震多維信息指標(，()，，)，形成了協同評價沖擊地壓的動靜態綜合指標體系。由于在統計時間區間內，各指標的計算值存在未知性。因此，首先通過歸一化的方法對各指標進行規范化處理，轉換為統一的正向指標。然后利用貝葉斯方法建立動靜態各指標參數的危險等級后驗概率關系，最后利用常規的加權平均模糊算法對建立的動靜態指標進行協同綜合評價。

貝葉斯概率的基本思想是：假定所研究的對象具有一定的認識，并且這些對象存在設定的先驗概率分布，然后通過新的樣本對設定的先驗分布進行不斷更新，從而得到后驗分布概率。本研究利用貝葉斯概率方法進行沖擊地壓多參量綜合危險性判別。貝葉斯概率計算公式為

(18)

式中，()為事件的先驗概率；(|)為事件的后驗概率；(|)()為調整系數。

將貝葉斯概率公式引入沖擊地壓評價時，式(18)改寫為

(19)

式中，為各評價指標計算值(=1,2,3,…,);為劃分的危險等級(=1,2,3,4)；為各指標設置的等級閾值。

則利用貝葉斯概率綜合評價沖擊地壓的計算步驟如下：

(1)計算()，即各評價指標參數值在沖擊危險性各等級的先驗概率。這里設定各指標計算結果處于每個風險等級的概率相等，即，(1)=(2)=(3)=(4)=14。

(2)計算(|)，概率采用幾何概率中距離計算值法，將各評價指標計算值與設定的危險等級評價標準之間距離的絕對值倒數進行計算可得

(20)

其中，=|-|。越小，則表示所測指標屬于對應等級的概率越大。

(3)根據式(19)，計算出(|)，即各評價指標參數值在沖擊危險性各等級的后驗概率。然后通過加權平均模糊綜合算法計算動靜態各指標在各危險等級下沖擊地壓綜合后驗概率，公式如下：

(21)

式中，為各動靜態指標的權重。

筆者基于“應力場-震動波場”監測預警方法，選擇微震多維信息、在線鉆孔應力及地震-CT波速異常系數作為應力場-震動場動態綜合指標體系，結合沖擊地壓靜態評價體系，形成沖擊地壓危險性動靜態協同綜合評價模型，其綜合評級流程如圖3所示。

圖3 沖擊地壓危險性等級綜合評價流程

3 微震事件定量預測數據驅動模型

3.1 模型基礎數據

筆者使用的微震數據來自某煤礦13230工作面。選取在2019-01-01—2019-12-31一年內工作面推進過程監測到的10 196次微震事件作為本文的基礎數據。微震事件預測是時間序列預測任務，模型的輸入是連續的微震事件序列，其中序列長度為模型超參數，筆者選取序列事件長度為12，即利用12個連續的微震事件預測未來的事件。模型的輸出是即將發生的微震事件序列，選擇模型要預測的事件序列長度為6。由于監測的微震事件在時間軸上的分布是不均勻的，因此，本研究為每個微震事件引入一個時移參量(單位為min)，表征其與前一個事件的時間差。MSNet模型輸入數據為長度為12的微震事件序列，每個事件包含空間位置(,,坐標)、能量、時移5個屬性，單個輸入樣本維度為12×5。MSNet模型輸出數據為長度為6的微震事件序列，事件屬性與輸入數據相同，單個輸出樣本維度為6×5。表3展示了工作面輸入微震事件時間、空間坐標、事件能量、時移的部分數據。

表3 工作面輸入微震事件屬性

3.2 數據驅動模型設計

微震事件在時間軸上的分布有長短程依賴特征，即一個微震事件或微震事件序列的發生不僅與之前短時間范圍內(如幾小時內)的微震事件有關，也會受到更長時間范圍(數天到數月)的微震事件影響。例如，在沖擊地壓發生前的一小段時間范圍會發生大量能量各異的微震事件，此為短程依賴；同時，頂板來壓等長周期事件也會影響微震事件的發生，此為長程依賴。因此，本研究建立可以同時獲取短程與長程的微震事件時序規律的MSNet模型，更精確地預測未來發生的微震事件。如圖4所示，MSNet的數據流向有兩大路徑，其中線性路徑是利用AR自回歸模型，輸入的微震事件序列直接預測未來的微震事件(序列)；非線性路徑則通過卷積層、循環神經網絡層、跳躍循環網絡等模塊，實現輸入的微震事件序列到未來微震事件序列的非線性映射。MSNet的最終預測結果是線性路徑結果與非線性結果的加和。下面依次介紹MSNet模型用到的幾個關鍵模塊，卷積模塊，循環神經網絡模塊，跳躍循環網絡模塊、自回歸模塊。

圖4 MSNet模型結構

3.2.1 MSNet非線性模塊

MSNet模型的非線性模塊主要包括卷積層、循環神經網絡及跳躍循環網絡3個模塊，通過這些模塊的綜合訓練實現微震事件的非線性預測。數據首先流入MSNet的卷積模塊，用于解析12個(單個序列長度)微震事件間的短程依賴關系以及事件5個屬性間的關系。如圖4所示，卷積核尺寸確定為6×5，卷積核數量為32個。卷積層激活函數選擇為ReLU(其表達式為ReLU()=max(0,))。式(22)展示了第個卷積核的運算過程，為卷積核矩陣，為輸入的時間序列矩陣，為偏置，*為矩陣的哈密頓算子。運算結果采用zero-padding保證輸入前后的維度一致。

=ReLU(*+)

(22)

(23)

MSNet設計一個跳躍循環網絡(skip-LSTM)，解決潛在的超長程依賴問題。流入LSTM單元的卷積層輸出不是連續時間序列，而是每隔一個特定周期(跳躍步長)選取一個數據點，即把當前隱藏單元與相鄰周期處于同一位置的隱藏單元鏈接起來作為LSTM的輸入，以擴展信息流的時間跨度。引入跳躍步長后，skip-LSTM的數據流向如式(24)所示。循環神經網絡模塊及跳躍循環網絡模塊的輸出拼接在一起形成一個長向量，然后通過一個隱藏層的全連接網絡得到最終的非線性路徑輸出。

(24)

式中，-和-分別為前時刻的隱藏狀態和記憶單元。

3.2.2 MSNet線性模塊

為了充分考慮事件之間的潛在線性關系，MSNet設計了一條平行于非線性路徑的自回歸線性路徑。假設各微震事件屬性向量之間具有線性關系，利用前次的微震事件預測第(+1)次事件。本研究中自回歸模型如式(25)所示，和為模型參數；為待預測事件的屬性向量。MSNet中的自回歸模塊利用單隱藏層的全連接網絡實現。非線性路徑結果與線性路徑結果進行向量加和，得到最終的模型預測結果。

(25)

3.3 模型訓練

深度學習模型通過訓練迭代更新模型參數，使損失函數(預測結果與真實值之間的差異)達到最小值。本研究中，模型的預測輸出是未來連續6次微震事件的屬性(維度為6×5的向量)，因此定義MSNet的損失函數是五維空間中2個長度為6的點序列之間的平均距離，損失函數如下

(26)

10 196次微震事件數據集人為劃分為訓練集與測試集，前者用于訓練模型使之達到最佳擬合效果，后者用于測試模型表現。考慮到數據樣本有限，數據集中不再單獨劃分驗證集用于模型超參數優化。選擇數據集中2019年1—10月的9 863次微震事件進行模型訓練，11月與12月的333次微震事件進行模型測試。為了盡可能多地利用微震數據，模型的輸入微震序列之間取步長為1，即第1個訓練/測試樣本輸入為～(+11)的12次連續微震事件，預測(+12)～(+17)的6次微震事件的坐標、能量、時移。本研究在構建訓練/測試樣本時，沒有樣本跨越劃分時間點(11月1日0時)，保證訓練樣本和測試樣本完全隔開，無數據無泄漏。

3.4 模型結果分析

在整個訓練過程中，采用小批量梯度下降法更新模型參數，訓練損失與測試損失最后都達到收斂，模型沒有出現過擬合現象。把樣本輸出的連續6個微震事件屬性分別看作6個長度為5的向量，利用余弦相似度來衡量屬性的預測值與真實值之間的相似性(值越接近于1，表示2個向量越相似，預測精度越高)。其計算方式如式(27)所示，其中為微震事件真實的屬性向量；為模型預測的屬性向量；為向量長度。

(27)

模型結果表明，模型對微震事件坐標與坐標的預測效果最好，測試樣本的平均余弦相似度分別為0.997與0.995。模型對微震事件坐標的預測結果稍差，明顯低于與的余弦相似度，其平均相似度為0.939。模型對事件時移的預測出現了明顯波動，其平均的預測余弦相似度為0.523，表明其預測值與真實值之間仍存在不錯的相關性，MSNet模型挖掘出了微震序列的時間規律。模型對微震事件能量幾乎沒能做出正確預測，超過一半測試樣本的余弦相似度徘徊在0附近，表明預測值與真實值之間幾乎沒有關聯。綜合分析，MSNet模型對微震事件的位置預測最準確，時間預測精度次之，對微震事件的能量預測精度較差。

圖5展示了模型訓練數據(1—10月的9 863次微震事件)的5個屬性進行[0,1]標準化之后的分布箱型圖。

圖5 模型訓練數據分布箱型圖

圖5中綠色菱形表示數據異常值，藍色圓點表示數據平均值，長方形上下邊緣分別代表屬性的上四分位數與下四分位數，紅色長線表示數據中位數。從圖5中可以發現微震事件的位置數據(,,)分布較為理想，中位數與平均值均處于數據集中區域，且與平均值相差不大，異常值并未過多影響數據分布。微震事件的時移數據分布稍差，數據集中在一個較小區域，有部分異常值超過了數據內限，影響了數據分布，但其中位數與平均值仍處于數據集中區域。微震事件能量分布出現了明顯的極端異常值，其平均值與中位數處于數據集中區，但是數據能級出現大范圍空白，直接影響了數據分布。這是由于工程現場絕大部分微震事件是小能量事件，超過1×10J的大能量事件占比極少。MSNet模型訓練輸入的絕大部分能量數據位于小能級，因此模型對大能量事件預測結果帶有明顯的小能量偏向性(bias)，能量數據不平衡(imbalance)極大影響了模型對事件能量的預測能力。

4 基于“雙驅動”模型的沖擊地壓智能預警平臺開發

4.1 平臺三維幾何模型

基于礦井采掘工程平面圖、礦井地貌圖、地質說明圖、早期勘探及測繪相關資料，利用3Dmax建模軟件通過虛擬三維空間構建三維幾何模型，能夠對回采工作面、掘進工作面、巷道布置、圍巖體等元素進行展示，準確表征斷層、巖性、煤層厚度變化等特征信息。如圖6所示，通過Unity 3D平臺技術建立了煤礦地質三維模型，包括整個礦區地上地下的三維空間結構、各采區的相對位置、巷道的地質構造及巖性等特征信息，實現了較為準確的三維表征。

圖6 采區三維幾何模型

4.2 “雙驅動”預警模塊實現流程

平臺搭建預警模塊的步驟如下：第1步，數據轉換。對監測系統獲得的數據進行數據轉換，把信息轉換成機器能夠接收的形式。第2步，數據分組。將數據按照不同類型進行分組，對煤巖數據、地質數據、回采數據、微震數據等進行分類。第3步，數據組織。把分好類的數據放到一個容器中，主要有數組、鏈表、列表等，筆者選擇的是列表容器。第4步，數據計算、存儲。通過C#語言對評價模型算法進行代碼編寫，實現對接入數據的計算，并對數據進行保存。第5步，數據顯示。根據需求將容器當中需要進行界面顯示的計算結果通過Unity3D平臺進行展示。

對于微震事件預測模塊，在訓練完畢以后進行在線部署，利用實時獲取的微震事件序列對未來事件的屬性進行動態預測。采用Intel公司的OpenVINO深度學習套件進行模型的在線部署。首先從Keras平臺導出TensorFlow凍結圖模型文件(*.pb文件)，使用OpenVINO套件中的Model Optimizer工具優化導出的凍結圖模型，獲得描述MSNet網絡拓撲結構的*.xml文件與儲存模型參數的*.bin文件，調用套件里面的推理引擎(Inference Engine)完成AI的推理計算。MSNet做預測時輸入的實時數據從AraWin微震數據處理軟件自帶數據庫中獲取，模型的預測結果以文本形式本地保存并被實時調用。

在數據處理和代碼化的一系列工作之后，實現最終的三維預警平臺需要進行如下操作：對所用到的應用程序接口、編寫的代碼、3D模型、3D模型材質以及應用程序需要展示界面的圖片進行打包集成，通過unity導出形成一個三維透明化的應用程序。通過數據收集及處理，得到評價可用數據，將數據按動靜態類別整理形成動靜態指標數據庫。預警軟件平臺調用數據庫數據，利用嵌入的評價模型對沖擊地壓綜合危險性進行等級評價并進行區域與時間的預警。圖7展示了具有數據處理及評價功能的三維預警平臺。三維平臺界面右邊主要顯示微震事件及應力數據，左邊顯示曲線變化。

圖7 三維預警平臺危險性等級展示

5 工程應用

現場應用地點為河南省義馬市耿村煤礦，其示范工作面煤層厚度為0～20.32 m，平均厚度為9.86 m，傾角8°～15°，平均10°，工作面斜長245 m，走向長為745 m，平均采深為630 m。近水平層理，裂隙發育，厚度10.7～17.8 m，煤層結構復雜，含夾矸3～5層，夾矸為灰黑色泥巖或細砂巖，夾矸累厚1.5 m，為較穩定煤層，煤層變化系數為12%。直接頂為泥巖，厚4.04～42.64 m，深灰色，致密，塊狀構造，含植物化石。煤層直接底為泥巖，厚度約2.8 m，基本底為泥巖砂巖交替層，厚度35 m以上，致密堅硬。基本底為泥巖、細～中砂巖和礫巖，厚0～35.7 m，平均13.24 m，煤層單向抗壓強度8.58 MPa，煤的彈性能指數2.54。

采用井下監測系統和預警系統對工作面沖擊地壓災害進行監測預警，如圖8所示。在工作面前方250 m左右區域的上、下巷外幫埋設應力計，淺孔埋深12 m，深孔埋深18 m，淺孔與深孔每組應力計的間距為0.5～1.0 m，每組20～25 m。微震拾震器布署在整個工作面上下巷道，間距為250 m左右。隨著工作面的推進，鉆孔應力計及微震拾震器不斷向前移動。通過應力數據、微震及結合地震CT形成的波速異常系數對工作面沖擊危險性進行實時動態評價。

圖8 工作面應力-微震在線監測設計

(1)應力監測結果。通過地震CT主動爆破的探測方式，獲得工作面監測區域相對應力情況，再利用產生的微震事件作為被動震源進行三維波速反演，形成三維速度梯度。根據反演的三維數據進行三維建模，形成三維數據云圖。同時結合監測區域鉆孔應力的變化特征，實現對工作面采動區域應力集中情況透明化顯示，如圖9所示。

圖9 工作面三維應力云圖

(2)微震監測結果。圖10展示了2019年監測的微震事件及部分微震多維指標的變化結果，其中圖10(c)，(d)顯示了2019-08-01—2019-08-07微震多維指標的變化曲線。從微震事件和指標對比來看，微震出現較為明顯的低值異常，與微震事件能量具有負相關。()和等效能級參數指標的變化規律具有相似性，同累計差量指標具有明顯高值異常。

圖10 監測的微震事件及微震多維指標變化曲線

(3)沖擊危險性等級預警結果。當平臺接入評價指標數據庫后，某個指標參數的改變都會使其區域災害危險性的結果及時更新。評價指標參數一旦發生改變，平臺則立刻對指標進行重新計算，實時更新綜合等級結果。

通過預警平臺對示范工作面2019年的數據進行現場測試，對比現場實際危險情況，確定預警模型的準確率。根據示范工作面的現場防沖日志報告，在2019年，共有25 d具有危險性，均屬于弱危險性。根據示范工作面2019年現場實際的動力現象和工人勘測結果及管理統計來看，沒有特別明顯的支架壓力增大、卡鉆等中等及強沖擊危險性現象，頂底板、支架未見明顯變形，時有輕微煤爆或煤炮聲等。根據《防治煤礦沖擊地壓細則》并結合現場專家意見，認為該示范工作面2019年的防沖日志具有參考意義和可比性。

通過平臺的預警結果表明，有22 d的預警結果同礦上的評價等級及實際危險情況是相符合的，有3 d出現了不一致，其預警平臺的評價結果均顯示為中等危險，都高于該礦的防沖日報結果，具體日期見表4。因此，在沖擊危險性預警方面，平臺的準確率為88%。需要說明的是，該模型準確率驗證的前提是防沖日志結果和實際現場情況是吻合的。評價結果產生偏差的原因可能是各指標結果的綜合疊加計算引起的，可通過不斷優化各指標的權重比和指標的隸屬關系來提高準確率。

表4 沖擊危險性平臺預測結果與現場結果對比

(4)雙驅動預警結果。結合物理驅動與數據驅動的預測輸出，形成沖擊地壓等級、時間、區域的綜合預警結果。圖11為沖擊地壓“雙驅動”綜合預警結果，圖中黃色和淡藍色圓圈分別為平臺確定的沖擊地壓危險性等級區域和預測的下一時刻潛在的危險區域，此時的危險性等級為安全。隨著采掘的推進，危險性等級區域將動態化地涵蓋整個回采工作面。危險預測區域確定的方式是，以預測的6個微震事件中能量最大的事件位置為圓心畫圓，其邊緣通過距離圓心最遠的微震事件位置，且危險性沿徑向遞減。該平臺一直保持顯示最近12個連續微震事件，同時預測接下來6個微震事件。當微震系統監測到最新的微震事件時，平臺右側的數據將會動態更新，并同步獲得最新的預警結果。沖擊地壓時間預測的方式是，首先根據本礦井沖擊地壓歷史事故發生前監測到的微震事件最大能量值作為判斷閾值，一旦預測的微震事件能量超過該值，則觸發沖擊地壓時間預警。

圖11 沖擊地壓”雙驅動”預警三維平臺

6 結論與展望

(1)建立沖擊地壓多源多參量動態評級體系，對擾動區域局部近場危險性進行評價。利用貝葉斯組合概率方法，實現了沖擊地壓動靜態協同的遠近場聯合監測預警。通過嵌入評價模型的三維預警平臺進行測試對比，在一年的防沖日報表里，具有沖擊危險性的25 d中有22 d評價結果相同，另外3 d的評價結果偏高，其危險性評價模型的準確率為88%，但仍需進一步擴大示范應用地點和模型的應用時間，以提高模型適用性和精度。

(2)以微震事件定量預測為切入點，量化沖擊地壓的危險性區域與時間。研究構建了一個普通卷積模塊、循環神經網絡模塊與AR自回歸模型相結合的深度學習網絡MSNet，利用深度學習技術首次實現了未來微震事件的時間、地點、能量預測。預測結果表明，MSNet模型對微震事件位置的預測最為準確，,,坐標的預測余弦相似度都超過了0.93，對事件發生時間的預測結果稍差，平均余弦相似度為0.523。由于訓練樣本中能量屬性的不平衡，MSNet對事件能量的預測具有明顯的小能量偏向性。

(3)開發了沖擊地壓綜合智能預警3D可視化平臺，實現了現場監測數據的實時采集、同步傳輸、智能預警、動態展示，直觀地展示了沖擊危險性等級，危險性區域及潛在的災害時間。

該研究的核心是構建貝葉斯動靜態協同評價模型進行沖擊地壓危險性等級預測，同時利用MSNet模型預測未來微震事件，并以此為依據判定的沖擊地壓危險區域及時間。MSNet的預測精度直接決定了沖擊地壓量化預測的準確性。MSNet雖然在位置預測方面精度很高，但在能量及時間預測上仍有提升空間。以下列出的幾個方面可以成為進一步研究的重點：

(1)調整模型的損失函數結構。該研究計算MSNet模型損失時(式(26))，把獨立微震事件作為一個長度為5的向量，并計算真實值與預測值的向量距離。這其中隱含了一個假設，即向量的5個維度(、、、能量、時移)是等權重的。如果工程現場更關注預測能量準確性，可以將向量維度中的時移權重調小，從而讓模型損失的大部分由其他4個屬性的損失構成，保證能量預測精度。

(2)考慮微震事件的全波形信息。微震事件的全波形是壓電傳感器電壓隨時間變化的曲線，其包含了震源位置、介質屬性等，信息量豐富。而本研究中對波形特征只考慮了微震事件能量(通過波形的時間積分求得)與時間(通過AraWin軟件拾取到時)，忽略了其他大量波形信息。未來研究中，可以把微震事件全波形(獨立的超高維向量)作為模型的一個輸入屬性，這將大幅提升模型信息量，有望提升模型能量與時間的預測精度。

(3)沖擊危險性區域與時間確定。本研究中確定危險性區域的方式并沒有考慮6個微震事件在能量與時間上的梯度變化，因此確定的沖擊區域并無危險性高低的進一步精準劃分。另一方面，該研究給出的沖擊地壓時間的預測方法只考慮了歷史災害數據中能量最大的微震事件作為觸發時間預警的閾值。未來研究中可以更進一步建立歷史沖擊地壓發生前的微震事件的能量序列分布，從數據分布中提取特征值或者特征向量作為判斷閾值，增強預測的魯棒性。

動力災害的智能監測預警是煤礦智能化的重要組成部分。筆者利用物理驅動與數據驅動相結合的方法解決沖擊地壓量化預測的難題，既考慮了沖擊地壓發生的物理機制，又充分利用了海量的災害前兆數據。研究提出的“雙驅動”智能預警架構，給出了沖擊地壓的危險性等級、區域及發生時間的預測方法。同時，研究給出了“雙驅動”預警架構未來的研究方向，為徹底解決沖擊地壓的量化預測難題奠定基礎。