不同級配矸石壓實特性數值模擬研究

程立朝，秦義嶺，李新旺，趙新元

(1.河北工程大學礦業與測繪工程學院，河北 邯鄲 056038；2.河北工程大學河北省高校煤炭資源開發與建設應用技術研發中心，河北 邯鄲 056038)

綜合機械化固體充填采煤技術是一種綠色開采技術[1-2]，具有方便、快捷、廉價的優點，其中矸石是一種極好的固體充填材料。矸石作為充填材料，自身的密實度是影響充填效果的主要參數[3-4]，影響矸石密實度的主要因素是矸石顆粒的粒徑和配比。針對不同粒徑和級配下矸石材料的壓實特性，涂強等[5]研究了矸石散體不同級配與壓縮量的關系，得到了壓縮過程中的軸向應力-應變曲線；王文等[6]研究了粒徑級配、軸向應力對矸石散體壓實度的影響；胡炳南等[7]采用大容器、大粒徑、大載荷和數值模擬方法進行了矸石壓縮試驗；姜振泉等[8]重點探討了煤矸石粒度分布的級配與壓密性之間的關系；其他專家學者[9-13]也從不同方面對矸石的壓實特性進行了研究，但實驗室對于矸石的壓實特性研究具有工作量大、材料耗費多、實驗設計復雜的缺點。

本文在前人的基礎上使用EDEM離散元軟件對矸石的壓實特性進行了研究，尋找適合矸石充填最佳粒徑級配的有效方法，并在實驗室重復本實驗，驗證EDEM離散元軟件的有效性，為現場矸石充填的應用提供參考。

1 離散元法的原理和散體級配模型

1.1 離散元法基本原理

離散元法是將一個復雜的稠密顆粒集劃分為相互獨立的顆粒單元，對每個顆粒單元之間的相互作用使用牛頓運動學定律進行單獨研究，并通過循環迭代的方法更新和確定每一個顆粒的位置來觀測整個顆粒集運動狀態的方法。 顆粒的運動方程見式(1)。

(1)

采用中心差分法對式(1)進行積分，以迭代時間步長的中點表示計算更新速度，以此計算下一時間步長的位移，見式(2)。

(2)

式中：N為對應的時間；Δt為時間步長。

將式(2)的計算結果帶入到接觸模型中可求得所需要的接觸力，再經過數次迭代計算，則可以分析出顆粒在某個時間段的運動狀態。

1.2 散體級配與泰波公式

富勒(Fuller)理論認為，固體顆粒可以按照粒徑的大小，不同的粒徑顆粒間相互搭配可以得到密度大，孔隙率低的顆粒集合體，顆粒級配曲線越接近拋物線，則其密度越大，孔隙率越小。泰波(Tabol)在富勒理論的基礎上進行了修正，提出了泰波公式，見式(3)。

(3)

式中：P為矸石各級粒徑的通過率，%；d為矸石各級粒徑，mm；D為矸石的最大粒徑，mm；n為泰波公式系數，是表征矸石散體中粒徑分布的參數，其大小與矸石散體的孔隙率、密實度和顆粒的排列形式有關，而與矸石顆粒自身的特征無關，此次試驗中取0.3～0.7。

根據泰波公式，n取值0.3～0.7時可得到連續級配矸石配比表(表1)。由表1可知，n值的增大會導致矸石顆粒級配中粒徑小的顆粒占比相對降低。

表1 連續級配矸石配比表Table 1 Continuous grade gangue proportion table

2 矸石散體壓縮物理模擬試驗矸石散體壓縮數值模擬試驗

2.1 試驗裝置

試驗采用的加載設備為長春新科YA-600型微機控制全自動壓力試驗機，試驗系統可以自行采集試件壓縮過程中的應力應變等參數，準確地反映試件的壓縮破壞全過程。由于本次試驗中矸石材料最大粒徑為40 mm，根據能源行業標準《固體充填材料壓實特性測試方法》(NB/T 51019—2014)要求：試驗設計剛性筒內徑不小于固體充填物料最大直徑的3倍，因此本試驗矸石壓實模具形狀設計為圓柱形鋼筒，內徑120 mm，壁厚10 mm，高度200 mm。試驗裝置實物圖如圖1所示。

圖1 實驗裝置圖Fig.1 The experimental platform

2.2 試驗材料選取與試驗結果

選取邢邯地區某礦的原生矸石作為試驗的材料，對完整的塊狀矸石取巖芯制備直徑50 mm、高100 mm的標準試件，使用三軸液壓伺服試驗機對標準試件進行加載，獲得矸石材料的本征參數。對初步破碎的矸石使用網孔直徑為10 mm、20 mm、30 mm、40 mm的篩網進行篩分，按照泰波系數的指導來配置五組質量為2 kg的矸石材料。對五組試驗材料分別進行加載試驗，加載方式采用位移加載方式，加載速率為0.001 m/s，數據每秒記錄一次。繪制連續級配矸石壓縮過程的應力應變曲線圖，如圖2所示。

