環(huán)境濕度影響下粉晶白云巖的力學特征與數值模擬

陳 偉,萬 文,王衛(wèi)軍,馮 濤,趙延林,3,吳秋紅,3,周 彧

(1. 湖南工程學院 建筑工程學院,湖南 湘潭 411104;2. 湖南科技大學 資源環(huán)境與安全工程學院,湖南 湘潭 411201;3. 湖南科技大學 南方煤礦瓦斯與頂板災害預防控制安全生產重點實驗室,湖南 湘潭 411201)

在非金屬礦山的開采中,為了高效開采和節(jié)約成本,往往采用房柱法開采[1-2]。開挖時所留設的保安礦柱中通常含有伊利石,高嶺石,蒙脫石,綠泥石,云母等親水性礦物,這些礦物成分易水解發(fā)生軟化[3-5]。一直以來,許多學者在研究水化學溶液對巖石力學性能的改造作用方面取得了巨大成就。ZHAO等[6]發(fā)現,節(jié)理巖體在水-巖作用下其強度與剛度的劣化會嚴重影響庫岸邊坡穩(wěn)定性。丁梧秀等[7]利用不同的水化學溶液對灰?guī)r進行侵蝕,發(fā)現灰?guī)r產生損傷劣化程度受控于溶液pH與凝結核豐度。MIAO等[8]研究了花崗巖在各pH以及流速條件的化學溶液中浸泡后,其單軸、三軸的變形規(guī)律和強度損傷,并借助電子能譜和掃描電鏡從微細觀角度總結了其化學動力反應原理。韓鐵林等[9]對砂巖進行了水化學溶液浸泡并開展了一系列力學試驗,發(fā)現砂巖的劣化參數與強度呈線性遞減關系。JIANG等[10]探討了單軸壓縮下砂巖經AMD溶液腐蝕后的損傷本構關系,并構建相關模型,該模型能較好的呈現砂巖在AMD蝕化下的劣化規(guī)律。ATKINSON等[11]及DUNNING等[12]發(fā)現了無機水化學溶液對石英的斷裂力學性能的影響規(guī)律,并用Griffith強度理論對巖石受水化學環(huán)境的作用機制作出理論分析。

上述的諸多研究主要是針對各種水化學溶液以浸泡方式對巖石力學特性產生腐蝕作用而展開,并未涉及水化學流體以氣態(tài)方式侵入巖體。近年來,深部開采進行得如火如荼,由于地溫的升高,結合地下水發(fā)育,蒸發(fā)量加大,深部礦房的環(huán)境濕度平均值會超過80%RH[13]。且礦柱在長期的地質構造運動過程中其內部形成了斷層、解理、裂隙等缺陷,缺陷的存在給氣態(tài)水分子提供了滲透通道,為水汽侵蝕礦柱創(chuàng)造了有利條件[14-17]。因此,關于礦柱在高濕環(huán)境下的損傷力學特性亟待研究。

針對上述情況,筆者對湖南郴州茶山嶺煤礦-850 m房柱法開采段進行了實地踏勘,其水文地質報告顯示,該區(qū)常年地下水發(fā)育并呈弱酸性,相對濕度處于90%RH這一高位。為了對比不同環(huán)境濕度值對礦柱力學性能的削弱程度,筆者對現場礦柱取心制備標準試件并模擬3種高濕環(huán)境,通過MTS-815型試驗系統對高濕處理60 d后的試件展開單軸壓縮試驗,總結白云巖微細觀結構在各濕度環(huán)境下的變化特征,最后利用PFC離散元軟件的二次開發(fā)虛擬實現同等條件下的濕度場水汽侵蝕數值試驗,結合模擬試樣的微裂紋發(fā)育特征與能量演化規(guī)律,揭示在濕度擴散效應影響下白云巖能量演化的微細觀機制。

1 高濕侵蝕礦柱試樣原理

1.1 試驗系統組成

1.2 試樣制備

對茶山嶺煤礦-850 m的3號礦房5號礦柱沿沉積方向取心后加工成50 mm×100 mm的圓柱試樣[23],如圖2所示。試件制備結束后選取外觀較為良好且無瑕疵的試樣,以期獲得可靠的試驗結果。為了避免風化,將試件密封并以塑料薄膜封裝。利用煮沸法及AiniMR-60核磁共振分析儀分別測得試件平均含有3.42%的飽水率及4.93%的孔隙率。

本次采樣點礦房溫度處于25～30 ℃,且濕度對試樣的力學損傷是此次研究重點,故以恒定25 ℃為試驗溫度。表1為試驗所選用白云巖試樣的基本參數。

圖2 巖心采取與試樣制備Fig.2 Core collection and specimens preparation

表1 白云巖試樣基本參數

1.3 加濕與加載程序

以濕度控制箱對試樣進行加濕處理:先將酸性Na2SO4溶液按1 ℃/min速率加熱到25 ℃;再利用加濕系統分別設定80%RH、90%RH、100%RH試驗濕度,并持續(xù)工作60 d;最后,擦干試樣表面,并遵循控制變量原則[24],在DZF-2060型真空干燥箱中以120 ℃干燥24 h。

2 試驗結果分析

2.1 粉晶白云巖物理特性變化規(guī)律

