張慧梅，袁 超，慕娜娜，張 嬋，路亞妮

(1.西安科技大學(xué) 理學(xué)院，陜西 西安 710054；2.西安科技大學(xué) 建筑與土木工程學(xué)院，陜西 西安 710054；3.湖北工程學(xué)院 土木工程學(xué)院，湖北 孝感 432000)

0 引 言

“一帶一路”倡議的提出，為中國(guó)煤礦行業(yè)提供了新的發(fā)展契機(jī)，隨著開(kāi)采水平的提高，開(kāi)采深度不斷增加[1]。相應(yīng)地，季節(jié)性?xún)鋈诘挠绊懮疃群陀绊懛秶搽S之增加，對(duì)巷道工程施工，尤其是砂巖巷道的支護(hù)提出了更高的要求。因此，開(kāi)展凍融循環(huán)作用下煤礦砂巖的力學(xué)特性研究對(duì)指導(dǎo)西部寒區(qū)煤礦安全開(kāi)采至關(guān)重要。

巖石凍融破壞為寒區(qū)工程建設(shè)中的常遇災(zāi)害，針對(duì)此問(wèn)題，國(guó)內(nèi)外學(xué)者進(jìn)行較為全面的研究，如徐光苗等通過(guò)進(jìn)行力學(xué)特性試驗(yàn)，分析不同宏觀力學(xué)參數(shù)與凍融循環(huán)次數(shù)之間的關(guān)系[2]；張慧梅等對(duì)巖石在凍融循環(huán)作用下的損傷劣化規(guī)律進(jìn)行研究，提出凍融荷載耦合的損傷模型[3-4]。任曉龍等對(duì)四川芙蓉白膠煤礦的煤樣進(jìn)行不同溫度下的單軸壓縮試驗(yàn)，研究溫度對(duì)煤巖物理力學(xué)性質(zhì)的影響[5]。LUO等通過(guò)以輝綠巖為研究對(duì)象，進(jìn)行不同次數(shù)的凍融循環(huán)試驗(yàn)，得到其質(zhì)量和P波波速隨凍融作用的變化規(guī)律[6]。KHANLARI等以伊朗的典型紅砂巖為研究對(duì)象，進(jìn)行室內(nèi)以及現(xiàn)場(chǎng)的凍融循環(huán)試驗(yàn)，分析其物理力學(xué)性質(zhì)隨凍融次數(shù)的變化情況[7]。ZHANG等通過(guò)真實(shí)紅砂巖預(yù)制裂隙的方法，研究不同裂隙幾何特征對(duì)巖體強(qiáng)度及裂紋擴(kuò)展演化規(guī)律的影響[8]。

但僅從宏觀角度出發(fā)研究，難以解釋試樣在凍融過(guò)程中產(chǎn)生的層狀破裂和片落破壞[9]。為揭示凍融破壞的細(xì)觀機(jī)理，細(xì)觀尺度的研究受到研究人員的重視，如翟成等利用核磁共振，分析凍融作用下煤體孔隙結(jié)構(gòu)特征的變化規(guī)律[10]；李杰林等進(jìn)行核磁共振試驗(yàn)，研究了巖石在凍融作用下內(nèi)部孔隙結(jié)構(gòu)的損傷特性，分析核磁共振圖譜與凍融循環(huán)次數(shù)之間的關(guān)系[11]。聞名等基于SEM圖像，分析水分及凍融環(huán)境對(duì)巖石動(dòng)態(tài)抗壓強(qiáng)度的影響[12]。在眾多細(xì)觀試驗(yàn)技術(shù)中，計(jì)算機(jī)層析掃描技術(shù)(CT)，因其可視化程度高，使用最為廣泛，國(guó)內(nèi)楊更社較早地開(kāi)展相關(guān)試驗(yàn)，并提出了一種采用CT數(shù)描述的凍融損傷模型，實(shí)現(xiàn)對(duì)凍融損傷的定量描述[13]；劉慧等運(yùn)用CT圖像直方圖技術(shù)，定量分析巖石在凍結(jié)情況下的損傷信息[14]。國(guó)際上，ARGANDONA等利用CT檢測(cè)技術(shù)對(duì)凍融巖樣進(jìn)行斷面掃描，研究不同凍融條件下巖樣的孔隙結(jié)構(gòu)變化規(guī)律[15]；PARK和KOCK等則分別利用CT技術(shù)研究火成巖和沉積巖在凍融循環(huán)作用下孔隙擴(kuò)展過(guò)程[16-17]。

