錦屏一級(jí)左岸巖體滲透特性的數(shù)值模擬研究

田洪圓，張 亮，李 宏

(1.大連理工大學(xué) 海岸和近海工程國(guó)家重點(diǎn)實(shí)驗(yàn)室， 遼寧 大連 116024；2.中國(guó)電建集團(tuán) 昆明勘測(cè)設(shè)計(jì)研究院有限公司， 云南 昆明 650041)

錦屏一級(jí)左岸巖體滲透特性的數(shù)值模擬研究

田洪圓1，張 亮2，李 宏1

(1.大連理工大學(xué) 海岸和近海工程國(guó)家重點(diǎn)實(shí)驗(yàn)室， 遼寧 大連 116024；2.中國(guó)電建集團(tuán) 昆明勘測(cè)設(shè)計(jì)研究院有限公司， 云南 昆明 650041)

基于RFPA2D-Flow軟件，并結(jié)合現(xiàn)場(chǎng)壓水試驗(yàn)，從表征單元體積(REV)及滲透張量入手，對(duì)裂隙巖體的滲透特性進(jìn)行數(shù)值模擬研究。以錦屏一級(jí)水電站左岸巖體為例，基于巖體裂隙幾何特征，運(yùn)用蒙特卡羅方法生成巖體裂隙網(wǎng)絡(luò)模型；模擬計(jì)算不同尺寸不同方向上的等效滲透系數(shù)，依據(jù)類常量與類張量特性初步確定裂隙巖體滲透表征單元體積(REV)及滲透張量，再結(jié)合現(xiàn)場(chǎng)壓水試驗(yàn)結(jié)果修正滲透張量。隨著模型尺寸的增大，巖體的滲透性逐漸趨于穩(wěn)定，表征單元體尺寸的大小約為跡長(zhǎng)的4倍，表明基于RFPA2D-Flow數(shù)值模擬并結(jié)合壓水試驗(yàn)確定裂隙巖體滲透張量和REV的方法合理有效。

裂隙網(wǎng)絡(luò)模型；表征單元體積；滲透張量；蒙特卡羅方法；RFPA2D-Flow

裂隙巖體中不同產(chǎn)狀、跡長(zhǎng)、隙寬等幾何特征，使巖體滲透特性在空間上存在客觀的各向異性。大型巖體工程的滲流分析中等效連續(xù)介質(zhì)模型應(yīng)用最為廣泛[1]，雖然通過等效連續(xù)介質(zhì)模型并不能對(duì)裂隙局部滲流進(jìn)行詳細(xì)描述，但它能夠闡述清楚裂隙巖體宏觀滲流特點(diǎn)，且計(jì)算量相對(duì)較小。

Sagar Budhi等[2]考慮了裂隙幾何特征以及裂隙幾何參數(shù)引起的滲透系數(shù)的誤差，提出二階統(tǒng)計(jì)分析方法；Oda Masanobu[3-4]基于裂隙空間分布的隨機(jī)性，通過統(tǒng)計(jì)理論對(duì)滲透張量進(jìn)行推導(dǎo)，同時(shí)針對(duì)荷載條件作用下巖體滲透性進(jìn)行了研究；Zhang X等[5]采用多孔介質(zhì)模型和離散裂隙網(wǎng)格模型，考慮裂隙幾何特征和應(yīng)力，數(shù)值模擬得到二維裂隙巖體滲透張量；Chen S H等[6]的復(fù)合單元法結(jié)合蒙特卡羅方法生成裂隙網(wǎng)絡(luò)，計(jì)算不同尺寸巖體試件不同方向上的等效滲透系數(shù)過程中考慮了巖塊與裂隙間的相互滲透作用；Rong Guan等[7]基于節(jié)理網(wǎng)絡(luò)模擬技術(shù)，應(yīng)用滲流能量疊加原理，推導(dǎo)得到計(jì)算節(jié)理巖體滲透張量公式，提出明確裂隙巖體滲透表征單元體積之法，并討論了裂隙幾何特征對(duì)巖體滲透特性的影響；王培濤等[8]通過離散裂隙網(wǎng)絡(luò)模型，構(gòu)建平面滲流分析法，對(duì)節(jié)里巖體中的滲透張量進(jìn)行了探討，利用對(duì)滲流定向性系數(shù)進(jìn)行定義，分析巖體滲流定向特征；楊金保等[9]基于裂隙幾何特征的基礎(chǔ)上，還考慮了應(yīng)力作用對(duì)滲透系數(shù)張量的影響；車富強(qiáng)等[10]把裂隙樣本法和現(xiàn)場(chǎng)測(cè)試法相結(jié)合，用裂隙樣本法初步確定裂隙巖體滲透張量，采用現(xiàn)場(chǎng)測(cè)試法修正裂隙巖體滲透張量。王俊奇等[11]利用離散裂隙網(wǎng)絡(luò)的管單元模型研究了裂隙巖體滲透張量，通過滲透橢圓擬合確定。