圖2 物理實驗連續級配矸石應力應變曲線圖Fig.2 Stress-strain curves of physical experiment of continuously-graded gangue

3 矸石散體壓縮數值模擬試驗

3.1 顆粒創建與顆粒接觸模型

利用傅里葉級數分析方法對矸石顆粒形狀的幾何參數(長度、寬度和高度等)進行定量分析，最終通過形狀指數和形狀系數對矸石形狀進行分類，得出矸石顆粒主要有錐形、棱柱形、橢球形、片形和圓柱形五種形狀，再對各粒徑中不同形狀的矸石顆粒數量進行統計，提高數值模擬試驗的準確性。采用矸石顆粒的典型形狀圖形在ProE軟件中建立三維模型，再將其三維模型文件在EDEM中加載并使用多球體組合的方式創建典型矸石顆粒形狀模型。典型矸石實物形狀與對應顆粒形狀模型如圖3所示。

圖3 典型顆粒形狀模型Fig.3 Typical shape model of gangue

由于矸石顆粒間黏聚力基本為0，因此本次數值模擬試驗顆粒接觸模型采用Hertz-Mindlin(No Slip)模型，模型基于Hertzian和Middlin的理論：顆粒單元法向力和切向力都具有阻尼分量，阻尼系數和碰撞恢復系數具有一定的關系，切向摩擦力遵循庫倫摩擦定律，滾動摩擦力通過接觸獨立定向恒轉矩模型實現。

法向力滿足式(4)。

(4)

式中：E*為當量楊氏模量；R*為當量半徑；δn為法向重疊量。

切向力滿足式(5)。

(5)

式中：St為切向剛度；G*為當量剪切模量；δt為切向重疊量。

滾動摩擦力滿足式(6)。

τi=-μrFnRiωi

(6)

式中：μr為滾動摩擦系數；Ri為接觸點到質心的距離；ωi為接觸點處的單位角速度矢量。

3.2 試驗材料參數獲取

1) 材料本征參數。本次試驗中使用的鋼材為Q235鋼，試驗中的矸石本征參數通過三軸液壓伺服試驗機對巖石標準試件加載試驗獲得，矸石和鋼材的本征參數見表2。

表2 試驗材料本征參數表Table 2 The intrinsic parameters of experiment materials

2) 接觸參數獲取。使用離散元法進行數值模擬試驗時，接觸參數的準確性直接決定了數值模擬試驗結果的準確性[14-15]。本實驗通過使用圓筒式法測量矸石自然安息角獲得接觸參數[16]，當室內物理試驗的安息角和數值模擬的安息角相等時，接觸參數可用于數值模擬試驗[17-18]。使用直徑120 mm、高200 mm的鋼桶，隨機在煤矸石堆中鏟取矸石物料，倒入鋼桶中。將鋼筒平穩放于鋪有薄紙的平臺上，緩慢豎向抬起鋼桶，使得煤矸石物料從鋼桶的底部緩慢放出，待煤矸石散落穩定后開始測定矸石物料的自然安息角，試驗過程如圖4所示。

圖4 實驗室安息角的測定Fig.4 Experimental of the rest angle

重復做3組試驗，每次試驗隨機選取三個不同方向測量矸石安息角，試驗數據見表3。最終得出本次試驗的矸石自然安息角角度為34.67°。

表3 實驗室安息角測定結果Table 3 Experimental determination results of the rest angle

利用EDEM中自帶的圖形編輯器創建鋼桶模型，鋼桶模型的參數與實際鋼桶尺寸相同，為直徑120 mm、高200 mm。在鋼桶模型的上方創建顆粒工廠用于生成矸石顆粒。生成的矸石顆粒級配與自然狀態下的矸石散體粒徑級配一致，各粒徑中不同形狀矸石顆粒數量分布遵循統計結果，仿真中矸石顆粒全部為隨機生成，可以視作最終生成的顆粒已經攪拌均勻。鋼桶內生成矸石顆粒完成之后，靜置0.5 s，使矸石顆粒間由于墜落產生的碰撞擠壓停止。在0.5 s后，使鋼桶以0.01 m/s的速度勻速豎直向上運行，使矸石顆粒從圓筒底部瀉出形成安息角。仿真結果如圖5所示。使用EDEM后處理中的Protractor功能對仿真的矸石安息角進行測量，最終實驗結果見表4。

表4 仿真安息角測定結果Table 4 Simulation experimental determination results of the rest angle

