李愛華,朱 江

(中國石油大學(xué)(華東)石油工程學(xué)院,山東 青島 266580)

巖石裂隙中水流脈動(dòng)壓力傳播的數(shù)值模擬

李愛華,朱 江

(中國石油大學(xué)(華東)石油工程學(xué)院,山東 青島 266580)

為了分析網(wǎng)絡(luò)結(jié)構(gòu)對(duì)裂隙水體內(nèi)脈動(dòng)壓力傳播規(guī)律的影響,引入蒙特卡羅模擬方法,通過數(shù)值模擬得到與天然裂隙網(wǎng)絡(luò)在統(tǒng)計(jì)上相似的仿真裂隙網(wǎng)絡(luò),以此作為裂隙水的賦存空間,建立一維瞬變流模型,對(duì)脈動(dòng)壓力在裂隙網(wǎng)絡(luò)中的傳播進(jìn)行數(shù)值模擬,得出裂隙水中脈動(dòng)壓力傳播的非定常時(shí)間過程,以及裂隙網(wǎng)絡(luò)中壓力脈動(dòng)強(qiáng)度的空間分布云圖。模擬結(jié)果表明，巖體裂隙網(wǎng)絡(luò)中脈動(dòng)壓力的傳播具有顯著的各向異性,裂隙中各點(diǎn)的瞬時(shí)脈動(dòng)壓力是由該點(diǎn)所銜接(尤其是直接銜接)的裂隙線元的水力瞬變作用以及裂隙入口干擾源、裂隙末端反射源三者共同作用的結(jié)果。

脈動(dòng)壓力;裂隙網(wǎng)絡(luò);蒙特卡羅方法;瞬變流模型;數(shù)值模擬

沖擊射流的脈動(dòng)壓力在縫隙中的傳播是導(dǎo)致溢洪道襯砌板塊揭底破壞的主要原因。1992年Fiorotto等[1]提出板塊縫隙中脈動(dòng)壓力傳播的瞬變流模型,現(xiàn)已被廣泛應(yīng)用于脈動(dòng)壓力傳播的數(shù)值研究中[2-7]。Jia等[8-10]對(duì)縫隙中的脈動(dòng)壓力傳播進(jìn)行了試驗(yàn)研究。但無論是試驗(yàn)研究還是數(shù)值研究均局限于簡(jiǎn)單一維規(guī)則裂隙或者二維規(guī)則裂隙面內(nèi),并沒有涉及實(shí)際基巖的復(fù)雜結(jié)構(gòu)特性。Annandale[11]在1995年力圖從水力和地質(zhì)兩個(gè)方面研究工程土石料的抗沖特性,定義了與土石料的強(qiáng)度、顆粒尺度、不連續(xù)性等相關(guān)的抗沖刷系數(shù)。但歸根結(jié)底,抗沖刷系數(shù)是由經(jīng)驗(yàn)系數(shù)確定的,同樣沒有涉及巖體的內(nèi)部網(wǎng)絡(luò)結(jié)構(gòu)。

事實(shí)上,巖體是由巖塊和分割巖塊的裂隙網(wǎng)絡(luò)所組成的結(jié)構(gòu)體。巖體裂隙網(wǎng)絡(luò)可以作為水體的賦存空間,因此裂隙網(wǎng)絡(luò)的結(jié)構(gòu)特性——如巖塊結(jié)構(gòu)面密度、粗糙度、開度、填充情況等,在很大程度上決定了巖體裂隙介質(zhì)中水流脈動(dòng)壓力傳播的特殊性和復(fù)雜性。巖體復(fù)雜裂隙網(wǎng)絡(luò)的提出始于對(duì)地下水滲流的研究,用于解決工程中的裂隙巖體滲流問題。國內(nèi)潘別桐等[12-13]在統(tǒng)計(jì)學(xué)和概率論基礎(chǔ)上提出了基于蒙特卡羅(Monte-Carlo)模擬技術(shù)的結(jié)構(gòu)面二維、三維裂隙網(wǎng)絡(luò)模擬方法,該法可以綜合反映裂隙系統(tǒng)的極不均勻性、滲透水流的各向異性和局部非連續(xù)性,還可以描述裂隙中水流瞬時(shí)變化的特征,因此可以應(yīng)用裂隙網(wǎng)絡(luò)模擬技術(shù)來研究巖體裂隙中壓力波傳播等問題。蒙特卡羅模擬技術(shù)現(xiàn)已被廣泛應(yīng)用于裂隙滲流研究領(lǐng)域[14-15]。