本文以錦屏一級(jí)水電站左岸巖體為例，根據(jù)巖體地質(zhì)勘探資料，生成隨機(jī)裂隙網(wǎng)絡(luò)模型，將模型數(shù)據(jù)導(dǎo)入RFPA2D-Flow軟件；選取巖體不同尺寸、不同方向的分析域進(jìn)行滲流數(shù)值模擬，得到等效滲透系數(shù)；基于類常量與類張量特性[12]來確定裂隙巖體滲流表征單元體積和滲透張量；再結(jié)合現(xiàn)場(chǎng)原位壓水試驗(yàn)數(shù)據(jù)，對(duì)模擬結(jié)果進(jìn)行修正。

1 裂隙巖體滲透張量及表征單元體積分析的基本原理

確定裂隙巖體滲透REV、滲透張量應(yīng)用的方法是：按照現(xiàn)場(chǎng)實(shí)測(cè)區(qū)域中內(nèi)節(jié)理幾何參數(shù)，通過Monte-Carlo法構(gòu)建出裂隙網(wǎng)絡(luò)生成域，然后于該域里的中心位置取出一系列尺度不一的正方形作為REV試算域。

一般情況下，REV的尺寸為裂隙平均跡長(zhǎng)的3～4倍[13]，而分析域的尺寸必須大于REV尺寸，才能使計(jì)算有意義。所以根據(jù)裂隙跡長(zhǎng)的統(tǒng)計(jì)數(shù)據(jù)大致確定出最大分析域的尺寸，按大于分析域的2倍大小來確定裂隙網(wǎng)絡(luò)模型生成域的大小，然后于該域里的中心位置取出一系列尺度不一、不同方向的正方形作為REV試算域。

1.1 裂隙巖體的平面滲流模型

二維裂隙網(wǎng)絡(luò)模型里，為了便于操作，設(shè)定節(jié)理是短直線相互交錯(cuò)所形成。根據(jù)現(xiàn)場(chǎng)試驗(yàn)資料來看，節(jié)理傾角以及間距還有跡長(zhǎng)和密度等幾何特征做了研究分析，獲得了某種概率密度分布，然后再通過Monte-Carlo法生成節(jié)理網(wǎng)絡(luò)，可參考圖1。

圖1 REV分析模型

根據(jù)圖1，定水頭邊界是MM′、NN′，而MM′、NN′間滲透系數(shù)如式(1)，不透水邊界是MN、M′N′。

(1)

在這個(gè)公式里：Δq表示研究區(qū)域中產(chǎn)生的總流量，m3/s；ΔH為M′M、NN′二者間水頭差，m；K則是MN方向等效滲透系數(shù)，m/s；MN、M′M是模型邊長(zhǎng)，m，二者相等長(zhǎng)。

1.2 滲透張量的計(jì)算

滲透張量與其他的二階張量一樣，也具有坐標(biāo)轉(zhuǎn)換特性。圖2中假定整體坐標(biāo)系xoy以原點(diǎn)o為圓心，逆時(shí)針旋轉(zhuǎn)角度形成局部坐標(biāo)系x′oy′，x′軸與x,y軸的方向余弦記為lx′x=cosθ，lx′y=sinθ，y′與x,y軸的方向余弦記為ly′x=-sinθ，ly′y=cosθ。

則[K′]各元素Ki,j與[K]中Kij間滿足轉(zhuǎn)換關(guān)系：

Ki′j′=Kij·li′i·lj′j

(2)

(3)

整體坐標(biāo)系下x軸和y軸分別是主滲透系數(shù)的方向，則局部坐標(biāo)系下軸方向的滲透系數(shù)為：

Kθ=K1·cos2θ+K2·sin2θ

(4)

極坐標(biāo)系(l,θ)下，橢圓方程可以表示為：

(5)

圖2 坐標(biāo)系相對(duì)位置

圖3 橢圓關(guān)系曲線

1.3 滲透張量的修正