隨著巖石細(xì)觀試驗(yàn)技術(shù)的成熟和試驗(yàn)裝置的普及，合理的圖像后處理成為決定試驗(yàn)的準(zhǔn)確性關(guān)鍵因素，王超等提出一種基于LBP和GLCM的煤巖圖像特征提取方法，可以實(shí)現(xiàn)巖石表面幾何特征的快速識(shí)別[18]。針對(duì)CT圖像，張嘉凡等基于改進(jìn)的K均值聚類(lèi)算法對(duì)砂巖CT圖像進(jìn)行閾值分割，實(shí)現(xiàn)凍融巖石橫截面孔隙分布的有效識(shí)別[19]；張慧梅等則考慮孔隙結(jié)構(gòu)的縱向空間分布，實(shí)現(xiàn)凍融巖石孔隙分布的三維重構(gòu)[20]。

常用的圖像處理技術(shù)多數(shù)以傳統(tǒng)歐氏幾何為基礎(chǔ)，僅通過(guò)二維孔隙面積或者三維孔隙體積來(lái)描述巖石的凍融損傷程度，無(wú)法對(duì)巖石孔隙的結(jié)構(gòu)變化特征進(jìn)行描述，巖石分形理論的提出為砂巖孔隙特征定量描述提供新的思路[21]，史宏財(cái)對(duì)煤樣進(jìn)行各個(gè)溫度預(yù)損傷和三軸聲發(fā)射試驗(yàn)，探討高溫預(yù)損傷后煤巖的力學(xué)行為[22]。

張志鎮(zhèn)等利用壓汞法研究高溫處理后花崗巖樣品的孔隙特征和分形結(jié)構(gòu)特征[23]。王秀娟等利用分形理論，計(jì)算得到煤滲流孔的分形維數(shù)，探討分形維數(shù)與滲透率的關(guān)系[24]。通過(guò)引入分形維數(shù)來(lái)描述圖像復(fù)雜性特征的方法受到廣泛認(rèn)可，該方法也被用于巖石孔隙幾何特征的研究，但目前已有研究多數(shù)分析特定CT圖像的較多，分析砂巖整個(gè)縱向高度分形特征的較少；定性分析的較多，定量描述的較少。

以新鮮采集的煤礦砂巖為研究對(duì)象，通過(guò)CT掃描試驗(yàn)，獲取掃描試樣在不同凍融循環(huán)次數(shù)下的斷面掃描圖像，經(jīng)預(yù)處理后，得到不同截面分形維數(shù)隨凍融循環(huán)次數(shù)的變化規(guī)律。并通過(guò)分形維數(shù)定義一種描述巖石孔隙結(jié)構(gòu)損傷的新型損傷變量，定量地分析不同凍融循環(huán)次數(shù)下煤礦砂巖的損傷演化規(guī)律。

1 試驗(yàn)概況

受地理環(huán)境的影響，中國(guó)礦區(qū)分布的煤礦砂巖受凍融災(zāi)害影響較多，考慮到研究對(duì)象的典型性，筆者依托陜西彬縣和長(zhǎng)武縣內(nèi)的煤礦礦井建設(shè)項(xiàng)目，采集一批新鮮的煤礦砂巖樣品。試樣制備完成后，初篩剔除外表明顯缺陷的試件，然后利用聲波檢測(cè)儀篩選波速相近的試件進(jìn)行相關(guān)試驗(yàn)，以便于盡可能的排除自身差異性。

1.1 凍融循環(huán)試驗(yàn)

依照《工程巖體試驗(yàn)方法標(biāo)準(zhǔn)》(GB/T50266—2013)[25]，凍融試驗(yàn)前先將巖樣置于(105±1)℃的恒溫箱中烘干24 h，然后采用真空飽和儀對(duì)巖樣進(jìn)行0.1 MPa條件下的強(qiáng)制飽水。