本文通過引入蒙特卡羅方法,模擬生成與天然巖體裂隙網(wǎng)絡(luò)在統(tǒng)計(jì)上相似的仿真裂隙網(wǎng)絡(luò),并以此作為水體賦存空間,建立一維瞬變流模型,對(duì)裂隙水體中脈動(dòng)壓力的傳播過程和分布規(guī)律進(jìn)行數(shù)值模擬和分析。在對(duì)裂隙水流脈動(dòng)壓力波傳播規(guī)律的研究中,引入蒙特卡羅方法,對(duì)于考察巖體裂隙結(jié)構(gòu)對(duì)脈動(dòng)壓力傳播的影響具有重要意義。

1 巖體裂隙網(wǎng)絡(luò)的蒙特卡羅模擬

大自然中的巖體裂隙可劃分為以下4個(gè)等級(jí):①一級(jí),真實(shí)裂隙網(wǎng)絡(luò)——規(guī)模較大的斷層、斷裂和軟弱夾層;②二級(jí),隨機(jī)裂隙網(wǎng)絡(luò)——分布具有隨機(jī)性的次一級(jí)裂隙網(wǎng)絡(luò);③三級(jí),等效連續(xù)介質(zhì)體系——規(guī)模較小、數(shù)量眾多且分布密集的小裂隙;④四級(jí),連續(xù)介質(zhì)體系——巖塊內(nèi)的孔隙網(wǎng)絡(luò)。在水利工程中,對(duì)一級(jí)真實(shí)裂隙需要進(jìn)行特殊處理;三、四級(jí)裂隙規(guī)模較小,在滲流研究中通常作為等效的連續(xù)體或連續(xù)體處理。因此,本文研究基巖裂隙中脈動(dòng)壓力的傳播以二級(jí)隨機(jī)裂隙網(wǎng)絡(luò)作為研究對(duì)象。

由于天然巖體在漫長(zhǎng)的地質(zhì)歷史時(shí)期中遭受多期地質(zhì)構(gòu)造運(yùn)動(dòng),裂隙網(wǎng)絡(luò)具有特殊的隨機(jī)性和不確定性。近些年發(fā)展起來的巖體裂隙網(wǎng)絡(luò)的計(jì)算機(jī)模擬是巖體結(jié)構(gòu)定量研究的最佳途徑。裂隙網(wǎng)絡(luò)的模擬過程與現(xiàn)場(chǎng)實(shí)測(cè)統(tǒng)計(jì)過程恰好相反，現(xiàn)場(chǎng)測(cè)量統(tǒng)計(jì)是根據(jù)巖體中裂隙的分布形式求出表征裂隙的各種幾何參數(shù)的分布函數(shù),而計(jì)算機(jī)模擬過程是根據(jù)統(tǒng)計(jì)得出的裂隙幾何參數(shù)的分布函數(shù)來推求服從這些分布規(guī)律的裂隙網(wǎng)絡(luò)圖形,因此,模擬是實(shí)測(cè)的逆過程。蒙特卡羅法是根據(jù)現(xiàn)場(chǎng)結(jié)構(gòu)面統(tǒng)計(jì)測(cè)量得出的參數(shù)分布來生成各參數(shù)的(偽)隨機(jī)數(shù),進(jìn)而產(chǎn)生一個(gè)與真實(shí)巖體結(jié)構(gòu)在統(tǒng)計(jì)上相似的結(jié)構(gòu)面網(wǎng)絡(luò)圖形。

