張 亞,彭樂文

(河海大學地球科學與工程學院,江蘇 南京 211100)

Monte-Carlo隨機模擬技術(shù)是通過一定的隨機數(shù)生成方法,生成服從某一隨機變量的概率分布形式的隨機數(shù)序列。由于裂隙的分布具有很大的隨機性、不確定性,可以把裂隙幾何參數(shù)當作一種隨機變量,并且大量研究表明裂隙幾何參數(shù)的分布服從一定概率分布形式,因此可以將Monte-Carlo隨機模擬技術(shù)運用到裂隙巖體中來研究裂隙分布特征以及裂隙網(wǎng)絡模型。離散裂隙網(wǎng)絡(DFN,discrete fracture network)將巖體切割成很多由相鄰結(jié)構(gòu)面接觸而組合的塊體,塊體的變形及運動均由結(jié)構(gòu)面控制,并且裂隙網(wǎng)絡為地下水流唯一的滲流通道。對裂隙巖體來說,離散裂隙網(wǎng)絡模型比等效連續(xù)介質(zhì)模型更能刻畫地下水實際的滲流情況[1]。

隨著計算機技術(shù)的發(fā)展,隨機裂隙網(wǎng)絡的模擬也成為當今研究的熱點之一。王晉麗等[2]利用Monte-Carlo隨機模擬技術(shù)生成二維裂隙網(wǎng)絡滲流模型,研究了裂隙網(wǎng)絡中水流在不同邊界條件下的流動特征;李新強等[3]根據(jù)結(jié)構(gòu)面的分布特征利用計算機模擬技術(shù)建立了裂隙網(wǎng)絡模型,并對該模型進行了等效巖體滲透張量的求解;張彥洪等[4]通過人工模擬方法建立隨機裂隙網(wǎng)絡模型,研究了裂隙開度和滲透壓力對裂隙巖體滲流應力耦合特性的影響。在有限的地質(zhì)測量基礎上,馮學敏等[5]運用Monte-Carlo方法建立了與實際巖體裂隙同分布的隨機裂隙網(wǎng)絡;敖雪菲等[6]在研究壩基裂隙巖體灌漿時,基于Monte-Carlo方法建立了裂隙巖體三維網(wǎng)絡模型,結(jié)果表明該法能有效反應裂隙真實特征并獲得貼近實際的灌漿結(jié)果;劉日成等[7]通過隨機模擬技術(shù)建立了2種DFN模型,考慮2種邊界條件研究了DFN的非線性滲流特征;劉泉聲等[8]基于裂隙網(wǎng)絡模擬技術(shù)Monte-Carlo法,通過Fish語言編制DFN模型生成程序DFN-GEN,研究了應力對等效滲透系數(shù)的影響。DFN建模技術(shù)應用領域比較廣泛,不僅在裂隙巖體滲流方面得到應用,而且郎曉玲等[9]、劉廣峰等[10]也將DFN模型建模技術(shù)運用到了油田開發(fā)等方面。

根據(jù)大量的實踐經(jīng)驗,可以發(fā)現(xiàn)裂隙分布非常復雜,并且裂隙幾何參數(shù)難以確定,但許多研究表明在裂隙網(wǎng)絡模擬中假設結(jié)構(gòu)面為薄圓盤形狀是比較符合實際情況的[11]。因此,研究利用Monte-Carlo隨機模擬技術(shù)建立裂隙為薄圓盤形狀的三維裂隙網(wǎng)絡模型,并對模型進行驗證、統(tǒng)計窗裂隙交線分析以及裂隙網(wǎng)絡連通率計算。

1 Monte-Carlo隨機模擬原理

Monte-Carlo隨機模擬方法的基本原理是利用[0,1]區(qū)間標準均分隨機數(shù),根據(jù)由已知裂隙幾何參數(shù)的密度分布函數(shù)求得的抽樣公式來獲得服從給定裂隙幾何參數(shù)分布形式的隨機變量。

1.1 標準均勻分布隨機數(shù)的產(chǎn)生

目前,在工程應用中多運用數(shù)學方法通過計算機編程來產(chǎn)生所要的隨機數(shù)。同余法是產(chǎn)生標準均勻分布隨機數(shù)最為常用的數(shù)學方法之一,它包括乘法線性同余法、混合線性同余法、二次同余法和加法同余法(見表1),其中研究使用的是混合線性同余法。

表1 常用的同余法

1.2 標準均勻分布隨機數(shù)檢驗

