基于巖體結構面統(tǒng)計與分維數的巖體質量評價及其在礦山工程中的應用

鞏城城，馬鳳山，張亞民，王善飛，劉自成

（1.中國科學院 工程地質力學重點實驗室 中國科學院 地質與地球物理研究所，北京 100029；2.山東黃金礦業(yè)股份有限公司 三山島金礦，山東 萊州 261400）

0 引言

大量的工程實例證明［1-3］，巖體結構面往往控制著巖體工程的穩(wěn)定性，巖體結構面分布狀況是評價巖體質量好壞的重要指標之一。研究巖體結構面的特征，并且對巖體質量做出評價，對礦山安全具有重要的意義。

以往評價巖體質量的方法主要有RQD法、CSIR法、巴頓巖體分類 Q法以及 RMQ法［4］，這些方法都是定性的或者半定量的評價方法，而巖體結構面的分維數可以作為巖體質量評價的一個定量指標。由Mandelbrot［5］于1982年發(fā)展起來的分形幾何是一門新的數學分支，其主要概念即是分維數，主要特征為自相似性。近年來，國內外學者對分形理論的認識正在逐漸加深［6-8］，該理論的應用范圍也越來越廣。Lee等測定了標準巖石節(jié)理面粗糙系數JRC值范圍的輪廓線分維值［9］；謝和平從統(tǒng)計的角度出發(fā)，得到了JCR與分維數的經驗公式［10］，并對分形損傷力學和分形斷裂力學做了初步總結［11］。此外，還有諸多學者將分形應用于研究巖體結構面的形態(tài)分布特征［12-14］，并將分形與物理實驗［15］和計算機模擬結合起來［16］。

目前，三山島礦區(qū)進入深部開采階段，在開采過程中遇到了坍塌、突水等問題，這就需要對礦井的巖體質量進行分析與評價。選取礦區(qū)-510中段和-555中段為例，根據實測數據，以 MATLAB7.1為平臺，首先對節(jié)理裂隙數據統(tǒng)計、分析，求得結構面特征的統(tǒng)計模型，然后做出Monte-Carlo網絡模擬圖，進而求得圍巖的分維數。從定性和定量的角度對巖體質量進行分析與評價，對礦山的安全生產具有一定的指導意義。

1 礦山地質概況

三山島礦區(qū)位于萊州灣濱海平原，區(qū)內構造以斷裂為主，影響礦區(qū)生產安全的主要是 F1和F3斷裂（圖1）。礦體主要賦存于F1下盤，北西向構造發(fā)育，斷裂帶及附近巖石受擠壓而破碎，掘進時易產生掉塊和塌方。F3斷裂走向北西，產狀近直立，是礦山開拓工程必須通過的主要構造，極易發(fā)生突水、塌方等災害。根據已有的礦山資料，F3斷裂發(fā)育寬度在15～36m，發(fā)育深度很大，在-600m中段依然存在。

2 巖體結構面特征統(tǒng)計分析

以實測數據為基礎，借助 Matlab對巖體結構面特征進行分布模型擬合，然后進行結構面的Monte-Carlo網絡模擬，求得網絡模擬圖。

圖1 研究中段相對位置Fig.1 Relative location of study levels

2.1 結構面產狀的統(tǒng)計模型

采用測窗法和裂隙量測法［17］在現場獲得544組節(jié)理裂隙數據，其中，-510中段320組，-555中段224組。在研究過程中發(fā)現，實測節(jié)理數據受F3斷裂影響較大。因此，以F3斷裂為界，分別研究各中段F3以南和F3以北節(jié)理產狀特征，主要為傾向、傾角、節(jié)理裂隙間距和跡長四個因素。

圖2～圖5為各中段分區(qū)的節(jié)理極點密度圖。經計算：

圖2 －510中段F3以北節(jié)理極點密度圖Fig.2 Contour plot of－510 level，north of F3