巖體裂隙的基本幾何參數(shù)包括裂隙密度、裂隙空間產(chǎn)狀(裂隙走向、裂隙傾向、裂隙傾角)、裂隙間距以及裂隙跡線長(zhǎng)等。裂隙幾何特性見圖1，裂隙的產(chǎn)狀為NαE∠SEβ。

圖1 巖體裂隙幾何特性示意圖

通常約定裂隙密度服從均勻分布,裂隙的傾向和傾角服從正態(tài)分布,裂隙跡線長(zhǎng)服從正態(tài)分布、對(duì)數(shù)正態(tài)分布或者負(fù)指數(shù)分布中的一種。本文對(duì)滿足表1裂隙統(tǒng)計(jì)分布的裂隙網(wǎng)絡(luò)給出由蒙特卡羅模擬得出的二維仿真裂隙網(wǎng)絡(luò)(圖2),并對(duì)此仿真裂隙網(wǎng)絡(luò)中的脈動(dòng)壓力傳播特性進(jìn)行瞬變流數(shù)值分析。

表1 模擬巖體10 m×10 m范圍內(nèi)裂隙統(tǒng)計(jì)分布規(guī)律

注:μ為平均值,σ為方差。

圖2 仿真裂隙網(wǎng)絡(luò)圖形及其節(jié)點(diǎn)、線單元編號(hào)

由表1可見,模擬巖體10 m×10 m范圍內(nèi)有兩組不同幾何特征的裂隙,巖體頂端受射流水體沖擊,左右兩端為透水邊界,底端設(shè)為不透水邊界。

由圖2可見,裂隙網(wǎng)絡(luò)中有17個(gè)節(jié)點(diǎn)、16個(gè)線單元、0個(gè)回路,分形理論中稱之為“樹”。裂隙網(wǎng)絡(luò)生成的同時(shí),實(shí)現(xiàn)了:①去除網(wǎng)絡(luò)中的孤立裂隙。若裂隙沒有最終連通到干擾波源,即連通到控制體頂邊界,則裂隙中不可能有干擾波的存在以及傳播。這種裂隙稱為孤立裂隙或者死裂隙,在研究中不予考慮。②裂隙網(wǎng)絡(luò)中線單元及節(jié)點(diǎn)自動(dòng)編號(hào)。③裂隙網(wǎng)絡(luò)的矩陣表示:銜接矩陣——定義了裂隙節(jié)點(diǎn)與裂隙線單元之間的銜接關(guān)系和裂隙方向；幾何參量矩陣——存放各條裂隙的厚度、粗糙系數(shù)以及堵塞系數(shù)等信息。

2 余弦干擾波在裂隙網(wǎng)絡(luò)中的傳播

以仿真裂隙網(wǎng)絡(luò)作為水體賦存空間,研究其中脈動(dòng)壓力的傳播規(guī)律。沿各條有向裂隙建立一維瞬變流模型:

(1)

圖3 不同節(jié)點(diǎn)的脈動(dòng)壓力時(shí)間過程

在裂隙網(wǎng)絡(luò)節(jié)點(diǎn)1處,施加單一的余弦干擾波h=cos(2πft),其中f為干擾波頻率,其最大脈動(dòng)幅值為1.0 m水頭。水工模型實(shí)驗(yàn)表明,在水墊塘底部,導(dǎo)致底板破壞的脈動(dòng)壓力波屬于低頻波,頻帶區(qū)位于0～10 Hz之間,這里取f=5 Hz。節(jié)點(diǎn)12處未感應(yīng)到此干擾波。在此干擾源作用下,取空間步長(zhǎng)為Δx=1 m,時(shí)間步長(zhǎng)為Δt=Δx/a=0.001 s,用特征線-有限差分法離散求解各裂隙中的瞬變流雙曲線方程。

圖3分別給出節(jié)點(diǎn)1、節(jié)點(diǎn)9、節(jié)點(diǎn)5、節(jié)點(diǎn)13的脈動(dòng)壓力時(shí)間過程。圖4給出了整個(gè)裂隙網(wǎng)絡(luò)中,脈動(dòng)壓力強(qiáng)度(即方差)的分布云圖。