圖5 圓筒式仿真模型Fig.5 Simulation model of cylinder

實驗室試驗結果與仿真結果的相對誤差見式(7)。

(7)

式中：ER為相對誤差；Ea為絕對誤差，是分析結果與數次平行測定結果的算術平均值之差；T為真值。

經計算求得相對誤差ER為1.1%，遠小于5%的誤差允許范圍，故數值仿真中矸石的接觸參數符合要求，接觸參數見表5。

表5 試驗材料的接觸參數表Table 5 The contact parameters of experiment materials

3.3 數值模擬結果與分析

按照實驗室物理試驗裝置的尺寸和各部位間的位置關系建立的本次模擬的試驗裝置模型，簡化的試驗裝置模型由直徑120 mm、高200 mm的無蓋鋼桶和直徑120 mm的壓實板組成。試驗中矸石顆粒模型與上文參數校準試驗中的顆粒形狀完全相同，直接采用上述的顆粒模型繼續運算。試驗裝置與試驗材料的本征參數見表2，接觸參數采用圓筒法測量安息角校準后的參數見表5。按照泰波系數n取0.3～0.7時的粒徑級配占比，各粒徑中不同形狀的矸石顆粒數量分布嚴格遵循統計結果，在鋼桶模型中生成2 kg不同級配的矸石散體。然后使用壓實板對不同級配的矸石散體進行壓縮試驗，壓實板加載方式與實驗室物理試驗一致，采用位移加載方式，加載速率為0.001 m/s。試驗結果如圖6所示。

圖6 數值模擬應力應變曲線圖Fig.6 Stress-strain curves of numerical simulation of continuously-graded gangue

對圖6數值模擬應力應變曲線進行分析可知。

1) 應力位于0～1.8 MPa階段，泰波系數n取0.3和0.4時，兩者的應變基本未發生變化，而泰波系數較大的試驗中應變迅速增加，是由于泰波系數較小時，粒徑為20～30 mm、30～40 mm的矸石顆粒較少，自然混合情況下矸石顆粒間結合緊密、孔隙較少，說明大顆粒少的情況下連續級配矸石散體密實度較高，需要更大的應力才能使矸石散體進一步密實。

2) 當應力為10 MPa，n取0.3和0.4時矸石散體的應變基本相等，約等于0.14，但隨著n取值的增大，應變逐步降低至0.19，可以看出矸石散體在同一應力條件下應變隨著n值的增大呈現先降低、后增高的變化規律，拐點為n取0.4時。說明粒徑大的矸石顆粒在矸石散體中起到了一定支撐作用，但是隨著大粒徑矸石顆粒比例的提高，矸石顆粒間的孔隙率也逐步增大，嚴重影響矸石散體的密實度。

3)n取0.4時的連續級配矸石散體應力應變曲線斜率變化最小，說明該條件下的連續級配矸石散體壓縮率最小，密實度最高，可作為現場充填級配參考。連續級配矸石壓縮應力應變曲線由快速增加到緩慢增加的拐點應力均約為3 MPa，此后隨著應力的增大應變增加幅度甚微。

4) 由圖2與圖6對比可知，數值模擬試驗中泰波系數n取0.4時試件的最終壓縮形變最小，最小應變為0.14，n取0.7時試件的最終壓縮形變最大，最大應變為0.19；實驗室物理試驗中n取0.4時試件的最終壓縮形變最小，最小應變為0.23，n取0.6時試件的最終壓縮形變最大，最大應變為0.27。可見，物理試驗與數值模擬試驗相同條件下最終應變存在0.1的誤差，主要原因在于：實驗室中，矸石在壓縮過程中，當達到一定的壓力時，矸石顆粒間已經達到極限平衡，進一步壓實，將引起大顆粒矸石發生破碎，導致最終變形較大。總體上仿真結果基本符合實驗室壓縮的趨勢，且泰波系數取0.4的矸石顆粒級配是最優級配。

4 結 論

1) 使用圓筒式安息角測量法，對矸石物料進行了安息角實驗測定，驗證了EDEM對應參數下數值模擬結果的準確性，表明測得的矸石本征參數與接觸參數是可靠的，這為以后物料參數的標定提供了一定的理論基礎。

2) 實驗室連續級配矸石壓縮應力應變曲線和仿真級配矸石應力應變曲線形態接近，應變隨應力增大呈先快速增加后逐漸變緩的趨勢，拐點約為3 MPa，泰波系數n=0.4的連續級配矸石，在所有泰波系數級配壓縮試驗中密實度最大，可以為現場試驗作為參考。

3) 實驗證明，EDEM離散元仿真是一種有效的驗證和獲取矸石最優配比、密實度等宏觀參數的方法，可用于實際生產快速尋找合適的矸石級配。