2.1.1 不同濕度作用下白云巖孔隙率的變化規(guī)律

氣態(tài)水分子侵入試樣后,水解、溶蝕部分親水性礦物,并通過外滲將反應物運移,巖石的孔隙率因其內部產生的溶洞連通為滲透通道而增加,其孔隙壓力和滲透特性因此受到影響[21],試樣結構疏松度、比表面積和次生孔隙等增大。利用NMR核磁共振技術獲取白云巖在不同濕度環(huán)境下60 d后的孔徑分布,如圖3所示。T2為飽和樣品孔隙中含水原子核的橫向弛豫時間,利用外部靜磁場所測得。鑒于試樣橫向弛豫時間T2與其孔隙尺寸成正比[25],其表達式為

圖3 4種濕度環(huán)境下放置60 d后的白云巖T2譜規(guī)律Fig.3 T2 spectra of dolomites under four humidity conditions for 60 d

(1)

式中,ρ為飽和樣品T2表面弛豫強度;S為巖石孔隙表面積;V為流體體積。

基于孔隙的半徑與喉道大小的正比關系,式(1)可表達為

(2)

式中,r為孔隙半徑;Fs為樣品孔隙喉道的幾何形狀因子。

令1/(ρFs)=N(轉換系數),式(2)可化簡為

T2=Nr

(3)