上述1.1節(jié)計(jì)算得到的滲透系數(shù)并不能代表某一工程巖體的平均滲透性，此處取幾何平均值作為巖體等效滲透系數(shù)：

(6)

式中：Kω為等效滲透系數(shù)；Ki為計(jì)算得到的滲透主值。取m=K0/Kω，對(duì)滲透張量進(jìn)行修正，其中K0為現(xiàn)場(chǎng)壓水試驗(yàn)得到的滲透系數(shù)。

1.4 表征單元體積的確定

裂隙巖體滲透REV即可以在宏觀方面等效表現(xiàn)出巖體滲透性相對(duì)應(yīng)的體積最小值，如果巖體尺寸比REV更大時(shí)，表明有二階對(duì)稱正定滲透張量。

把裂隙巖體等效為連續(xù)介質(zhì)，在此引入類張量和類常量特性[11]。當(dāng)所研究的巖體的滲透性參數(shù)在給定的誤差范圍內(nèi)同時(shí)滿足類張量和類常量特性，此時(shí)研究區(qū)域大小即為相應(yīng)的滲透REV。

在裂隙網(wǎng)格生成域內(nèi)以模型中心為基點(diǎn)，以2 m起始至28 m構(gòu)成14個(gè)研究域，每個(gè)對(duì)應(yīng)尺寸每隔30°旋轉(zhuǎn)一次，構(gòu)成12個(gè)方向，通過統(tǒng)計(jì)下側(cè)邊界流量求得裂隙巖體等效平均滲透系數(shù)。此處利用裂隙網(wǎng)絡(luò)對(duì)稱性，僅對(duì)0°～150°上6個(gè)方向進(jìn)行模擬計(jì)算，根據(jù)對(duì)稱性得到其余6個(gè)方向上等效滲透系數(shù)，確定滲透特性表征單元體REV大小。

2 工程實(shí)例

國(guó)內(nèi)現(xiàn)在已建成的拱壩中，最高的即錦屏一級(jí)水電站，該水電站的水庫(kù)正常蓄水位1 880 m，庫(kù)容77.6億m3。巖層的傾向在285°～295°范圍內(nèi)，與之相應(yīng)的傾角則在35°～45°范圍內(nèi)。拱壩左邊是反向坡，右邊是順向坡，整個(gè)區(qū)域中巖體結(jié)構(gòu)面明顯，同時(shí)受地質(zhì)構(gòu)造作用的影響強(qiáng)大[14]，壩址區(qū)節(jié)理裂隙劃分成5組：

(1) N15°—80°E，NW∠25°～45°；

(2) N50°—70°E，SE∠50°～80°；

(3) 近SN—N30°E，SE∠60°～80°；

(4) N30°—50°W，NE∠60°～80°；

(5) N60°W—EW，NE(SW)或S(N)∠60°～80°。

各級(jí)巖體巖體結(jié)構(gòu)、風(fēng)化卸荷以及所處岸坡不同部位情況，各級(jí)巖體的滲透系數(shù)如下：

Ⅱ級(jí)巖體：q<3 Lu，新鮮巖體；

Ⅲ1級(jí)巖體：3≤q<10 Lu，深部微卸荷巖體；

Ⅲ2級(jí)巖體：10≤q<100 Lu，淺表弱卸荷巖體；

Ⅳ1級(jí)巖體：q≥100 Lu，岸坡淺表部強(qiáng)卸荷巖體，卸荷裂隙發(fā)育且多張開。

Ⅳ2級(jí)巖體：q≥100 Lu，左岸深部裂縫發(fā)育段巖體，多松弛破碎張開呈空縫。

根據(jù)表1中裂隙的分布規(guī)律，可以發(fā)現(xiàn)跡長(zhǎng)均值為2.29 m，結(jié)合參考文獻(xiàn)[15]，可取分析域尺寸為10 m×10 m，則生成域應(yīng)大于分析域的2倍，此處取為25 m×25 m，由于幾何參數(shù)中缺乏對(duì)裂隙張開度的統(tǒng)計(jì)數(shù)據(jù)，所以此處參考文獻(xiàn)[2]，取裂隙張開度的估計(jì)值微0.2 mm?；诿商乜_方法生成裂隙網(wǎng)絡(luò)模型，得到裂隙的數(shù)據(jù)文件，然后將裂隙幾何信息導(dǎo)入RFPA2D-Flow軟件，邊界條件如圖1所示，通過數(shù)值模擬結(jié)算得到裂隙巖體的滲透特性變化如圖4所示。

