基于Weibull分布的非均質性對花崗巖聲發射與能量影響的離散元模擬

沈琴冰

(廣東粵源工程咨詢有限公司，廣州 510600)

1 概 述

巖石是天然產出的、由一種或多種礦物或其他物質構成的固態集合體。力學性質的非均質性是巖石基本的性質之一，巖石內部各種礦物成分均對巖石材料的宏觀力學性質有重要的影響[1]。巖石內部礦物組成成分的非均質性是造成巖石在外荷載下表現出宏觀力學性質離散性的原因之一[2]。

巖石的非均質性控制在巖石的宏觀力學響應，到目前為止，C.A.Tang等[3]對巖石的非均質性對巖石強度與變形特性的影響進行了大量的研究。結果表明，非均質性對于巖石強度的最大影響在于非均質性會引起加載過程中巖石內部的局部應力集中現象，從而使得巖石力學性質的弱化和破壞的快速進行[4]。由于在室內試驗中難以準確測量巖石的非均質性，因此目前研究巖石的非均質性的方法主要為數值模擬。宏觀力學參數的統計規律較為符合Weibull分布，因此許多學者將巖石的力學參數假定服從Weibull分布，以此來研究力學參數的非均質性對巖石強度的影響。如RFPA軟件中巖石單元的力學參數可以設置為服從Weibull分布，研究表明，非均質性對裂紋擴展和強度產生很大影響[3,5]。羅榮等[6-8]采用細胞元數值模擬方法研究了Weibull分布對于模擬多礦物組成的巖石的適用性，研究表明非均質巖石的力學特性由細胞單元力學特性和其非均質特性共同決定。

采用數字圖像處理技術，可以將巖石彩色照片識別出不同礦物組分的具體位置，從而可以基于真實巖石圖像建立非均質模型，實現細觀尺度的非均質巖石模型建模。如岳中琦等[9-10]、于慶磊等[11]以及李國瑩[12]等均通過對組成巖石模型的不同單元賦予不同的力學參數來定量模擬巖石的力學性質的非均質性。

花崗巖作為一種良好的建筑材料，具有優越的力學性能。花崗巖主要由斜長石、黑云母、石英等礦物成分組成。不同礦物成分之間的力學性質差別巨大，例如黑云母單軸壓縮強度很小，而石英單軸壓縮強度較高。因此，造成了花崗巖力學性質的非均質性。而PFC軟件[13]在研究花崗巖的力學性質方面表現出良好的性能，由于PFC軟件屬于微細觀建模，非常適于花崗巖的非均質性研究。本文通過北山花崗巖的室內試驗數據[14]進行對比，驗證本文所建模型的正確性。然后通過PFC2D中微觀參數服從Weibull分布，從而研究花崗巖力學參數的非均質性對巖石聲發射和能量的影響。

2 基于Weibull分布的非均質模型的建立

基于Weibull分布的力學參數賦值方法是數值模擬當中表征巖石材料非均質性時最廣泛的方法[6]。本文將PFC2D中的4個微觀參數(顆粒接觸模量、平行黏結模量、平行黏結抗拉強度和平行黏結黏聚力，具體數值見表1)服從Weibull分布，以此研究Weibull分布中形狀參數m對巖石聲發射和能量演化的影響。

接下來以平行黏結抗拉強度為150 MPa時，介紹PFC2D中如何實現微觀參數非均質賦值。Weibull分布的概率密度函數表達式和累積分布函數表達式為：

(1)

F(σt)=1-e-(x/σt0)mx>0

(2)

采用Monte-Carlo法生成隨機變量：

(3)

式中：R為(0,1]之間的均勻分布。

本文中，m取值分別為2、4、6、10。平行黏結抗拉強度為150 MPa時，平行黏結抗拉強度的概率密度函數圖見圖1。本文所建立的模型寬50 mm，高100 mm，高寬比為2，與室內實驗[14]保持一致。在PFC2D中顯示見圖2。

圖1 平行黏結抗拉強度概率密度圖

圖2 平行黏結抗拉強度數值分布 (圖中數值單位為Pa)

PFC2D中一般通過“試錯法”不斷調整微觀參數，使得所建立的非均質模型的宏觀力學性質和室內實驗中花崗巖宏觀力學性質基本一致。

首先通過均質模型標定微觀參數，即同一參數不同平行黏結接觸取同一值。最終選取的微觀參數見表1。宏觀力學性質對比見表2。表2證明本文所取微觀參數基本可以反映室內實驗中花崗巖的宏觀力學性質。

表1 PFC微觀參數表

表2 宏觀力學性質對比表

接下來，采用Weibull分布隨機賦值的方法，將形狀參數m設置為2、4、6、10，特征參數見表1，研究力學參數的非均質性對巖石聲發射和能量演化的影響規律。

3 結果分析

3.1 應力-應變曲線

本文所建立的5個模型(均質模型、m=2、m=4、m=6、m=10)的應力-應變曲線見圖3。宏觀力學性質(彈性模量、單軸壓縮強度和泊松比)見圖4。

圖3 應力-應變曲線

圖4 宏觀力學性質隨形狀參數變化圖