天然巖石具備多孔介質屬性,因此目前還未能形成統一的標準對其進行孔徑劃分,文獻[26]以10 μm為界,將巖石中孔徑區(qū)分為大孔(r>10 μm)、小孔(r≤10 μm);針對低滲透巖石,文獻[27]劃分為大孔(1～5 μm)、中孔(0.1～1.0 μm)、小孔(r<0.1 μm)。筆者參考毛細管壓力測量孔隙半徑分級方法,結合試驗用樣的孔徑集中分布的范圍,將孔徑劃分為小孔(r≤ 0.1 μm)、中孔(0.1 μm

T2=10r

(4)

由圖3結合式(4)可將粉晶白云巖的T2譜中的橫向馳豫時間分布在0～1 ms劃分為小孔隙、1～10 ms為中孔隙、10 ms以上為大孔隙。試件內不同孔徑的譜面積分布及占比統計見表2。

表2 試件內不同孔徑、不同環(huán)境濕度的譜面積分布及占比

干燥試樣主要分布在0.1～10.0 ms弛豫時間內,且振幅遠低于經過高濕處理后的樣本,其T2譜總面積為3 255.26。對比來看,80%RH環(huán)境下試樣的T2譜總面積為160 631.61,增幅達49倍之多,且在濕度逐漸接近100%RH過程中,波峰不斷向右偏移,表明在高濕環(huán)境下與試件內部親水化合物發(fā)生劇烈的物理、化學反應,孔隙結構持續(xù)劣化。值得說明的是,100%RH環(huán)境下試樣譜線出現3峰,這是由于試樣內部黏土礦物吸水膨脹內部擠壓,導致一部分中、大孔隙體積減小轉化成小孔隙,小孔譜面積增大,即形成左側第一波峰。具體反應如下所述。

粉晶白云巖中主要碎屑礦物為正長石和云母,在酸性氣態(tài)水分子侵蝕下產生溶蝕孔洞,形成黏土礦物:

(5)

(6)

生成的黏土礦物高嶺石(Al2Si2O5(OH)4),其粒徑不大于2 μm,具備較大表面能,包括強親水性與膨脹性。高嶺石吸水膨脹化學式為

(7)

測定干燥狀態(tài)及3種高濕環(huán)境作用下60 d并烘干后白云巖試樣的總孔隙率(包括開口與閉合兩部分孔隙),計算公式為

n=(ρp-ρd)/ρp×100%

(8)

式中,n、ρd、ρp分別為試件的總孔隙率、干視密度、真密度。

從統計結果,以及白云巖表面的高清照片(50倍率,圖4)可知,巖樣礦物顆粒間的孔隙結構隨環(huán)境濕度的升高而愈發(fā)松散。

具體來看,白云巖孔隙率的增幅呈現明顯的階段性變化趨勢,當試樣所處環(huán)境從干燥狀態(tài)轉為80%RH時,其平均孔隙率增幅較大,由4.96%增至7.41%,增幅高達49.39%。100%RH高濕作用后的白云巖平均孔隙率較80%RH濕度作用后累計增幅達8.09%,符合負指數函數增長。

2.1.2 不同環(huán)境濕度作用下白云巖質量的變化規(guī)律

高濕作用帶來白云巖的孔隙率增加,因此有必要對各濕度環(huán)境作用前后的質量差異進行分析[28]。稱得3種經高濕環(huán)境處理60 d并沖洗外表松散脫落物后烘干的白云巖試件的質量,計算出質量差。根據前后質量差,可求出試樣的相對質量變化率,即質量損失因子,定義為D,即

圖4 4種濕度環(huán)境下放置60 d后的白云巖試樣孔隙率Fig.4 Porosity of dolomite specimens under 4 humidity conditions for 60 d

(9)

其中,M0、M(h)、ΔM(h)分別為巖樣的初始質量、經各濕度環(huán)境作用60 d并烘干后的質量、經各濕度環(huán)境作用60 d前后質量差。不同高濕作用前后各巖樣的質量變化與質量損失因子分布如圖5所示。

圖5 3種高濕環(huán)境下放置60 d的白云巖試樣的質量變化對比Fig.5 Comparison of mass changes of dolomite specimens under 3 humidity conditions for 60 d

從圖5(a)可以看出,白云巖試件進行各濕度作用后其質量皆有不同的減少,濕度越大質量損失越嚴重。試樣的平均質量損失在100%RH環(huán)境下為5.23 g,是80%RH環(huán)境下的5.03倍及90%RH環(huán)境下的2.41倍。同時,觀察圖5(b)、(c)發(fā)現,環(huán)境濕度越大,試樣的平均質量損失因子越大,且增幅不斷上升,符合冪函數增長趨勢。試樣在100%RH濕度下其平均質量損失因子是0.95%,遠高于80%、90%RH濕度下的0.19%和0.39%。

2.2 粉晶白云巖力學特性劣化規(guī)律

本次試驗中20組試件的單軸壓縮力學參數見表3。

表3 白云巖單軸壓縮力學參數

圖6展示了各濕度環(huán)境處理60 d后所有試樣的單軸壓縮應力-應變曲線。

不同濕度條件下,試件的單軸壓縮應力-應變曲線可分為(Ⅰ)～(Ⅳ)階段,分別為裂紋閉合壓密階段(Ⅰ),彈性變形至裂紋穩(wěn)定擴展階段(Ⅱ),裂紋非穩(wěn)定擴展階段(Ⅲ)與峰后擴容階段(IV)[4]。分析各階段的曲線特征以便進一步解讀白云巖礦柱(巖心)的漸進破壞過程。

(Ⅰ)為裂紋閉合壓密階段,在此階段,白云巖內部的微裂隙結構出現了閉合。此時應力-應變曲線呈下凹型,試樣的變形模量逐漸增大但側向膨脹較小。90%RH和100%RH濕度處理過的試樣甚至出現了微弱的側向收縮。主要是由于這2種環(huán)境下,空氣中的水分子密度較高,水巖接觸界面將會產生濕度梯度,迫使氣態(tài)水分子由高質量濃度區(qū)躍遷至低質量濃度區(qū)。濕度梯度越大,水分子對內部礦物組分中的親水性礦物溶蝕隨之加大,形成空隙。加載初期,孔隙被壓實。隨著載荷的增加試樣體積逐漸減少,此時軸向壓縮占主體地位。隨著環(huán)境濕度的增加,其裂紋的閉合應力降低,由35.97 MPa(干燥)降為13.56 MPa(100%RH),降幅達62.31%(圖7(a))。

試樣內部裂紋閉合完畢后進入彈性變形至裂紋穩(wěn)定擴展階段(Ⅱ)。白云巖在此階段的應力-應變曲線近似直線段。參數泊松比μ和彈性模量E則用來反映試樣這一階段的變形特性。結合式(10)、(11)求得平均彈性模量和泊松比[29],結果如表3和圖7(b)所示。

(10)

(11)

式中,參考試樣應力-應變曲線近似直線段(圖6),起點A和終點B的應力分別為σA和σB;軸向應變分別為εhA和εhB;徑向應變分別為εdA和εdB。

應力-應變曲線上直線段的E越小說明水分子對白云巖的侵蝕程度越深[22]。由圖7(b)可知,白云巖在3種濕度環(huán)境處理過后的彈性模量均低于干燥狀態(tài),隨濕度值呈負指數函數形式降低并趨于穩(wěn)定。相對濕度為80%、90%、100%環(huán)境下60 d后,試樣平均彈性模量分別為10.39、9.34和8.41 GPa,與干燥試樣相比分別降低了54.08%、58.76%和62.88%。由此可見,化學作用引起的腐蝕軟化效應會隨著環(huán)境濕度的上升而增強。此外,白云巖的變形參數泊松比與環(huán)境濕度值呈正相關。與干燥狀態(tài)相比,在100%RH條件下的白云巖泊松比由0.314增至0.826,累積增加量為163.41%。

當軸向應力達到了損傷應力時,表明試樣進入了裂紋非穩(wěn)定擴展階段(Ⅲ)。彈性變形至裂紋穩(wěn)定擴展階段(Ⅱ)和裂紋非穩(wěn)定擴展階段(Ⅲ)以擴容起始應力點為界(損傷應力),即體積壓縮與體積膨脹的拐點。各濕度環(huán)境下試樣的損傷應力在31.25～63.22 MPa,為峰值應力的51.86%～58.49%。該階段的軸向應力-應變曲線改變直線狀轉而呈現上凸形態(tài),即進入塑性變形階段。干燥、80%RH與90%RH環(huán)境下,軸向應力-應變曲線在峰值處產生應力平臺并波動(圖6(b)、(c))。此階段,新生新裂紋持續(xù)不斷的擴展致使試件的側向應變呈非線性增加。試樣內部微裂紋相互交匯連通導致試件出現宏觀破壞。在每種濕度下放置60 d后,干燥狀態(tài)和100%RH環(huán)境下試樣的平均峰值應力分別為121.91和53.44 MPa,明顯降低了41.05%。相反,隨著環(huán)境濕度值的增加,其試樣峰值點的軸向和側向應變不斷增加:干燥狀態(tài)和100%RH環(huán)境下試樣的平均軸向應變分別為0.72×10-2和0.96×10-2,增幅達32.78%;此外,這2種狀態(tài)下的平均側向應變分別為0.43×10-2和1.03×10-2,增加了137.97%。

圖6 不同相對濕度作用60 d后白云巖的單軸壓縮應力-應變曲線Fig.6 Uniaxial compressive stress-strain curves of dolomite specimens under different humidity actions for 60 d

這里,分別定義試樣峰值點應力的劣化度Q(h)和峰值點軸、側向應變的增長率P(h)、R(h),用來定量分析環(huán)境濕度對白云巖強度變形的影響規(guī)律,分別為

(12)

(13)

(14)

式中,σc(h)為巖樣經不同高濕環(huán)境處理60 d后峰值強度;ε1c(h)、ε2c(h)分別為其峰值點軸向和側向應變;σc0、ε1c0和ε2c0為干燥環(huán)境巖樣峰值點強度、軸向應變和側向應變。

圖7(c)展示了60 d的處理時間過后,白云巖單軸壓縮下峰值點的應力、軸向和側向應變隨環(huán)境濕度值而發(fā)生的典型階段變化率。白云巖試樣在峰值點的強度的劣化度和軸向、側向應變的增長率均隨著環(huán)境濕度的升高而增加。與干燥狀態(tài)相比,同樣的濕度作用時間60 d后,經80%RH、90%RH和100%RH處理后的白云巖峰值強度平均劣化度隨濕度呈負指數函數上升,分別為42.37%、54.14%和60.28%,分布曲線的增幅具有顯著非均勻性,即由干燥環(huán)境進入80%RH環(huán)境后,兩者的劣化度增幅較大,由80%RH環(huán)境進入100%RH環(huán)境導致的兩者劣化度的增幅相對較小。80%RH、90%RH和100%RH處理后的白云巖較干燥試樣峰值點軸向、側向應變對應的平均增長率隨濕度值線性增加,軸向應變增長率分別為18.76%、21.52%和33.01%,側向應變增長率分別為45.58%、92.79%和140.32%。

圖7 不同相對環(huán)境濕度下60 d后白云巖的單軸壓縮力學參數Fig.7 Uniaxial compression mechanical parameters of dolomite under different humidity conditions for 60 d

應力達到峰值強度后,試樣進入峰后擴容階段(IV)。此時,軸向應力-應變曲線分別呈斷崖式跌落(圖6(a))和峰值后繼續(xù)發(fā)生部分殘余變形這2種狀態(tài)(圖6(b)～(d))。第1種下跌狀態(tài)是由于干燥試樣的結構具有致密性,呈現出較強的脆性,應力出現驟降現象的同時往往伴隨脆裂聲,且產生短時間的峰值軟化。另一種下跌中試樣峰后殘余變形會隨著濕度的增加而增加。這一階段,裂紋發(fā)展速度加快,側向應變增長加劇,且發(fā)生多次應力降,每次應力降都會導致側向應變急劇增加而呈延性破壞。因此,白云巖隨環(huán)境濕度的升高會由脆性破壞向延性破壞過渡,且加劇試樣峰后的累積殘余變形。

3 基于微細觀損傷的數值模型

3.1 SEM電鏡掃描分析

對各濕度環(huán)境處理后的白云巖側面進行SEM電鏡掃描拍攝(500倍率),選取部分典型圖片展示(圖8)。觀察可知,干燥的白云巖試樣在切割打磨后表面結構相對穩(wěn)定致密,基本無微裂隙發(fā)育,不同晶體之間存在層理面,礦物顆粒整齊排列,完全被膠結物覆蓋,形態(tài)呈現明晰片狀及塊狀結晶(圖8(a))。這說明在進入高濕環(huán)境前,白云巖試件的宏觀力學性能良好。經80%RH環(huán)境60 d后,原規(guī)則緊密的塊狀、片狀結晶結構部分膠結物在水化作用下溶解,展露出海綿狀的團絮形態(tài),層間隆起,形成陡峭不一的“山脊”(圖8(b))。90%RH環(huán)境下,外表面部分礦物顆粒由于水化侵蝕溶解了膠結物,削弱了顆粒間的粘結力,部分顆粒呈凸出狀。與此同時,試樣表面發(fā)育有蜂窩狀的次生孔隙(圖8(c))。100%RH環(huán)境下,表層礦物顆粒間的膠結物繼續(xù)水化分解,片狀晶體在溶蝕作用下產生孔洞。相較80%RH的試件,礦物顆粒間膠結物缺失嚴重且出現顆粒大面積掉落現象,孔隙通道交匯連通產生微裂隙,整體微細結構趨于疏松,規(guī)則排列的晶體顆粒消失殆盡,白云巖微細觀結構損傷愈發(fā)嚴重。

圖8 不同濕度環(huán)境下60 d白云巖表觀電鏡掃描Fig.8 Scanning electron microscopy of dolomite surfaces under different humidity conditions for 60 d

由于礦柱上、下端面與頂底板連接,并未與水蒸氣直接接觸。因此,為了得到濕度擴散引起的白云巖礦柱腐蝕深度變化規(guī)律,將試件上、下端面用薄膜覆蓋,僅側面用來接受氣態(tài)水分子的侵蝕作用。為了定量分析白云巖在高濕環(huán)境下的微細觀損傷規(guī)律,將微細觀損傷定義為電鏡觀測面中的微裂隙和微孔隙。對干燥和3種濕度環(huán)境中處理60 d后的白云巖試件的頂底面利用掃描電鏡進行拍攝。為了拍攝精度,已對儀器內置的載物臺進行了坐標網格劃分(即以載物臺中心為圓心,半徑為2.5 cm圓形)共計78張電鏡掃描圖,作為損傷原始面(圖9)。

圖9 載物臺坐標網格劃分Fig.9 Coordinate gridding of loading platform

載物臺的尺寸經后期廠家改裝為半徑5 cm的圓盤,可容納高度13 cm的試樣進入拍攝。對各試件頂、底面按圖9中的編號進行依次逐點掃描。根據前期嘗試的效果,以及采集清晰度的需求,最終選取1 000倍作為本次定點定倍的拍攝方式。在SEM圖中,相對于完整部位來說,損傷部位的灰度為0。表4中利用Coreldraw軟件對每張掃描圖進行二值化處理。根據軟件內置的柱狀圖工具得到最佳灰度分割閾值為97,從而得出原始面每個點的SEM的“損傷因子”O(jiān)[30]為

(15)

式中,以二值化圖片為參考,LD為數值為0的面積;L為總面積。

表4展示了各濕度環(huán)境下白云巖試樣上、下端面的掃描圖以及二值化圖片。

平均損傷因子選取試件上、下端面對應兩點的損傷因子平均值,圖10展示了試件上、下端面平均損傷因子的頻率分布?？紤]到干燥試件存在天然缺陷以及加工制備期間端面產生的初始損傷,因此其電鏡掃描圖的損傷因子不為0。觀察可知,每個試件上、下端面的平均損傷因子分布具有明顯差別,以圖10中的數據經三參數峰值函數擬合后來直觀對比試件置于不同濕度下的損傷程度,結果表明白云巖上、下端面的平均損傷因子頻率直方圖符合Gauss函數分布。式(16)為Gauss函數擬合公式:

表4 白云巖試樣上、下端面的掃描圖以及二值化圖片

(16)

式中,J為比例參數;w為形狀參數;xl為位置參數。

干燥試件損傷因子在0.02～0.04內占總拍攝點數40.08%,而100%RH環(huán)境下試件的損傷因子在0.04～0.06內占總拍攝點數29.21%,反映出高斯曲線的波峰會因環(huán)境濕度升高而右移(圖10)。可以看出,損傷因子分布逐漸由低損區(qū)往高損區(qū)過渡,環(huán)境濕度增加直接加劇了微觀損傷的形成。其內部會在長期的高濕環(huán)境下產生一定的細觀損傷累積,這種巖體的強度漸進劣化作用會影響深部巖石工程的穩(wěn)定性。

基于圖9中的原始面二值化圖片分布坐標網格,可得到各試件上、下原始面平均損傷因子分布(圖11)。由圖11可知,對比干燥試樣,經高濕侵蝕作用后,試件邊緣區(qū)多見高損傷區(qū),中心區(qū)偶見高損傷區(qū)。隨著濕度的增加,細觀損傷由邊緣側逐漸向中心發(fā)育,不斷形成的高損區(qū)以“損傷-聯結-侵入”往試樣中心延展,且濕度越高,損傷侵入的深度越深。