凍融循環(huán)試驗(yàn)由圖1所示的XMT605快速凍融試驗(yàn)機(jī)完成，試驗(yàn)過(guò)程中試樣處于無(wú)約束狀態(tài)，按溫度控制凍融循環(huán)時(shí)間，巖心溫度由+20 ℃降低到-20 ℃后再重新升溫到+20 ℃為一個(gè)凍融循環(huán)。為保證試樣中心溫度達(dá)到設(shè)定值，試驗(yàn)中設(shè)置有控制試樣，內(nèi)置有高精度熱敏電阻溫度傳感器。

圖1 XMT605凍融循環(huán)試驗(yàn)機(jī)Fig.1 XMT605 freeze-thaw cycle testing machine

通過(guò)對(duì)試驗(yàn)現(xiàn)象的觀測(cè)發(fā)現(xiàn)，凍融次數(shù)較少時(shí)，砂巖試樣整體完好，僅表面有零散的顆粒脫落；凍融次數(shù)達(dá)到5次時(shí)，大部分砂巖試樣仍能保持完整，有少量試樣表面出現(xiàn)淺裂紋；當(dāng)凍融次數(shù)增加到10次時(shí)，試樣表面顆粒開(kāi)始大量脫落，脫落體多為顆粒狀，表現(xiàn)為剝落破壞；當(dāng)凍融次數(shù)達(dá)到20次時(shí)，試樣的破壞形式變化不大，部分試樣表面開(kāi)始出現(xiàn)片狀脫落；當(dāng)凍融次數(shù)到達(dá)40次時(shí)，巖樣的破壞主要仍集中于兩端，主要為嚴(yán)重的層裂破壞，小部分出現(xiàn)了貫穿的橫向裂紋，試件沿著裂紋面逐漸滑移成2個(gè)部分，典型的破壞形式如圖2所示。

圖2 砂巖試樣典型的破壞模式Fig.2 Typical failure modes of sandstone specimens

1.2 計(jì)算機(jī)層析掃描(CT)試驗(yàn)

CT試驗(yàn)試驗(yàn)裝置為圖3的YTU225型工業(yè)X射線CT機(jī)，掃描精度為0.005 mm，掃描照片分辨率為96 DPI。

圖3 YTU225型工業(yè)CT Fig.3 YTU225 industrial CT machine

實(shí)驗(yàn)掃描層數(shù)預(yù)設(shè)為900層，由于吸漲作用和端面脫落，實(shí)際得到不同凍融次數(shù)下的掃描圖像數(shù)量有輕微差異，凍融次數(shù)n=0，5，10，20，40時(shí)第450層的CT圖像，如圖4所示。

圖4 不同凍融次數(shù)下第450層的CT圖像Fig.4 CT image of layer 450 under different numbers of freeze-thaw cycle

CT圖像是與巖石密度相關(guān)的灰度圖像，它是由0～255個(gè)灰度值所構(gòu)成，灰度值直接反映了物體的密度，圖4中高亮點(diǎn)代表的是密度大的巖石顆粒，暗些的點(diǎn)代表的則是密度較小的孔隙結(jié)構(gòu)。

2 CT圖像預(yù)處理

CT圖像雖然能夠更加直觀地體現(xiàn)孔隙變化規(guī)律，但由于肉眼分辨能力較低，不能準(zhǔn)確察覺(jué)灰度圖像中的細(xì)微變化，且獲取圖像的過(guò)程中會(huì)受到來(lái)自各個(gè)方面的噪聲影響，給后續(xù)凍融巖石的損傷劣化定量分析帶來(lái)了困難，因此還需進(jìn)一步對(duì)圖像進(jìn)行預(yù)處理，減少噪聲，豐富CT圖像中的力學(xué)信息。

2.1 中值濾波