-510中段F3以北有三組優(yōu)勢節(jié)理，第一組優(yōu)勢方位 52°∠76°，第二組優(yōu)勢方位 129°∠53°，第三組優(yōu)勢方位 300°∠63°。

-510中段F3以南有四組優(yōu)勢節(jié)理，第一組優(yōu)勢方位 43°∠81°，第二組優(yōu)勢方位 114°∠66°，第三組優(yōu)勢方位 224°∠71°，第四組優(yōu)勢方位 303°∠44°。

圖3 －510中段F3以南節(jié)理極點密度圖Fig.3 Contour p lot of －510 level，south of F3

圖4 －555中段F3以北節(jié)理極點密度圖Fig.4 Contour plot of－555 level，north of F3

圖5 －555中段F3以南節(jié)理極點密度圖Fig.5 Contour plot of－555 level，south of F3

-555中段F3以北有三組優(yōu)勢節(jié)理，第一組優(yōu)勢方位 300°∠57°，第二組優(yōu)勢方位 56°∠74°，第三組優(yōu)勢方位 118°∠61°。

-555中段F3以南有四組優(yōu)勢節(jié)理，第一組優(yōu)勢方位 46°∠77°，第二組優(yōu)勢方位 113°∠63°，第三組優(yōu)勢方位 227°∠75°，第四組優(yōu)勢方位 330°∠59°。

根據以往的經驗統(tǒng)計模型［18-19］，利用 Matlab，分別對-510中段、-555中段的F3以北和F3以南的節(jié)理裂隙數據進行模型擬合分析。研究結果表明，結構面傾向和傾角的多數為正態(tài)分布模型，結構面間距則主要是Weibull分布模型，裂隙跡長主要是對數正態(tài)分布模型（表1）。其中，各分布模型的概率密度函數如下：

均勻分布：

正態(tài)分布：

（k=1時為指數分布）；對數正態(tài)分布：

表1 各中段分區(qū)結構面特征的統(tǒng)計模型Table 1 Rock mass discontinuities’statistical model of both levels

結合各種模型的分布形式及計算所得數據（表1），可以得到 -510中段和 -555中段 F3以南的巷道圍巖節(jié)理發(fā)育較密集，巖體質量較F3以北的質量差，這與實際情況也是相符的，這從定性的角度對巖體質量進行了評價。

2.2 結構面網絡模擬

結構面網絡模擬一般采用實測素描和計算機模擬兩種方法。實測素描方法簡單，但是工作量大，特別是有上萬組數據的時候，繪制網絡圖需很長時間。因此，在研究中采用計算機模擬，這種方法省時省力，并且可以模擬出露頭上看不到的結構面。

在現場結構面調查的基礎上，根據上述所得結構面特征的統(tǒng)計模型，利用Matlab7.1編制Monte-Carlo模擬程序，用計算機自動生成隨機數代替實際節(jié)理裂隙數據，計算出各隨機數生成的裂隙網絡的控制點坐標，將其連線，生成各中段F3以北及以南的裂隙網絡圖。其中，模擬范圍用長×寬來表示。在模擬中，寬對應實際巷道的高度，分別取為5m、10m等進行模擬計算；長對應巷道的走向長度，分別取10m、20m、30m等進行模擬計算，以求最優(yōu)解。如圖6～圖8是模擬區(qū)域分別是 10m×5m、20m×10m、40m×5m時，-510中段F3以北的節(jié)理網絡模擬圖。

圖6 －510中段F3以北10×5（m2）范圍節(jié)理網絡模擬圖Fig.6 The fissures’network simulation of －510 level for 10×5（m2），north of F3

圖7 －510中段F3以北20×10（m 2）范圍節(jié)理網絡模擬圖Fig.7 The fissures’network simulation of －510 level for 20×10（m 2），nor th of F3

圖8 －510中段F3以北40×5（m 2）范圍節(jié)理網絡模擬圖Fig.8 The fissures’network simulation of －510 level for 40×5（m 2），north of F3

3 分維數的確定及分析