圖10 白云巖試樣上、下端面的平均損傷因子頻率直方圖與擬合曲線Fig.10 Average damage factor frequency histogram and fitting curves of dolomite specimen top and bottom end faces

3.2 高濕侵蝕深度

上述微觀表征顯示高濕環(huán)境下白云巖試樣產生了由外及內的侵蝕作用。對各濕度環(huán)境下處理60 d后的試樣沖洗表層脫落顆粒后,進行尺寸測定,圖12展示了各環(huán)境下5個試樣的平均直徑變化曲線。

圖12顯示,白云巖試樣經高濕作用60 d后由于表層受剝蝕脫落其直徑減小,然而隨著濕度值越高,下降速度隨之平緩。進入80%RH環(huán)境后,試件直徑呈最大降幅,為1.06 mm,縮減率達2.31%。100%RH對比80%RH環(huán)境下的試樣直徑累計減小量為0.15 mm,縮減率為0.31%。由圖12可知,80%RH環(huán)境中的試樣較干燥試樣而言,氣態(tài)水分子會對白云巖產生強烈侵蝕作用,致使其外層發(fā)生較為嚴重的剝蝕脫落,80%RH～90%RH侵蝕速率相對較緩。

3.3 濕度場數值模型的建立

根據濕度傳質理論可知,存在濕度差的2種物體一旦產生接觸面就會形成濕度梯度,迫使氣態(tài)水分子由高濃度組分移至低濃度組分,通過調節(jié)最終趨于平衡[31]。這個過程稱為濕度擴散,本次試驗中的白云巖與濕度環(huán)境通過擴散效應進行能量交換[32]。能量交換的主要方式是分子擴散和對流擴散:分子運動引起的能量交換被稱為分子擴散,主要存在于流體和固體之間;而對流運動引起的物質轉移被稱為對流擴散,主要存在于2種流體之間。2種傳遞方式相比,分子擴散的物質傳遞比較平穩(wěn)[33]。在深部礦房中,主要存在空氣和礦柱之間水分的遷移,以分子擴散進行濕度的傳遞?；趥髻|理論,空氣和白云巖礦柱之間的水分能量交換如式(17)所示,表示單位時間內交換的對流水量。交換對流水量的影響因素眾多,比如對流交換系數,能量交換時間和接觸面積等。其中,不穩(wěn)定因素是對流交換系數。