由圖3、圖4可知，力學參數的非均質性對巖石的宏觀力學響應影響非常大。隨著形狀參數m的增加，模型非均質性減弱，模型的力學性質越來越均一化，使得模型的單軸壓縮強度逐漸增加。當m=10時，單軸壓縮強度較均質模型減小4.6%，彈性模量減小2.5%，而泊松比增加2.8%。彈性模量和單軸壓縮強度隨形狀參數擬合成指數函數關系，泊松比隨形狀參數擬合成線性函數關系。

3.2 聲發射特征分析

在室內實驗中，花崗巖內部的礦物顆粒及其膠結物在加載過程中會發生破裂而釋放能量，產生聲發射事件。在PFC2D模擬中，通過顆粒間接觸斷裂來模擬聲發射事件[15-18]。圖5為單軸壓縮下每20個時步內聲發射數(接觸斷裂數)與應變關系曲線。

由圖5可知，5個模型聲發射特征基本類似，大致可以分為4個階段：初始壓密階段、彈性階段、峰前塑性階段和峰后破壞階段。在加載初期，聲發射數較少(圖5(b))或者基本沒有(圖5(a)、圖5(e))；隨著加載的進行，軸向應力-軸向應變曲線進入彈性階段，聲發射數逐漸增加，

圖5 聲發射、應力-應變關系

但是整體依然較少；到達峰前塑性變形階段，聲發射數突然增加，表明模型發生大量損傷，即將到達峰值強度點；峰值強度點至峰后破壞階段，聲發射數達到頂峰，表明模型形成大量宏觀可見的貫通裂紋，模型已經失去承載能力。

從圖5可知，非均質度對聲發射特征有一定的影響。其中最明顯的區別就是當形狀參數小(模型非均質度高，模型力學參數分布不均勻)的時候，初始壓密階段會隨機產生一定量的聲發射信號；而當形狀參數大的時候，初始壓密階段則基本沒有聲發射信號。原因是因為當模型力學參數分布不均勻的時候，少數接觸力學強度很小，在加載初期就會發生破裂。因此，可以認為相較于均質模型(圖5(a))，采用Weibull分布賦值的模型(圖5(b)-圖5(e))產生了一定初始損傷。

3.3 能量演化特征分析

在PFC2D中，可以測量出模型在加載過程中各種能量的變化[19]。總應變能為顆粒黏結能和應變能之和，研究證明，從PFC2D中細觀力學的角度研究能量演化特征極具優勢[20]。根據能量守恒定律，定義耗散能=邊界能-總應變能。不同種類的能量(邊界能、動能、耗散能、摩擦能、總應變能)隨形狀參數變化關系見圖6。

圖6 形狀參數為2時能量變化

由圖6可知，隨著加載的進行，外界輸入模型的邊界能逐漸增加，但是由于初期模型無裂紋產生，邊界能全部轉化成模型的應變能， 而耗散能、動能和摩擦能都基本為0。在彈性階段，微裂紋逐漸萌生，裂紋擴展所消耗的能量逐漸增加，到達峰前塑性階段時，耗散能急劇增加，表示模型微裂紋數量急劇增加，對應圖5中峰前塑性階段聲發射數猛增。此時，顆粒之間相互滑移，導致模型摩擦能和顆粒動能也迅速增加。軸向應力到達峰值強度時，總應變能也幾乎同步到達峰值，然后模型破壞失效，總應變能迅速釋放，耗散能持續增加。在整個模型加載過程中，顆粒動能占系統能量比例很小；總應變能隨軸向應變先不斷增加，后減小；耗散能不斷增加，且占系統能量比例不斷增加。

非均質性對系統能量演化的影響見圖7。

由圖7可知，力學參數的非均質性對巖石能量演化具有很大影響。隨著形狀參數m由2增加至10，達到同一軸向應變時，邊界能和總應變能逐漸增加，而顆粒動能、耗散能和摩擦能逐漸減小。當形狀參數小時(如m=2)，模型非均質度高，在加載過程中，模型更加容易產生應力集中現象，使得微裂紋的萌生與擴展更加容易，從而使得微裂紋表面能所耗散的能量較少，耗散能較其他模型(如模型m=10)較少，且模型儲存的總應變能也最小。因此，當模型中微裂紋擴展所需耗散能量小時，模型更加容易在加載過程中形成宏觀可見的裂紋，使得模型的單軸壓縮強度降低。

4 結 論

1) 力學參數的非均質性對巖石的宏觀力學響應影響非常大。形狀參數m越大，單軸壓縮強度越高。彈性模量和單軸壓縮強度隨形狀參數m擬合成指數函數關系，泊松比隨形狀參數m擬合成線性函數關系。

2) 聲發射大致可分為4個階段：初始壓密階段、彈性階段、峰前塑性階段和峰后破壞階段。非均質度對聲發射特征有一定的影響。當形狀參數m小時，模型力學參數分布不均勻，初始壓密階段會隨機產生一定量的聲發射信號；而當形狀參數大的時候，初始壓密階段則基本沒有聲發射信號。

3) 隨著形狀參數m增加，同一軸向應變時，邊界能和總應變能逐漸增加，顆粒動能、耗散能和摩擦能逐漸減小。