隨機數(shù)檢驗是對利用隨機數(shù)生成器產(chǎn)生的偽隨機數(shù)序列與真正的[0,1]均勻分布隨機數(shù)序列相似程度的檢驗。常用的檢驗方法包括:參數(shù)檢驗、均勻性檢驗、獨立性檢驗等,其中均勻性檢驗又包括χ2檢驗和K-S檢驗。統(tǒng)計檢驗不能完全肯定或否定某一假設,最好能夠多做幾種統(tǒng)計檢驗,以便有較大的把握保證使用的隨機數(shù)列有較好的統(tǒng)計性質(zhì)。研究所用的隨機數(shù)檢驗方法如下:

(1)

該統(tǒng)計量公式漸近地服從χ2(m-1) ,通過查χ2分布表,即可對頻率差異性做顯著性檢驗。

(2) 獨立性檢驗 獨立性檢驗是檢驗一組隨機數(shù)相互之間是否存在相關性。設u1,u2,…,ui為一組按先后順序排列的隨機數(shù)序列,前后距離為k的隨機數(shù)相關系數(shù)為

(2)

因此,統(tǒng)計量為

(3)

當在n充分大時,Vk近似服從N(0,1)分布。

1.3 隨機變量抽樣

在Monte-Carlo隨機模擬中,隨機變量是通過在標準均勻分布隨機數(shù)的基礎上經(jīng)過一定的概率分布函數(shù)公式變換處理得到的,并且服從不同分布形式,這一過程稱為隨機變量抽樣。抽樣的方法主要包括:直接抽樣法(又稱反函數(shù)法)、復合抽樣法、舍選抽樣法和變換抽樣法等,其中直接抽樣法最為常用和有效。

直接抽樣法的基本思路為:設隨機變量t服從分布p(t),累積分布函數(shù)P(t),則P(t)的值域為[0,1],所以可以均勻分布隨機數(shù)ui作為P(ti)的函數(shù)值。根據(jù)累積分布函數(shù)與分布函數(shù)的關系,在P(ti)和ui之間建立t與u對應關系:

(4)

求式(4)的反函數(shù)得

t=p-1(u),

(5)

將p(t)代入式(4),再將式(4)代入求反函數(shù)之后的式(5),即可得到抽樣表達式。幾種常見分布形式的抽樣表達式如表2所列。

表2 常見分布的抽樣表達式

2 結(jié)構(gòu)面分布規(guī)律及幾何特征

2.1 結(jié)構(gòu)面分布規(guī)律

雖然結(jié)構(gòu)面的發(fā)育具有隨機性,但是同組結(jié)構(gòu)面在巖體內(nèi)分布具有一定的規(guī)律可循,即結(jié)構(gòu)面具有等距性和韻律性。結(jié)構(gòu)面的分布規(guī)律不僅為研究巖體結(jié)構(gòu)力學和巖體水力學模型提供了依據(jù),還有助于研究結(jié)構(gòu)面的力學特性。

結(jié)構(gòu)面的等距性是指同組結(jié)構(gòu)面在一定尺度上、一定范圍內(nèi)具有相同的間距。無論是大尺度還是小尺度,這種等距現(xiàn)象都存在。結(jié)構(gòu)面的等距性主要受以下4種因素的影響:巖石力學性質(zhì)、應力強度、構(gòu)造作用和巖層厚度。除了等距性,結(jié)構(gòu)面還表現(xiàn)出韻律性,即裂隙發(fā)育密集帶和稀疏帶交替相見出現(xiàn),且?guī)c帶之間具有一定間距。結(jié)構(gòu)面韻律性分布的形成機制一直是學者研究的熱點,目前公認的2種理論模式為:飽和模式與傳播模式。

2.2 結(jié)構(gòu)面幾何特征

裂隙巖體中力學特征及水力特征受裂隙控制,而裂隙的空間分布又可以通過幾何特征來表征。結(jié)構(gòu)面的幾何特征主要是指結(jié)構(gòu)面的形態(tài)以及結(jié)構(gòu)面幾何參數(shù),其中結(jié)構(gòu)面幾何參數(shù)又包括方位、規(guī)模、張開度和粗糙度等。

(1) 結(jié)構(gòu)面形狀 結(jié)構(gòu)面形狀是指在忽略厚度的情況下,在延伸平面上的幾個形狀。由于地質(zhì)條件的復雜性,在不同形成條件下、結(jié)構(gòu)面的形狀完全不同,很難用一個統(tǒng)一的形狀來概括,因此許多學者提出了不同的概念模型,主要包括:正交模型(見圖1)、圓盤模型(見圖2)和多邊形模型等。

圖1 結(jié)構(gòu)面正交模型Fig.1 Structural surface orthogonal model

圖2 結(jié)構(gòu)面圓盤模型Fig.2 Structural surface disc model