Hf=mC(vw-vf)

(17)

式中,Hf為礦房中的濕交換量,kg/(m2·s);mC為濕度交換系數,m/s;vw為礦柱水巖接觸面的飽和水汽質量濃度,kg/m3;vf為礦房中對流層的水汽質量濃度,kg/m3。

水汽流動本質上屬于多孔介質滲流,氣體在多孔介質附近發(fā)生滲流的阻力公式[34]為

圖11 白云巖試樣上、下端面的平均損傷因子分布Fig.11 Average damage factor distribution images of dolomite specimen top and bottom end faces

圖12 試樣的直徑和直徑變化量隨環(huán)境濕度變化曲線Fig.12 Diameter and diameter variation curves of specimen with environment humidity

(18)

式中,Δp為壓降,N/m2;l為礦房中滲流帶長度,m;я為礦房中流體黏度,N/(m2·s);ο為礦柱離散元顆粒平均直徑,m;g為礦柱孔隙率,%;ρ0為礦房水汽質量濃度,kg/m3;v為礦房風流流速,m/s。

利用PFC2D模擬水汽傳輸,動量方程中的源項分別表示為黏性、慣性損失,如

(19)

式中,等號右側首項、次項分別為黏性、慣性損失系數項;Ki為x、y、z方向上礦房中多孔介質在動量方程中的源項;Aij為黏性損失的預定義矩陣;Bij為慣性損失的預定義矩陣;vj為x、y、z方向上的風流流速的分速度。

礦房水汽流動方向的黏性阻力系數q和慣性阻力系數I2為

(20)

(21)

式中,Ap為顆粒粒徑。

在上述由外及內的細觀損傷基礎上,引入應力腐蝕理論建立了高濕作用下礦柱的顆粒離散元接觸模型。在恒載地應力下,水分子持續(xù)腐蝕礦柱內部結構,礦柱的強度變形隨著時間的推移不斷劣化。編寫了FISH函數將PBM的模型進行改造,結合細觀侵蝕深度對礦柱試樣模型進行二次開發(fā),從而實現了應力腐蝕模擬:基于PBM模型和斷裂力學中KI型張拉裂紋和KⅡ型剪切裂紋(圖13)的亞臨界擴展速率的概念定義損傷速率v(計算見式(22)),不斷縮減顆粒的黏結半徑,得到試樣時效損傷變形的腐蝕數值模型。

圖13 濕度擴散引發(fā)顆粒間接觸力損傷的應力腐蝕模型Fig.13 Stress corrosion model for moisture diffusion-induced interparticle contact force damage

(22)

其中,rb為PBM線性平行黏結本構模型的黏結半徑;N1、N2為巖石材料的相關性系數;σb0、σb1、σb2分別為PBM線性平行黏結本構模型的法向應力、法向抗拉強度與起裂應力。采用PFC建模前,首先需要調試試樣的微觀力學參數,即依據室內試驗的應力-應變曲線,標定模型中的宏微觀參數。根據應力-應變曲線特征,確定干燥試樣的起裂強度σci為26.51 MPa左右。表5為最終選用的細觀力學參數。

值得說明的是,本次建模未能考慮非飽和滲流過程對巖樣力學行為的影響。我們通過簡化模型,直接以試樣在非飽和吸水狀態(tài)下的水汽入侵深度引入非飽和滲透系數。將粉晶白云巖試樣顆粒的潤濕性直接定義成微細觀損傷的損傷程度,認為顆粒之間的微張拉力和微剪切力的削弱是引起粉晶白云巖試樣微細觀時效損傷演化的內在驅動力[35-36]。

表5 數值模型的微觀和宏觀性質

結合圖11的濕度擴散導致的白云巖損傷速率,并通過自編程FISH語言來實現兩側加載板的濕度傳遞程序,計算得出3種環(huán)境濕度下試樣60 d內濕度擴散損傷區(qū)分布云圖(圖14),此處的大氣壓力是以標準大氣壓為參考。

由圖14可知,在礦柱(巖心)上、下隅角處濕度擴散速度最小,損傷最小,越靠近礦柱中段位置,腐蝕深度越深,且損傷分布面呈“扇型”對稱分布。由于水汽流的不斷涌入,礦柱中段容易接觸的水分子面積較大,且滲透壓力隨著垂直方向開始不斷減小。相同的作用時間下,更高濕度下水汽侵蝕的深度較深。