3.1 分維數的確定

分維數是定量描述分形幾何特征和復雜程度的重要參數。目前可用于巖體結構面分維數測定的方法眾多，主要有碼尺法，盒計維數法，自仿射分形法（方差法），冪律譜法和周長-面積關系法等。

在研究中采用盒計維數法，這種方法抓住了分形結構的特點?；驹砣缦拢菏褂眠呴L為r的正方形去覆蓋分形曲線（圖6～圖8），所需盒子的數目將隨r的值而變化。如果所量測的對象具有分形特性，則N（r）-r滿足以下關系：

N（r）=ar-D

式中a為常數。對上式兩邊取對數，得：

其中，D即是分維數。

利用Matlab編制程序，計算覆蓋的盒子數，覆蓋到的部分記為1，沒覆蓋到的盒子記為0。變化 r值，得到一組N（r）-r數據，然后繪制出雙對數 ln r-ln［N（r）］曲線，最后用最小二乘法求出 ln r-ln［N（r）］線性關系的斜率k（圖9），此斜率的絕對值即是所求的分維數D。

選取不同的模擬范圍5×5（m2）、10×5（m2）、10×10（m2）、20 ×5（m2）、20 ×10（m2）、20 ×20（m2）、30×5（m2）、30 ×10（m2）、30 ×20（m2）、30 ×30（m2）、40×5（m2）、40 ×10（m2）、50 ×5（m2）、60 ×5（m2）、70 ×5（m2）、80 ×5（m2），進行結構面的分形研究與分析。研究發(fā)現，相同模擬區(qū)域的分維數基本相同，如-510中段 F3以北，模擬區(qū)域為10×10（m2）時分維數是 1.46001，模擬區(qū)域為 20×5（m2）是1.46527；模擬區(qū)域為30×10（m2）的分維數是1.64619，模擬區(qū)域為 60×5（m2）的分維數是1.66018；模擬區(qū)域為20×20（m2）的分維數是1.66186，模擬區(qū)域為 40×10（m2）的分維數是1.67414，模擬區(qū)域為 80×5（m2）的分維數是1.69832。在此，認為相似的原因是結構面特征模型相同，而每次模擬產生的隨機數不可能完全相同，則導致相同模擬區(qū)域的分維值不完全一致。另外可以看出，模擬區(qū)域面積一定時，區(qū)域的長寬比增大，分維數隨之緩慢增加。

表2 研究區(qū)分維數計算結果Table 2 Fractal calculated value of study area

根據模擬結果可得，隨著模擬區(qū)域的增大，分維數呈現同步增大趨勢，即出現升維現象，并且增大到某一值后趨于穩(wěn)定。圖10為-510中段F3以北不同模擬范圍的分維值變化曲線圖，由圖可知，模擬范圍是40×5（m2）時，分維值是1.5981，模擬范圍是80×5（m2）時，分維值是1.69832，變化范圍很小。當模擬區(qū)域達到一定的尺寸后，認為是合理的計算尺寸。另外，模擬范圍相同時，長寬比越大，分維值越大，相應的巖石完整性越差［20］。因此，選取40×5（m2）作為合理的模擬區(qū)域。表2是研究區(qū)的計算分維值。

圖9 －510中段 F3以北 lnr－ln［N（r）］相關關系圖Fig.9 Correlation of lnr－ln［N（r）］of－510 level，north of F3

圖10 －510中段F3以北不同模擬范圍的分維數Fig.10 Fractal dimension of different simulation area of －510 level，north of F3

3.2 分維數與巖體質量關系

各中段分區(qū)的ln r-ln［N（r）］曲線相關性好（見表2），首先說明了研究區(qū)具有分形特征。以往學者在分形理論上已做過大量研究，但由于計算方法、參數選取以及研究區(qū)域的不同，提出了不同的巖體質量評價的分維數標準［21］。然而，相同的是，分維值越大，巖體質量越差。