圖4 裂隙網(wǎng)絡(luò)內(nèi)脈動(dòng)壓力強(qiáng)度分布云圖

由圖3可以看出,當(dāng)入口施加單一余弦干擾波時(shí),干擾波將在極短時(shí)間內(nèi)傳播至巖體內(nèi)各裂隙處。其時(shí)間過程曲線與干擾源時(shí)間過程曲線(圖3(a))不同,是一條復(fù)雜的波動(dòng)曲線。文獻(xiàn)[6]對(duì)單體板塊下縫隙中脈動(dòng)壓力傳播的研究認(rèn)為，這是由于縫隙中類似于水錘的水力瞬變作用導(dǎo)致的,而且水力瞬變的頻率與縫隙長(zhǎng)度密切相關(guān)。因此,在復(fù)雜裂隙網(wǎng)絡(luò)中,不同長(zhǎng)度裂隙可以導(dǎo)致不同頻率的水力瞬變,裂隙中某點(diǎn)的脈動(dòng)壓力與其所直接銜接或者間接銜接(通過其他裂隙連通至該節(jié)點(diǎn))裂隙的水力瞬變作用相關(guān),直接銜接裂隙的水力瞬變影響尤其顯著。

由圖4容易得出:其一,裂隙網(wǎng)絡(luò)中的脈動(dòng)壓力傳播具有很強(qiáng)的各向異性特性,這主要取決于作為脈動(dòng)壓力傳播介質(zhì)的水體,其賦存的空間即巖體裂隙網(wǎng)絡(luò)具有特殊的結(jié)構(gòu)特性和不連續(xù)特性。其次,裂隙內(nèi)某點(diǎn)的脈動(dòng)強(qiáng)度與各裂隙單元的端點(diǎn)有關(guān)。裂隙單元的端點(diǎn)或?yàn)殚_口端(即緊鄰?fù)杆畬踊蛄严堆由熘量刂企w外),或?yàn)槊ざ?即緊鄰不透水層或裂隙終止于巖體內(nèi)部)。若裂隙末端為開口端,則裂隙內(nèi)各處的脈動(dòng)壓力強(qiáng)度由此處與干擾源的等效距離決定。所謂等效距離是指沿裂隙結(jié)構(gòu)面連通到干擾源的最短距離。一般而言,某點(diǎn)與干擾源的等效距離越近,則其脈動(dòng)強(qiáng)度越強(qiáng);反之,若某點(diǎn)與干擾波源的等效距離越遠(yuǎn),則其脈動(dòng)強(qiáng)度越弱。反映在圖3上,節(jié)點(diǎn)9、節(jié)點(diǎn)5、節(jié)點(diǎn)13與干擾波源(節(jié)點(diǎn)1)的等效距離依次漸遠(yuǎn),各點(diǎn)的脈動(dòng)強(qiáng)度分別為0.43 m、0.37 m、0.26 m。若裂隙末端為盲端,當(dāng)干擾波傳至裂隙末端時(shí),因存在壁面反射,使得其附近各點(diǎn)的脈動(dòng)壓力強(qiáng)度較大，且到反射源等效距離越小,反射效應(yīng)對(duì)其影響越大。反映在圖4上,裂隙單元8、16與裂隙單元9相比,雖然前者到干擾源的等效距離大,由于其所在裂隙的末端為盲端,致其壓力的脈動(dòng)強(qiáng)度遠(yuǎn)大于后者。裂隙14與控制體左邊界相交,裂隙末端為開口端,因此其壓力脈動(dòng)強(qiáng)度急劇減小。

3 結(jié) 論

a. 復(fù)雜裂隙網(wǎng)絡(luò)中的脈動(dòng)壓力傳播具有很強(qiáng)的各向異性特性,與巖體裂隙網(wǎng)絡(luò)的結(jié)構(gòu)特性密切相關(guān)。