中值濾波就是把圖像中的每個(gè)點(diǎn)都用其領(lǐng)域窗口內(nèi)的中值代替，從而消除孤立的噪聲點(diǎn)，其好處是可以減少圖像中的非線性噪聲，有效地保護(hù)邊緣信息，方法簡(jiǎn)便，算法也容易實(shí)現(xiàn)。圖5為n=40時(shí)第450層圖像進(jìn)行處理前后的對(duì)比。

圖5 中值濾波前后對(duì)比Fig.5 Comparison before and after median filtering

從圖5可以看出，由于巖石邊緣被其鄰域的中值代替，所以巖石邊緣和顆粒間隙略顯模糊，但前后并沒(méi)有太大的灰度差異，仍能保持圖像中重要細(xì)節(jié)部分。

2.2 圖像直方圖均衡化

圖5(b)的灰度分布直方圖如圖6所示，可以看出圖像的灰分布范圍較窄且極不均勻，直方圖多密集靠近在一起且中間突出一個(gè)高峰，在50～150的像素點(diǎn)個(gè)數(shù)幾乎為0。

圖6 CT圖像灰度直方圖Fig.6 CT Image gray histogram

如若將直方圖分布在整個(gè)水平方向上壓縮，展開(kāi)成一個(gè)新直方圖，其清晰程度會(huì)提高，所需目標(biāo)信息會(huì)更突出。也就是說(shuō)，若圖像存在噪聲，直方圖均衡化會(huì)放大圖像中的噪聲。因此，若使用直方圖均衡化對(duì)圖像進(jìn)行增強(qiáng)，應(yīng)先對(duì)圖像進(jìn)行降噪處理。圖7為依次進(jìn)行中值濾波和直方圖均衡化處理后n=40時(shí)第450層CT圖像。

圖7 均衡化后的CT圖像Fig.7 CT image after equalization

對(duì)比圖5及圖7可知，經(jīng)中值濾波和均衡化處理后的CT圖像整體效果更加鮮明，對(duì)比度也更加明顯，巖石內(nèi)部細(xì)節(jié)部分處理的更好，使得巖石與孔隙對(duì)比更明顯，巖石孔隙更加清晰。

2.3 圖像二值化處理

為了進(jìn)一步觀察巖石顆粒與孔隙的分布特點(diǎn)，本文采用閾值法對(duì)CT圖像進(jìn)行二值化處理，其原理為：先由用戶指定或通過(guò)算法生成一個(gè)閾值，若圖像中某像素點(diǎn)的灰度值小于該閾值，則將該像素灰度值設(shè)置為0或255，否則設(shè)置為255或0。

二值化處理的變換函數(shù)表達(dá)式如下

(1)

式中T為閾值。閾值選擇是灰度圖二值化的關(guān)鍵步驟，常見(jiàn)的閾值選取方法有人工選取法、最大類(lèi)間方差法等。人工選取分割閾值容易受主觀因素影響，而迭代法受圖像目標(biāo)和背景比例差異影響較多，因此本文采用最大類(lèi)間方差法對(duì)圖像進(jìn)行二值化處理。圖8為n=40時(shí)第450層圖像二值化處理后的圖像。

圖8 二值化處理后的CT圖像Fig.8 CT image after binary operation

二值化處理后，圖中白色代表巖石顆粒，黑色代表巖石孔隙和背景，從圖8可以看出，巖石顆粒與孔隙黑白分明，極大地提高了原始CT圖像的直觀性。預(yù)處理完成后凍融次數(shù)n=0,5,10,20,40時(shí)第450層的二值化CT圖像如圖9所示。

圖9 預(yù)處理完成后的CT圖像Fig.9 CT image after all pre-operations

3 基于CT圖像的細(xì)觀特征分析

3.1 巖石孔隙率的分析

為了定量分析巖石孔隙大小的變化規(guī)律，基于最大類(lèi)間方差法對(duì)不同凍融循環(huán)次數(shù)下不同截面的孔隙率進(jìn)行了計(jì)算，并給出其不同凍融次數(shù)下的孔隙率均值，如圖10所示。

圖10 截面計(jì)算孔隙率Fig.10 Computational porosity of different cross-sections