結合本次研究，可知巖體結構面分布的分維數是衡量結構面分布的指標，大小能反映巖體質量的好壞程度。分維值D越大，結構面分布越密集，且優(yōu)勢結構面越多，巖體質量越差。表2給出了巖體質量的定量指標，從中可以得出各中段的巖體質量從好到差依次為：-510中段F3以北，-555中段F3以北，-555中段F3以南，-510中段F3以南。

分維值表征的巖體質量好壞與實際情況是相符合的，這與巖體結構面特征的統(tǒng)計模型所得結果一致，相互驗證。根據對巖體質量的定性和定量分析結果，為保證礦山的安全生產，對F3以南的巖體穩(wěn)定性的防護措施應強于F3以北的巖體防護措施。然而，對于采取何種級別的防范措施，還應當配合其它方法的綜合分析與評價。

4 結論

（1）以實測數據為基礎，對巖體結構面特征進行統(tǒng)計模型分析，結果顯示-510中段和-555中段F3以南的巖體質量較F3以北的巖體質量差，這給出了巖體質量的定性評價。

（2）巖體結構面分布的分維數大小能反映巖體質量的好壞程度。分維值D越大，結構面分布越密集，且優(yōu)勢結構面越多，巖體質量越差。研究區(qū)的圍巖質量的分維數為1.5981至1.6668，進行了巖體質量的定量分析。

（3）分維數與裂隙網絡模擬區(qū)域有密切聯系。相同的模擬區(qū)域得到的分維值基本相同，并且隨著模擬區(qū)域的長寬比增大而緩慢增加。另外，隨著模擬區(qū)域的增大，會出現升維現象。當模擬區(qū)域增大到一定范圍時，分維數趨于穩(wěn)定。

（4）為提供可靠的巖體質量評價，保證礦山的生產安全，還需配合其它方法，進行綜合分析與評價。

［1］盛永清，周創(chuàng)兵，陳益豐，等.平行 F42-9發(fā)育的結構面對錦屏水電站高邊坡穩(wěn)定性的影響［J］.巖土力學，2008，29（10）：2613-2618.SHENG Yongqing，ZHOU Chuangbing，CHEN Yifeng，et al.Effect of occurred texture planes parallel to fault F42-9on stability of high slope at Jinping Hydropower Station［J］.Rock and Soil Mechanics，2008，29（10）：2613-2618.

［2］師金峰，張應龍.超大斷面隧道圍巖的穩(wěn)定性分析［J］.地下空間與工程學報，2005，1（2）：227-230，241.SHI Jingfeng，ZHANG Yinglong.FEM analysis of rock stability on super section tunnel［J］.Chinese Journal of Underground Space and Engineering，2005，1（2）：227-230，241.

［3］樊啟祥，王義峰.向家壩水電站地下廠房緩傾角層狀圍巖穩(wěn)定性分析［J］.巖石力學與工程學報，2010，29（7）：1307-1313.FAN Qixiang，WANG Yifeng.Stability analysis of layered surrounding rock mass of large underground powerhouse of Xiangjiaba hydropower station ［J］.Chinese Journal of Rock Mechanics and Engineering，2010，29（7）：1307-1313.

［4］李智毅，楊峪云.工程地質學概論［M］.武漢：中國地質大學出版社，1994.LIZhiyi，YANG Yuyun.Overview of engineering geology［M］.Wuhan：China University of Geosciences，Press，1994.

［5］Mandelbrot B B.The fractal geometry of nature［M］.San Francisco：W H Freeman＆ Company，1982.

［6］Aguilar-Hernandez A，Ramirez-Santiago G.Self-similar and self-affine properties of two-dimensional fracture patterns in rocks［J］.Mathematical Geosciences，2010，42（8）：925-954.

［7］Park S L，Kim Y S，Ryoo C R，et al.Fractal analysis of the evolution of a fracture network in a granite outcrop，SE Korea ［J］.Geosciences Journal，2010，14（2）：201-215.