b. 當(dāng)瞬時(shí)壓力干擾波在裂隙網(wǎng)絡(luò)中傳播時(shí),不同長(zhǎng)度裂隙可以導(dǎo)致不同頻率的水力瞬變。裂隙網(wǎng)絡(luò)中某點(diǎn)的脈動(dòng)壓力時(shí)間過程為復(fù)雜的波動(dòng)曲線,與該點(diǎn)所銜接(尤其是直接銜接)的裂隙線元的水力瞬變作用密切相關(guān)。

c. 巖體裂隙網(wǎng)絡(luò)中某點(diǎn)的壓力脈動(dòng)強(qiáng)度不僅與此點(diǎn)離干擾波源的等效距離大小有關(guān),而且與各條裂隙的端點(diǎn)條件有關(guān),壓力脈動(dòng)強(qiáng)度受裂隙入口干擾源和裂隙末端反射源共同影響。到干擾源或者反射源的等效距離越小,受其影響越大。

[1] FIOROTTO V,RINALDO A.Fluctuating uplift and lining design in spillways stilling basins[J].Journal of Hydraulic Engineering,1992,118(4): 578-596.

[2] LIU Peiqing,LI Aihua.Model discussion of pressure fluctuations propagation within lining slab joints in stilling basins[J].Journal of Hydraulic Engineering,2007,133(6): 618-624.

[3] LIU Peiqing,LI Aihua.Fluctuating uplift acting on rock blocks at the bottom of river bed and estimation of the limiting scour depth[J].Journal of Hydraulic Research,2007,45(4): 478-485.

[4] LI Aihua,LIU Peiqing.Mechanism of rock-bed scour due to impinging jet[J].Journal of Hydraulic Research,2010,48(1): 14-22.

[5] 李愛華,劉沛清.脈動(dòng)壓力在消力池底板縫隙傳播的瞬變流模型和滲流模型統(tǒng)一性探討[J].水利學(xué)報(bào),2005,36(10): 1236-1240.(LI Aihua,LIU Peiqing.Models for simulating the propagation of fluctuating pressure in gaps of stilling basin slab[J].Journal of Hydraulic Engineering,2005,36(10): 1236-1240.(in Chinese))

[6] 李愛華,王騰,劉沛清.溪洛渡壩區(qū)巖石河床沖刷過程數(shù)值模擬[J].水力發(fā)電學(xué)報(bào),2012,31(5): 154-158.(LI Aihua,WANG Teng,LIU Peiqing.Numerical simulation of rock bed scour behind the dam of Xiluodu Hydropower Station[J].Journal of Hydroelectric Engineering,2012,31(5): 154-158.(in Chinese))

[7] 李愛華,劉沛清.脈動(dòng)壓力在板塊縫隙中傳播衰變機(jī)理研究[J].水利水電技術(shù),2006,37(9): 33-37.(LI Aihua,LIU Peiqing.Study on mechanism of pulsating pressure propagation disintegration in joints of bottom slab[J].Water Resources and Hydropower Engineering,2006,37(9): 33-37.(in Chinese))

[8] JIA Y,KITAMURA T,WANG S S Y.Simulation of scour process in plunging pool of loose bed material[J].Journal of Hydraulic Engineering,2001,127(3): 219-229.

[9] 劉安富,楊敏,李會(huì)平,等.寬尾墩消力池透水底板脈動(dòng)壓力試驗(yàn)研究[J].水科學(xué)進(jìn)展,2012,23(2): 243-248.(LIU Anfu,YANG Min,LI Huiping,et al.Experimental study of the fluctuating pressure on pervious slab in stilling basin under jet flow from tail-flaring piers[J].Advances in Water Science,2012,23(2): 243-248.(in Chinese))

[10] 辜晉德,彭秀芳.消力池底板縫隙水流脈動(dòng)壓力頻譜分析研究[J].水力發(fā)電學(xué)報(bào),2013,32(6): 177-182.(GU Jinde,PENG Xiufang.Study on spectrum of fluctuating pressure under bottom slabs of stilling basins[J].Journal of Hydroelectric Engineering,2013,32(6): 177-182.(in Chinese))