表1 巖體結(jié)構(gòu)面幾何參數(shù)

圖4 錦屏一級(jí)左岸巖體不同尺寸滲透特性變化

給定誤差限15%，根據(jù)圖4，由類常量特性圖4(a)可見，當(dāng)尺寸達(dá)到11 m×11 m以上，Kxx、Kyy、Kxy滿足類常量特性；由類張量特性圖4(b)可見，當(dāng)尺寸達(dá)到10 m×10 m，變異系數(shù)在允許誤差范圍內(nèi)；再通過圖4(c)可以進(jìn)一步驗(yàn)證，當(dāng)尺寸達(dá)到12 m×12 m以后，滲透主值以及兩比值系數(shù)亦趨于穩(wěn)定，所以取REV尺寸為12 m×12 m，相當(dāng)于最大裂隙跡長(zhǎng)的4倍，與文獻(xiàn)[15]的結(jié)論相符合。

上文裂隙張開度取為估計(jì)值來確定REV大小無可厚非，但以此得到的滲透張量并不能用來表述整個(gè)巖體的滲透特性，因?yàn)閹r體各巖層、各分區(qū)的滲透性亦不同，為了更好的反應(yīng)巖體滲透特性，結(jié)合現(xiàn)場(chǎng)試驗(yàn)數(shù)據(jù)，將各巖層滲透系數(shù)相近者歸為一類，取其平均滲透系數(shù)近似推算出裂隙開度均值，取分析域尺寸為15 m×15 m，通過數(shù)值模擬分別計(jì)算各研究域滲透特性，得到初始滲透張量如表2所示。

表2 初步確定滲透張量

由于隨機(jī)裂隙網(wǎng)絡(luò)具有隨機(jī)性，并不能完全反應(yīng)出真實(shí)的裂隙分布以及裂隙連通情況，所以再次根據(jù)現(xiàn)場(chǎng)壓水試驗(yàn)數(shù)據(jù)(見表3)進(jìn)行修正，修正后的滲透張量如表4所示。

表3 各巖層水壓試驗(yàn)水力傳導(dǎo)系數(shù)[16]單位：m/s

表4 修正后各分區(qū)的滲透張量

3 結(jié) 論

(1) 結(jié)合數(shù)值試驗(yàn)法和現(xiàn)場(chǎng)單孔壓水試驗(yàn)，本文以錦屏一級(jí)水電站壩區(qū)左岸巖體為例，對(duì)滲透張量及滲流REV進(jìn)行了研究。

(2) 在研究域內(nèi)，根據(jù)巖體內(nèi)節(jié)理的分布特征，運(yùn)用蒙特卡羅方法生成隨機(jī)裂隙網(wǎng)絡(luò)，再導(dǎo)入RFPA2D-Flow進(jìn)行滲流計(jì)算，能同時(shí)考慮節(jié)理的分布及連通情況，此法特別適用于裂隙巖體二維滲透張量的研究。由于節(jié)理裂隙分布的隨機(jī)性，數(shù)值試驗(yàn)中采用的裂隙產(chǎn)狀及幾何參數(shù)也不能完全的反映巖體內(nèi)節(jié)理的實(shí)際分布，由此計(jì)算得到的滲透張量存在理論誤差。特別地，通過現(xiàn)場(chǎng)壓水試驗(yàn)結(jié)果來修正模擬結(jié)果，保證數(shù)值和工程要求相一致。

(3) 該研究通過滲透張量計(jì)算方法，考慮了節(jié)理的實(shí)際分布及連通情況、不同巖層及巖體等級(jí)對(duì)滲透張量的影響。該法具有計(jì)算方便、易于理解的優(yōu)勢(shì)，非常適合用于研究裂隙巖體滲透性。

[1] Ansari M A, Deodhar A, Kumar U S, et al. Isotope Hydrogeological Factors Control Transport of Radon-222 in Hard Rock Fractured Aquifer of Bangalore,Karnataka[M]. Berlin: Springer International Publishing, 2016.

[2] Sagar Budhi, Runchal Akshai. Permeability of fractured rock: effect of fracture size and data uncertainties[J]. Water Res. Research, 1982,18(2):266-274.

[3] Oda Masanobu. Permeability tensor for discontinuous rockmass[J]. Geotechnique,1985,35(4):483-495.

[4] Oda Masanobu. An equivalent continuum model for coupled stress and fluid flow analysis in jointed rock masses[J]. Water Resources Research, 1986,22(13):1845-1856.