［8］Feng Z C，Zhao Y S，Zhao D.Investigating the scale effects in strength of fractured rock mass［J］.Chaos Solitons＆ Fractals，2009，41（5）：2377-2386.

［9］Lee Y J，Carr J R，Barr D J，et al.The fractal dimensions a measure of the roughness of rock discontinuity profiles ［J］.International Journal of Rock Mechanics and Mining Sciences＆Geomechanics Abstracts，1990，27（6）：453-464.

［10］謝和平，W G Pariseau.巖石節(jié)理粗糙系數（JRC）的分形估計［J］.中國科學（B 輯），1994，24（5）：524-530.XIE Heping.W G Pariseau.Rough coefficient of rockjoints fractalestimate［J］.Science in China（Series B），1994，24（15）：524-530.

［11］謝和平.分形 -巖石力學導論［M］.北京：科學出版社，1995.XIE Heping，Fractal-rock-mechamics introduction［M］.Beijing：Science Publishing Company，1995.

［12］徐光黎.巖體結構面幾何特征的分形與分維［J］.水文地質工程地質，1993，20（2）：20-22.XU Guangli.Fractal ananysis for rock mass joint geometry［J］.Hydrogeology and Engineering Geology，1993，20（2）：20-22.

［13］Amitava Ghosh，Jaak JK Daemen.Fractal characteristics of rock discontinuities ［J］.Engineering Geology，1993，34（1）：1-9.

［14］Carpinteri A，Lacidogna G，Niccolini G.Fractal analysis of damage detected in concrete structural elements under loading［J］.Chaos Solitons ＆ Fractals，2009，42（4）：2047-2056.

［15］謝和平，于廣明，楊倫，等.采動巖體分形裂隙網絡研究［J］.巖石力學與工程學報，1999，18（2）：147-151.XIE Heping，YU Guangming，Yang Lun，et al.Research on the fractal effects of crack network in overburden rock stratum ［J］.Chinese Journal Rock Mechanics and Engineering，1999，18（2）：147-151.

［16］張向東，徐崢嶸，蘇仲杰，等.采動巖體分形裂隙網絡計算機模擬研究［J］.巖石力學與工程學報，2001，20（6）：809-812.ZHANG Xiangdong，XU Zhengrong，SU Zhongjie，et al.Research of the imitation about the fractal crack in overburden rock stratum by computer ［J］.Chinese Journal Rock Mechanics and Engineering，2001，20（6）：809-812.

［17］賈洪彪，唐輝明，劉佑榮，等.巖體結構面三維網絡模擬理論與工程應用［M］.北京：科學出版社，2008.JIA Hongbiao，TANG Huiming，LIU Yourong，et al.Theory and engineering application of 3-d network modeling of discontinuities in rockmass［M］.Beijng：Science Press，2008.

［18］伍法權.統(tǒng)計巖石力學原理［M］.武漢：中國地質大學出版社，1993.WU Faquan.Principles of statistical mechanics of rock masses［M］.Wuhan：China University of Geosciences Press，1993.

［19］Priest S D，Hudson J A.Estimation of discontinuity spacing and trace length using scanline surveys［J］.International Journal of Rock Mechanics and Mining Sciences and Geomechanics Abstracts，1981，18（3）：183-197.

［20］范留明，黃潤秋.一種估計結構面跡長的新方法及其工程應用［J］.巖石力學與工程學報，2004，23（1）：53-57.FAN Liuming，HUANG Runqiu.New method for estimating trace length of discontinuities and its engineering application ［J］.Chinese Journal Rock Mechanics and Engineering，2004，23（1）：53-57.

［21］龐大鵬，陳劍平，王丹微.節(jié)理網絡分形在隧道超前地質預報中的應用［J］.巖土力學，2009，30（5）：1415-1420.PANG Dapeng，CHEN Jianping，WANG Danwei.Application of fractal description of rock mass joint and cracks network to advanced geological forecast of tunnels［J］.Rock and Soil Mechanics，2009，30（5）：1415-1420.