(2) 結(jié)構(gòu)面參數(shù) 結(jié)構(gòu)面的幾何參數(shù)一般包括:產(chǎn)狀、間距、跡長、粗糙度、張開度、充填物和連通性等,其中間距是指相鄰結(jié)構(gòu)面的垂直距離,用來反映結(jié)構(gòu)面發(fā)育的密集程度,另外根據(jù)所測的平均間距可將巖體進行描述(見表3);跡長是指結(jié)構(gòu)面與露頭面的交線,用來描述和評價結(jié)構(gòu)面的連續(xù)性,國際巖石力學學會提出的分級標準見表4;張開度是指裂隙兩壁間的垂直距離,是表征滲透性重要的幾何參數(shù)之一。

表3 結(jié)構(gòu)面間距分級表

表4 結(jié)構(gòu)面連續(xù)性分級表

3 裂隙網(wǎng)絡模型的建立

3.1 結(jié)構(gòu)面分布的概率分析

裂隙發(fā)育具有很大的隨機性,即裂隙幾何參數(shù)具有隨機性,屬于隨機變量,因此可以用相應的概率分布來描述。目前,常見的裂隙幾何參數(shù)的概率分布函數(shù)包括:均勻分布、負指數(shù)分布、正態(tài)分布和對數(shù)正態(tài)分布等[12]。

根據(jù)大量的野外資料表明,結(jié)構(gòu)面的產(chǎn)狀、跡長、間距及張開度等幾何要素服從于某種隨機分布,其中均勻分布一般適用于產(chǎn)狀的分布形式,負指數(shù)分布一般適用于跡長、間距和張開度分布形式,如圖3所示;正態(tài)分布一般適用于產(chǎn)狀和跡長分布形式,如圖4所示;對數(shù)正態(tài)分布一般適用于跡長、間距和張開度分布形式。常見的幾何要素概率分布規(guī)律見表5[13]。

圖3 間距的負指數(shù)分布形式Fig.3 Negative exponential distribution of spacing

3.2 DFN的建立

野外巖體裂隙分布是復雜的,大量的工程實踐證明,裂隙的分布特征服從一定的規(guī)律。因此,可以通過野外露頭、平洞等裂隙的統(tǒng)計分析,建立裂隙分布的概率模型,并通過Monte-Carlo隨機模擬方法,在計算機上建立反映服從該概率分布規(guī)律的三維裂隙網(wǎng)絡模型。裂隙網(wǎng)絡模擬過程是根據(jù)現(xiàn)場裂隙統(tǒng)計資料獲得裂隙幾何參數(shù)的概率分布規(guī)律,應用Monte-Carlo隨機模擬方法,對已知的密度函數(shù)進行隨機抽樣,進而得到與實際分布函數(shù)相應的隨機變量,進而推算出每條裂隙的空間位置坐標,并結(jié)合裂隙幾何參數(shù)在計算機上生成形象直觀的三維裂隙網(wǎng)絡。

圖4 傾向的正態(tài)分布形式Fig.4 Normal distribution of tendency

實際上,隨機裂隙網(wǎng)絡建模的過程就是根據(jù)露頭得到有限的裂隙信息推求服從這些概率分布特征的巖體內(nèi)部以及更大區(qū)域的裂隙分布,建立能夠反映現(xiàn)場實際情況的裂隙網(wǎng)絡模型[12]?；镜娜S裂隙網(wǎng)絡模擬流程如圖5所示。

4 實例分析

4.1 工程概況

根據(jù)文獻[1]中某水電站壩址區(qū)裂隙的幾何參數(shù)統(tǒng)計資料,并結(jié)合Monte-Carlo模擬技術(shù)建立三維裂隙網(wǎng)絡模型。裂隙的幾何參數(shù)如表6所列,其優(yōu)勢方向為北東向和北西向。對裂隙的幾何分布進行擬合,部分擬合結(jié)果如圖6、圖7所示。

表5 結(jié)構(gòu)面幾何要素經(jīng)驗概率分布

圖5 三維裂隙網(wǎng)絡模擬流程Fig.5 3D crack network simulation flowchart

圖6 間距(NE)頻數(shù)直方圖(負指數(shù)分布)Fig.6 Frequency histogram (Negative exponential distribution) of crack spacing (NE)

圖7 傾向(NE)頻數(shù)直方圖(正態(tài)分布)Fig.7 Frequency histogram (Normal distribution) of tendency (NE)

4.2 巖體三維裂隙網(wǎng)絡的模擬

設置裂隙生成區(qū)域為200×200×200 m3,正北向為X軸,正東向為Y軸,垂直向為Z軸。由于該區(qū)域發(fā)育2組優(yōu)勢裂隙,故將NE向、NW向分布命名為DFN01和DFN02。2組裂隙的數(shù)量大小由結(jié)構(gòu)面體密度控制,而求得所需要生成的2組裂隙數(shù)量的方法可參考文獻[14]中所述方法,主要思路如下:

由線密度λd和結(jié)構(gòu)面間距d的定義可知兩者互為倒數(shù),則有

(6)

(7)

表6 某水電站壩址區(qū)裂隙幾何參數(shù)統(tǒng)計

(8)

因此,模擬區(qū)域V內(nèi)的結(jié)構(gòu)面數(shù)量n為

(9)

經(jīng)計算可知,NE向裂隙組DFN01所需生成裂隙數(shù)量為307條,NW向DFN02所需61條,由此生成的三維裂隙網(wǎng)絡模型如圖8所示,其中DFN01標記為藍色,DFN02標記為綠色。

圖8 三維裂隙網(wǎng)絡模型Fig.8 3D crack network model

4.3 模型的驗證與分析

(1) 模型的驗證 根據(jù)裂隙網(wǎng)絡模型求得所生成DFN的幾何參數(shù)信息,結(jié)合現(xiàn)場裂隙統(tǒng)計得到的裂隙幾何特征對通過Monte-Carlo隨機模擬技術(shù)建立的三維裂隙網(wǎng)絡進行驗證。模型驗證一般可以考慮對比模型數(shù)據(jù)與實測數(shù)據(jù)之間結(jié)構(gòu)面平均傾角傾向、平均跡長、線密度等參數(shù)。在模型中設置4條測線來計算所建DFN的線密度,如圖9所示。對比模型數(shù)據(jù)與實際數(shù)據(jù)之間線密度的差異性,并計算相對誤差百分比,結(jié)果見表7。

圖9 測線及統(tǒng)計窗布置示意圖Fig.9 Measuring line and statistical window layout map

由表7裂隙網(wǎng)絡線密度模擬值與實測值對比來看,相對誤差百分比均≤10%,認為該模型模擬結(jié)果比較貼近實際情況。

(2) 模型分析 研究對模型的分析主要包括統(tǒng)計窗裂隙交線分析和裂隙連通率計算分析。統(tǒng)計窗設置如圖9中紅色平面所示,統(tǒng)計窗坐標分別為(-200,-200,0)(-200,200,0)(200,200,0)(200,-200,0)。統(tǒng)計窗裂隙交際線如圖10所示。統(tǒng)計窗上藍色裂隙交線為DFN01與統(tǒng)計窗的交線,綠色為DFN02與統(tǒng)計窗的交線。

由圖10可知,DFN01和統(tǒng)計窗相交而成的交線與DFN02和統(tǒng)計窗交線的夾角約為70°,這與表6中所示的2組裂隙走向的差值大致一致。DFN裂隙連通率主要是通過裂隙間的交線的長度來反映的,即交線長度越長,裂隙連通率越大,裂隙連通率計算如圖11所示。由圖11可知,大部分裂隙交線的走向為NE向,這是因為DFN01數(shù)量較多,約為DFN02數(shù)量的5倍;裂隙連通率較大的交線是2組優(yōu)勢裂隙中半徑較大的裂隙相交而形成。

表7 裂隙網(wǎng)絡線密度模擬值與實測值對比

圖10 統(tǒng)計窗上裂隙交線分布Fig.10 Distribution of the crack intersection on the statistical window

圖11 裂隙連通率計算Fig.11 Crack connectivity calculation

5 結(jié)語

研究采用Monte-Carlo隨機模擬方法建立了三維裂隙網(wǎng)絡模型,主要通過混合線性同余法生成隨機數(shù),并且對隨機數(shù)進行χ2檢驗和獨立性檢驗,結(jié)合實測裂隙統(tǒng)計與分析得到的裂隙幾何參數(shù)概率分布的抽樣公式對隨機變量進行抽樣,得到與實際分布函數(shù)相應的隨機變量,并且假設結(jié)構(gòu)面為薄圓盤形狀,進而建立DFN模型。

在裂隙網(wǎng)絡模擬之后,需要進行模型的驗證,通過對比實測數(shù)據(jù)與模擬數(shù)據(jù)之間的線密度參數(shù)的差異性來檢驗模型的可靠性,驗證結(jié)果表明模型擬合度較高。建立模型之后,除了對模型統(tǒng)計窗裂隙交線分析和裂隙連通率計算分析,還可以近似地求得RQD、滲透系數(shù)張量等參數(shù)。利用Monte-Carlo隨機模擬技術(shù)對三維裂隙網(wǎng)絡模型的研究,不僅為裂隙巖體的建模提供了一種很好的思路,并且為后續(xù)對裂隙巖體滲流的研究工作提供了技術(shù)支持。