[5] Zhang X, Sanderson D J, Harkness R M, et al. Evaluation of the 2-D permeability tensor for fractured rock masses[J]. International Journal of Rock Mechanics and Mining Science & Geomechanics Abstracts, 1996,33(1):17-37.

[6] Chen S H, Feng X M, Isam S. Numerical estimation of REV and permeability tensor for fractured rock masses by composite element method[J]. International Journal for Numerical and Analytical Methods in Geomechanics, 2008,32(12):1459-1477.

[7] Rong Guan, Peng Jun, Wang X J, et al. Permeability tensor and representative elementary volume of fractured rock masses[J]. Hydrogeology Journal, 2013,21(7):1655-1671.

[8] 王培濤,楊天鴻，于慶磊，等.基于離散裂隙網(wǎng)絡(luò)模型的節(jié)理巖體滲透張量及特性分析[J].巖土力學(xué),2013,34(S2):448-455.

[9] 楊金保,李 蕾,謝 妮.考慮應(yīng)力狀態(tài)的裂隙巖體滲透系數(shù)張量計(jì)算[J].人民黃河,2013,35(3):131-134.

[10] 車富強(qiáng),陳益峰，等.基于等效連續(xù)介質(zhì)壩區(qū)裂隙巖體滲控效應(yīng)分析[J].中南大學(xué)學(xué)報(bào)(自然科學(xué)版),2014,45(6):1967-1974.

[11] 王俊奇，董 曄，李鵬宇，等.三維離散裂隙網(wǎng)絡(luò)管單元模型確定巖體滲透張量的嘗試[J].水利與建筑工程學(xué)報(bào),2016,14(3):47-54.

[12] 張宜虎.巖體等效水力學(xué)參數(shù)研究[D].武漢：中國(guó)地質(zhì)大學(xué)(武漢)，2006.

[13] 張貴科,徐衛(wèi)亞.裂隙網(wǎng)絡(luò)模擬和REV尺度研究[J].巖土力學(xué),2008,29(6):1675-1670.

[14] 李楊鵬,薛新華，楊興國(guó)，等.錦屏一級(jí)水電站降雨入滲條件下左岸邊坡穩(wěn)定性分析[J].水利與建筑工程學(xué)報(bào),2014,12(1):16-20.

[15] 安玉華,王 清.基于三維裂隙網(wǎng)絡(luò)的裂隙巖體表征單元體研究[J].巖土工程,2012,33(12):3775-3780.

[16] 荔 強(qiáng).錦屏一級(jí)左岸邊坡霧化滲流穩(wěn)定性問題研究[D].北京：清華大學(xué)，2010.

Numerical Simulation of Permeability of Rock Mass on the Left Bank of Jinping-I

TIAN Hongyuan1, ZHANG Liang2, LI Hong1

(1.StateKeyLaboratoryofCoastalandOffshoreEngineering,DalianUniversityofTechnology,Dalian,Liaoning116024,China;2.PowerChinaKunmingEngineeringCorporationLimited,Kunming,Yunnan650041,China)

Based on the numerical simulation and finite element analysis software RFPA2D-Flow and single-hole packer tests, the two-dimensional representative element volume and permeability tensor was analyzed for fractured rock mass. According to the geometric characteristics of the joints, such as the length, the density and the yield of the joints. The random fracture networks of Jinping I hydropower station are generated by using Monte Carlo method, numerical simulation of RFPA2D-Flow is adopted to calculate the equivalent permeability coefficient of different size in different directions of fractured rock masses, in order to obtain the representative elementary volume and permeability tensor. Then they are further modified by the single-hole packer test results. The results show that as the model size increases, the permeability of rock gradually begins to stabilize, the size of REV is about 4 times the length of trace. The method is reasonable and effective to determine the permeability tensor and REV of fractured rock masses combining single-hole packer tests and numerical simulation of RFPA2D-Flow.

fracture network model; representative element volume; permeability tensor; Monte Carlo method; RFPA2D-Flow

10.3969/j.issn.1672-1144.2017.03.013

2017-02-14

2017-03-10

國(guó)家重點(diǎn)基礎(chǔ)研究“973”計(jì)劃項(xiàng)目(2011CB013503)

李 宏(1969—)，男，遼寧遼陽人，教授，主要從事儲(chǔ)層流動(dòng)數(shù)值模擬、邊坡節(jié)理巖體滲流研究。E-mail:hong.li@dlut.edu.cn

TV139.14

A

1672—1144(2017)03—0067—05