從圖10可以看出，凍融0次時(shí)巖石的孔隙率為13.97%，與試驗(yàn)測(cè)得的孔隙率14.3%較為接近，說(shuō)明了處理方法的可靠性。紅砂巖孔隙率均值隨凍融次數(shù)的增加整體呈上升趨勢(shì)，但不同凍融次數(shù)下試樣孔隙率在軸向的分布比較穩(wěn)定。然而在凍融10次與20次時(shí)，孔隙率增加的并不明顯，且曲線有較長(zhǎng)的交錯(cuò)段。然而同時(shí)進(jìn)行的力學(xué)試驗(yàn)表明[4]，凍融10次與20次試樣的彈性模量及抗壓強(qiáng)度有明顯差異，因此僅對(duì)孔隙率大小變化規(guī)律進(jìn)行研究無(wú)法準(zhǔn)確描述凍融巖石細(xì)觀損傷劣化規(guī)律。

3.2 巖石分形維數(shù)的計(jì)算

分形理論作為一種新方法新理論，正在許多領(lǐng)域應(yīng)用探索，它跳出一維線、二維面、三維體的傳統(tǒng)藩籬，采用連續(xù)維度的思維來(lái)描述物體的屬性與狀態(tài)，進(jìn)而表達(dá)出物體的多樣性與復(fù)雜性。分形維數(shù)是分形理論及其應(yīng)用中最重要的概念之一，是描述物體復(fù)雜度的重要參數(shù)[26]。根據(jù)其不同的定義方法可以分為Hausdorff維數(shù)DH、相似維數(shù)DS、信息維數(shù)Di等等。其中，許多等價(jià)Hausdorff維數(shù)的維數(shù)被提出來(lái)，例如文中使用的計(jì)盒維數(shù)。

3.3 計(jì)盒維數(shù)的計(jì)算方法

采用不同半徑的盒子將圖像中的巖石部分全部覆蓋，當(dāng)選取的盒子尺寸不同，其數(shù)量也不相同。設(shè)小盒子的半徑為ε，盒子的數(shù)量為Nε，選取不同的ε值，假設(shè)兩者滿足

Nε=AεB

(2)

式中A,B為參數(shù)。

對(duì)(2)式變形，可以得到

B=DB+k

(3)

圖像二值化處理完成后，選取不同的ε對(duì)圖像分區(qū)，在每一種網(wǎng)格劃分下，計(jì)算出所有覆蓋圖像中含巖石區(qū)域的網(wǎng)格的數(shù)目，記為網(wǎng)格覆蓋數(shù)，根據(jù)ε取值的不同，可得到一系列尺寸不同的“網(wǎng)格”和相應(yīng)“網(wǎng)格覆蓋數(shù)”的數(shù)據(jù)對(duì)，即Nε-ε數(shù)據(jù)對(duì)，然后在雙對(duì)數(shù)坐標(biāo)系下畫(huà)出數(shù)據(jù)對(duì)，利用最小二乘法擬合為一次函數(shù)，得到一次函數(shù)的斜率k，其相反數(shù)即為圖像的計(jì)盒維數(shù)。

3.4 計(jì)算結(jié)果分析

在CT掃描試驗(yàn)中，經(jīng)歷0次、5次、10次、20次、40次凍融循環(huán)后得到的掃描圖像分別為 884張、899張、775張、925張、905張，由于試件兩端掃描圖像不完整，故選取100～800層的圖像進(jìn)行計(jì)算，在MATLAB軟件中通過(guò)FOR循環(huán)語(yǔ)句依次讀取文件夾中所保存的CT圖像，按照盒維數(shù)的計(jì)算原理編寫(xiě)程序，分別計(jì)算出每一張CT圖像的分形維數(shù)，最后將其保存到統(tǒng)一的XLSX文件中，進(jìn)行ORIGIN作圖，結(jié)果如圖11所示。

圖11 截面計(jì)盒維數(shù)Fig.11 Box dimension of different cross-sections