圖14 不同濕度環(huán)境下60 d白云巖濕度擴散分布云圖Fig.14 Moisture dispersion distribution clouds of dolomite under different humidity conditions for 60 d

3.4 數值模擬結果分析

圖15對比了室內試驗和PSC模型所得到的單軸壓縮應力-應變曲線。由于模型所賦參數是基于試驗數據的平均值,而不是單個試樣,因此2者的曲線之間存在一些誤差。但仍可以注意到,模擬試樣的應力-應變曲線和破裂形態(tài)與試驗結果基本一致。

單軸壓縮下,試樣只受軸向應力,因此巖樣的側向變形往往沒有約束。在本試驗中,試樣的最終破壞主要包括2種類型:拉伸破壞和混合剪切/拉伸破壞(圖15)。

圖15 白云巖受壓數值模型與室內試驗結果對比:應力-應變曲線與破裂形態(tài)Fig.15 Comparison of numerical modeling with laboratory results of dolomite under compression:Stress-strain curves and failure modes

干燥狀態(tài)下的白云巖主要呈現拉伸破壞,并且破壞形態(tài)會隨著環(huán)境濕度的增加轉變?yōu)榛旌霞羟?拉伸破壞。天然巖樣破裂面較順直,起伏不明顯,用手觸摸摩擦感不強,且破裂面上有少量巖屑分布;100%RH處理后的巖樣破裂面凹凸不平、起伏明顯,用手觸摸摩擦感較強且有明顯濕潤感,破裂面上巖屑分布多、顆粒大。由外及內的結合水膜的潤滑作用以及孔隙水壓力使得顆粒間膠結強度減弱,降低了巖樣外層原生裂隙尖端的應力強度因子。在單軸壓縮下,試樣外層最先發(fā)生起裂并呈斜對角線連接和貫通,最終導致宏觀剪切裂隙數量增加[37],致使白云巖宏觀破壞由拉伸破壞發(fā)展為剪切/拉伸混合破壞。

此外,觀察圖15與未示出的試樣發(fā)現,巖樣破壞時其表面發(fā)育新裂紋的數量會隨濕度值的增大而增多。這是由于干燥試樣受載前在長期的地應力作用下,內部缺陷已經形成穩(wěn)定的力學狀態(tài),因此荷載作用下表面只萌生導致試樣形成最終破壞的宏觀裂紋。另一方面,水分子侵蝕改變了試樣的內部結構的穩(wěn)定性。因此,在相同的加載速率下,新生裂紋由于松散結構的重組充分擴展,試樣表面的宏觀裂紋急劇增多。試樣數值模擬結果也顯示兩側的裂紋數量隨濕度的升高而不斷增加,與室內試驗結果基本吻合,證明所構建的腐蝕離散元顆粒模型(PSC)是可以反映現場高濕環(huán)境對白云巖礦柱力學特性由外及內的侵蝕損傷。

4 討 論

4.1 微裂紋發(fā)育特征

對PSC模型進行加載時,一旦其顆粒黏結力低于內部接觸強度,黏結力鏈即發(fā)生斷裂,導致平行黏結接觸轉化為線性接觸。這標志著PSC模型內部開始萌生微裂紋。試樣微裂紋的數量與類型極大地影響了顆粒的黏結半徑與顆粒平均半徑的比值。當黏結力小于接觸面顆粒的抗拉強度時,會產生張拉微裂紋。反觀黏結剪切強度大于黏結力時,剪切微裂紋隨之發(fā)育。圖16比較了4種不同濕度狀態(tài)下模擬試樣的微裂紋發(fā)育數量特征。

圖16 模擬試樣微裂紋發(fā)育特征Fig.16 Microcrack developmental characteristics of samples

可見,當試樣發(fā)生破壞時,內部微裂紋總數量隨著環(huán)境濕度的升高而增加,但張拉微裂紋數的增幅遠低于剪切微裂紋的增幅。如干燥狀態(tài)下的拉伸微裂紋和剪切微裂紋數分別是124和98,以拉伸微裂紋為主,而100%RH環(huán)境下的兩者數量基本接近,分別為689和736,漲了4.5和6.5倍。試驗過程中白云巖微裂紋的發(fā)育特征與數值模擬結果基本類似,說明在濕度擴散效應下由外及內的漸進侵蝕作用劣化了試樣的強度。同時,濕度越高,白云巖試樣壓縮破壞后更為破碎,呈現明顯的剪切特性和較弱的張拉特性。根據模擬結果可以得出,不同濕度狀態(tài)下試樣的微裂紋集中萌生于塑性變形段并在荷載達到峰值強度后數量陡增。這也進一步說明了峰后宏觀裂縫斷面的貫通和摩擦滑移促使礦物顆粒間的黏結斷裂,導致裂紋數量劇增。通過試驗所得出的耗散能演化特征與PSC模擬所得出的微裂紋發(fā)展規(guī)律相一致,揭示了高濕環(huán)境對白云巖侵蝕后受載所產生裂紋的能量機制。

4.2 能量演化的微細觀機制

單軸壓縮下PSC離散元模型的數值模擬,所有顆粒之間的能量關系為

Eboundary=Estrain+Eslip+Epbstrain+Ekinetic

(23)