[11] ANNANDALE G W.Erodbility[J].Journal of Hydraulic Research,1995,33(4): 471-494.

[12] 潘別桐,井蘭如.巖體結(jié)構(gòu)概率模擬和應(yīng)用[M].沈陽: 東北工學(xué)院出版社,1989: 1-79.

[13] 陳劍平.巖體隨機(jī)不連續(xù)面三維網(wǎng)絡(luò)數(shù)值模擬技術(shù)[J].巖土工程學(xué)報(bào),2001,23(4): 397-402.(CHEN Jianping.3D network numerical modeling technique for random discontinuities of rock mass[J].Chinese Journal of Geotechnical Engineering,2001,23(4): 397-402.(in Chinese))

[14] 劉曉麗,王恩志,王思敬,等.裂隙巖體表征方法及巖體水力學(xué)特性研究[J].巖石力學(xué)與工程學(xué)報(bào),2008,27(9): 1814-1821.(LIU Xiaoli,WANG Enzhi,WANG Sijing,et al.Representation method of fractured rock mass and its hydraulic properties study[J].Chinese Journal of Rock Mechanics and Engineering,2008,27(9): 1814-1821.(in Chinese))

[15] 劉曉麗, 王恩志, 王思敬.裂隙巖體精細(xì)結(jié)構(gòu)描述及工程特性數(shù)值試驗(yàn)[J].巖石力學(xué)與工程學(xué)報(bào), 2008, 27(增刊2): 3935-3940.(LIU Xiaoli, WANG Enzhi, WANG Sijing.Fine structure description for fractured rock mass and numerical tests on its engineering properties [J].Chinese Journal of Rock Mechanics and Engineering, 2008, 27(Sup2): 3935-3940.(in Chinese))[16] 束龍倉，張春艷， 閔星，等.不同隙寬條件下L形裂隙水頭損失試驗(yàn)[J].河海大學(xué)學(xué)報(bào)：自然科學(xué)版,2013,41(5):383-388.(SHU Longcang,ZHANG Chunyan,MIN Xing,et al.Experimental study on water head loss of L-shaped fractures with different aperture widths[J].Journal of Hohai University：Natural Sciences,2013,41(5):383-388. (in Chinese))

Numerical simulation on propagation of fluctuating pressure of flow within rock fractures//

LI Aihua, ZHU Jiang
(SchoolofPetroleumEngineering,ChinaUniversityofPetroleum,Qingdao266580,China)

We used the Monte-Carlo method to analyze the influence of network structure on the propagation law of fluctuating pressure of flow within fractures. A replica of the fracture network was created by numerical simulation, which has resembled statistical characteristics with natural fracture network. Taking the replica of fracture network as main space of fissure water, a one-dimensional transient flow model was established to numerically simulate the propagation of fluctuating pressure. And then, the unsteady time process of fluctuating pressure’s propagation of fissure water and the spatial distribution nephogram of fluctuating intensity were generated. We concluded that the propagation of fluctuating pressure in rock fracture network has significant anisotropy. And the transient fluctuating pressure at any node within fractures is determined by the combination of three factors: the hydro-transient effect of fracture line elements connected with the node(especially of those directly connected with the node), the interfering source at fracture entrances, and the reflection source at fracture ends.

fluctuating pressure; fracture network; Monte-Carlo method; transient flow model; numerical simulation

中央高?；究蒲袠I(yè)務(wù)費(fèi)專項(xiàng)(12CX04024A,13CX02049A)

李愛華(1978—),女,山西祁縣人,副教授,博士,主要從事流體力學(xué)及計(jì)算流體力學(xué)研究。E-mail: aiwali0524@163.com

10.3880/j.issn.1006-7647.2015.04.006

P641.2

A

1006-7647(2015)04-0024-04

2014-04-18 編輯：鄭孝宇)