從圖11可以看出，計(jì)盒維數(shù)在試樣軸向的分布具比較穩(wěn)定，但存在有小幅震蕩，為了挖掘計(jì)盒維數(shù)變化的統(tǒng)計(jì)學(xué)規(guī)律，應(yīng)用移動(dòng)平均法對(duì)圖11中的數(shù)據(jù)進(jìn)行平滑處理，移動(dòng)平均期數(shù)取值為50時(shí)的結(jié)果如圖12所示。

圖12 移動(dòng)平均處理后不同截面分形維數(shù)Fig.12 Box dimensions of different cross sections after moving average operation

從圖11和圖12可以看出，雖然各截面的計(jì)盒維數(shù)不同，但整體上與凍融循環(huán)次數(shù)負(fù)相關(guān)，說(shuō)明巖石內(nèi)部不同位置具有不同的孔隙特征，且隨著凍融循環(huán)次數(shù)的增加，巖石內(nèi)部復(fù)雜的微小孔隙結(jié)構(gòu)逐漸擴(kuò)展并聚合為分布簡(jiǎn)單的較大孔隙。不同凍融次數(shù)下的計(jì)盒維數(shù)平均值如圖13所示。

圖13 孔隙分形維數(shù)均值隨凍融次數(shù)的變化Fig.13 Variation of mean value of pore fractal dimension with freeze-thaw times

從圖13可以看出，0～5次凍融階段，孔隙分形維數(shù)下降速度較慢，僅下降0.18%，孔隙率增加也較小，說(shuō)明凍融初期巖石內(nèi)部孔隙發(fā)展緩慢，凍融損傷增長(zhǎng)不明顯；隨著凍融循環(huán)次數(shù)的增加，巖石內(nèi)部損傷開(kāi)始加速演化，表現(xiàn)為孔隙分形維數(shù)的大幅度降低和孔隙率的快速增加，即內(nèi)部孔隙增長(zhǎng)的同時(shí)也不斷發(fā)生連通和擴(kuò)展，當(dāng)凍融次數(shù)由5次增長(zhǎng)為10次時(shí)巖石內(nèi)部萌生許多新的微孔隙，與此同時(shí)原有的微孔隙也不斷地聚集聯(lián)通，損傷快速演化；隨后凍融次數(shù)由10次增長(zhǎng)到20次的階段，孔隙率的增長(zhǎng)放緩，而孔隙分形仍保持較為高速的降低，說(shuō)明這一階段試樣的凍融損傷主要為孔隙重組；當(dāng)凍融超過(guò)20次以后，全截面孔隙分維均值下降趨勢(shì)也開(kāi)始變緩，僅下降原來(lái)的0.6%，是因?yàn)閹r石孔隙空間已發(fā)展到足夠容納水相變時(shí)的體積膨脹，無(wú)法使孔隙進(jìn)一步擴(kuò)展，凍融循環(huán)作用對(duì)紅砂巖孔隙結(jié)構(gòu)的影響逐漸減弱。

4 結(jié) 論

1)運(yùn)用中值濾波法和直方圖均衡化手段，對(duì)煤礦砂巖的CT圖像進(jìn)行預(yù)處理，極大的降低圖像噪聲，提高圖像的分辨率，為后續(xù)計(jì)盒維數(shù)的計(jì)算奠定基礎(chǔ)。

2)通過(guò)計(jì)算圖像分形維數(shù)，發(fā)現(xiàn)不同凍融次數(shù)下煤礦砂巖的截面分形維數(shù)具有相似性，說(shuō)明初始損傷對(duì)巖石凍融損傷具有持續(xù)性的影響；且隨著凍融次數(shù)的增加，孔隙率逐漸增加，孔隙分形維數(shù)降低，說(shuō)明巖石孔隙發(fā)育的同時(shí)其復(fù)雜度逐漸降低。

3)將煤礦砂巖細(xì)觀損傷特征與凍融荷載聯(lián)系，發(fā)現(xiàn)凍融損傷是由孔隙引起的有效面積損傷和分形維數(shù)表征的結(jié)構(gòu)性損傷共同構(gòu)成。從內(nèi)部細(xì)觀角度揭示巖石的破壞機(jī)制，為指導(dǎo)煤礦安全開(kāi)采提供一定的理論基礎(chǔ)。