式中,Eboundary為顆粒邊界能;Ekinetic為顆粒動能;Eslip為顆粒間摩擦滑移能;Estrain為顆粒間應變能;Epbstrain為顆粒間黏結能。

PSC模型顆粒間的能量形式主要有2種,分別是動能和接觸能。顆粒動能主要是模型受外力作用所吸收的邊界能,類似于總吸收能U。接觸能主要是儲存線性彈簧中的應變能以及平行黏結彈簧中的黏結能,類似于彈性應變能Ue;其次為耗散于微裂紋形成與剪切滑移中的摩擦能,與耗散能Ud類同。所以,受載時PSC模型的邊界能約等于應力達到峰值前的顆粒應變能、摩擦能和黏結能的總和。

為了方便對比,圖17給出了干燥狀態(tài)和100%RH環(huán)境下白云巖PSC模型在單軸壓縮下的能量演化監(jiān)測曲線?？芍?第1,加載過程中所做的功大部分轉化為顆粒的應變能。當試樣中的微裂紋開始出現時,應變能開始釋放并轉化為顆粒耗散能和局部阻尼能。第2,軸向荷載達到峰值以后,應變能釋放并轉化為摩擦斷裂能,表明斷裂后顆粒的摩擦作用直接導致峰后殘余強度的形成。圖17還顯示,在應力達到峰值之前,應變能和黏結能的總和在邊界能中的占比超過了90%。當總值超過能量累積的極限,儲存在模型內部的能量會隨試樣破壞而釋放,這個儲能極值受控于巖石材料的強度劣化度。此外,這與試驗中高濕作用降低白云巖的儲能極限相吻合。隨著微裂隙界面摩擦作用的加劇,摩擦能也在不斷地攀升。與此同時,在應力達到峰值時,宏觀裂隙發(fā)生貫通,導致摩擦能出現陡增的趨勢。這與室內試驗中摩擦耗散能隨巖樣塑性變形的發(fā)展而增加,峰值后其能量值急劇增加相一致。

圖17 PSC模型試樣的能量演化曲線Fig.17 Energy evolution curves of PSC model specimens

圖17中能量演化特征的計算結果表明:高濕環(huán)境下白云巖的邊界能、黏結能和應變能遠小于干燥狀態(tài)下試樣所對應的能量值。由此可知,環(huán)境濕度的升高降低了白云巖的總吸收能U與彈性應變能Ue,即高濕作用使得試樣的吸能性能減弱。對比干燥狀態(tài)和高濕環(huán)境下試樣的摩擦能,白云巖在高濕侵蝕后的能量耗散略有增加,說明濕度擴散效應導致內部裂隙發(fā)育,其微裂紋數量的增加提升了顆粒間摩擦滑移耗散能。

在單軸壓縮條件下,根據白云巖模型在干燥和3種濕度狀態(tài)下的模擬試驗數據,分別從能量演化角度進一步揭示白云巖受環(huán)境濕度影響的微細觀機制:① 白云巖在軸向應力的持續(xù)做功下不斷吸收能量,并以顆粒黏結能與應變能結合的形式儲存和積累,形成邊界能;高濕侵蝕迫使顆粒間的黏結發(fā)生斷裂,因此進入塑性階段時首先需要克服的是顆粒間的黏結能,然后微裂紋才能隨著應變能的增加而繼續(xù)擴展;摩擦能對裂紋滑移有著決定性的作用,其占邊界能的比例隨著裂紋擴展和滑移程度的加劇而增加。② 應力達到峰值后,儲存的黏結能和應變能迅速釋放;高濕處理后試件的微裂紋數量以及發(fā)育、貫通程度都在增加,使摩擦能在總能量中占據主導地位;然后,剪切微裂紋迅速擴展并貫通,成為宏觀裂紋的主要表現形式。

5 結 論

(1)茶山嶺煤礦粉晶白云巖的力學特性受控于環(huán)境濕度。對比高濕環(huán)境與干燥狀態(tài)下的試樣,80%～100%RH環(huán)境下試樣的承載力、抗變形能力、強度、脆性變形都在減弱,試樣受載破壞時表面的宏觀剪切裂紋數量增加。

(2)環(huán)境濕度的升高促使白云巖由密實片狀結晶向團絮狀態(tài)發(fā)展,框架結構趨于松散、層間分界逐漸消失、微細孔隙裂隙大量發(fā)育、礦物顆粒剝蝕脫落,最終致使微觀致密結構演化為多孔結構。高濕環(huán)境中白云巖試樣的損傷不斷以“損傷-聯結-侵入”由外沿往中心發(fā)展,且濕度越高,損傷侵入的深度越深。

(3)基于濕度擴散效應以及濕度傳質理論建立了濕度場二維數學模型。白云巖礦柱的濕度擴散速度受控于環(huán)境濕度,侵蝕損傷分布規(guī)律呈“扇形”對稱分布。在垂直方向中部侵蝕深度較深,上、下隅角處數值較小。

(4)對顆粒離散元軟件PFC2D中的PBM模型進行FISH語言的二次開發(fā),虛擬完成了濕度侵蝕PSC數值模擬。從微裂紋發(fā)育以及能量演化的角度,揭示了高濕環(huán)境下粉晶白云巖受載破壞時能量耗散的微細觀機制。