杜子建 高永濤 鄧代強(qiáng) 韓浩亮

(北京科技大學(xué)土木與環(huán)境工程學(xué)院,北京 100083)

目前,國內(nèi)外學(xué)者對天然巖體的水力學(xué)特性研究較多,而對于擾動(dòng)后的巖體滲流特性卻研究較少,例如,對于目前比較常見的地下采空區(qū),礦層開采之后破壞了圍巖的應(yīng)力狀態(tài),造成應(yīng)力的重新分布,導(dǎo)致圍巖變形、破壞,改變了天然巖體的裂隙分布,巖體的滲透性也會隨之改變,采空區(qū)范圍的滲流場與天然巖體必然有較大區(qū)別[1,2].特別當(dāng)采空區(qū)位于水體滲流區(qū)域,如本文將提到的水庫下伏采空區(qū),水庫水體會沿著新產(chǎn)生的采動(dòng)裂隙流動(dòng),并在采空區(qū)裂隙帶形成水力通道,水體對裂隙體滲透作用,以及對裂隙縫中破碎體以及充填物長期沖蝕,會造成裂隙體的強(qiáng)度進(jìn)一步削弱,引起巖體應(yīng)力場的環(huán)境變遷,導(dǎo)致裂隙巖體的滲透變形,增大巖體的滲透性,在兩者長期相互作用下,如不采取措施必然會導(dǎo)致較嚴(yán)重的后果.

1 擾動(dòng)巖體應(yīng)力場對滲流場影響

擾動(dòng)巖體重新分布的應(yīng)力場改變了裂隙的位移及結(jié)構(gòu)面的隙寬,進(jìn)而影響地下巖體整體滲透特性,應(yīng)力場對滲流特性的影響主要表現(xiàn)為巖體滲流場滲透系數(shù)k的改變[3].

假設(shè)裂隙巖體介質(zhì)統(tǒng)計(jì)為規(guī)則和均質(zhì)的單純裂隙巖體介質(zhì),對于布有N個(gè)方向裂隙組成的導(dǎo)水裂隙網(wǎng)絡(luò),根據(jù)滲透性的可疊加原理,將每組滲透裂隙面的滲透張量進(jìn)行疊加,則有幾何參數(shù)表達(dá)的總滲透張量公式[4-6],裂隙巖體滲透張量K的數(shù)學(xué)模型為

式中,ei為第i組裂隙的等效水力隙寬;λi為第i組裂隙的平均線密度;nxi、nyi、nzi為第i組裂隙面法向的方向余弦.

在三向應(yīng)力作用下,垂直于裂隙面方向的正應(yīng)力對滲流特性起著至關(guān)重要的作用[7-8],Snow(1968)、Louis(1974)、Kilsall(1984)等通過試驗(yàn)研究,相繼得出了裂隙巖體滲透系數(shù)k隨正應(yīng)力的經(jīng)驗(yàn)關(guān)系式.以此為基礎(chǔ),結(jié)合大量勘探經(jīng)驗(yàn)成果,得出擾動(dòng)巖體應(yīng)力場改變環(huán)境下裂隙巖體滲透系數(shù)k隨正應(yīng)力σ變化的經(jīng)驗(yàn)關(guān)系式:

(1)當(dāng)正應(yīng)力σ變化范圍比較小時(shí),經(jīng)驗(yàn)關(guān)系式為直線衰減形式

式中,k′為σ=σ0時(shí)裂隙巖體的滲透系數(shù),λ為待定系數(shù),由試驗(yàn)確定.

(2)當(dāng)正應(yīng)力σ變化范圍比較大時(shí),經(jīng)驗(yàn)關(guān)系式為負(fù)指數(shù)衰減形式

式中,k′為σ=σ0時(shí)裂隙巖體的滲透系數(shù),λ為待定系數(shù),由試驗(yàn)確定.

2 擾動(dòng)巖體滲流場對應(yīng)力場影響

目前水荷載的計(jì)算方法忽略了滲流場的影響,而以靜水壓力的形式表示水荷載.實(shí)際上,在任何透水介質(zhì)中,水荷載應(yīng)以滲透體積力和滲透壓力(面力)表示,是和透水介質(zhì)中滲流場的分布密切相關(guān)的,一定的滲流場分布對應(yīng)一定的水荷載分布;而滲流場分布的變化也將引起水荷載分布的變化.所以,滲流場對應(yīng)力場的影響是通過改變巖體外荷載(水荷載)而改變巖體的應(yīng)力場分布的[9].

若已知某時(shí)刻覆蓋層內(nèi)的水頭分布函數(shù)h(x,z),則可求出某作用面上的滲透水壓力u分布為

進(jìn)而可求出覆蓋層內(nèi)的單位土體所受的滲透體積力f分布為

式中,fx、fz為滲透體積力f在x、z方向上的分力.采用有限元數(shù)值方法計(jì)算時(shí),將上述求得的滲透水壓力和體積力轉(zhuǎn)化為單元的結(jié)點(diǎn)外荷載而加以考慮.

3 耦合分析方法及主要步驟

滲流場與應(yīng)力場耦合模型的建立,關(guān)鍵是兩場耦合關(guān)系的正確建立,耦合的非線性在于滲透系數(shù)是孔隙率的非線性函數(shù).巖體應(yīng)力場對滲流場的影響在于應(yīng)力場改變了介質(zhì)本身的結(jié)構(gòu),介質(zhì)滲透性能發(fā)生改變;滲流場對應(yīng)力場的影響在于介質(zhì)中單元體所受的滲透力和外部水荷載在不斷改變,引起巖體介質(zhì)中應(yīng)力的不斷調(diào)整.

運(yùn)用物理環(huán)境耦合方法進(jìn)行分析,均質(zhì)巖體滲流場與應(yīng)力場耦合分析的有限元數(shù)值解法及迭代步驟[10]如下:(1)首先假設(shè)應(yīng)力場為0,求出均質(zhì)巖體的滲透系數(shù),讀取流體物理環(huán)境,運(yùn)行FLOT RAN流體模塊求解滲流場,得出水頭分布h;(2)將水頭分布h代入式(4)、(5),求解出滲透水壓力分布u,滲透體積力f;(3)讀取結(jié)構(gòu)物理環(huán)境,加載滲透水壓力u及滲透體積力f,按有限元法求解應(yīng)力場,得出應(yīng)力場的分布以及位移場分布;(4)將應(yīng)力場載入模型,代入式(2)或式(3),可得出擾動(dòng)巖體應(yīng)力場改變環(huán)境下裂隙巖體滲透系數(shù)k′分布;(5)回到第1步,重復(fù)1～4步,迭代直至滿足精度要求,最后得出耦合條件下滲流場的水頭壓力和孔隙水壓力.

ANSYS實(shí)現(xiàn)物理環(huán)境耦合模擬分析的方法流程圖如圖1所示.

圖1 物理環(huán)境法耦合分析流程

4 模型建立

4.1 模型的假設(shè)條件

①巖體是各向同性彈塑性體,屈服滿足Drucker-Prager強(qiáng)度準(zhǔn)則,為等效連續(xù)介質(zhì)[11];②孔隙水滲流服從Darcy定律;③巖體的滲透系數(shù)k是動(dòng)態(tài)值,滲透系數(shù)與孔隙率呈一定函數(shù)關(guān)系[12].

4.2 數(shù)值計(jì)算模型尺寸設(shè)計(jì)

計(jì)算模型取地表面為模型頂面,忽略表層少量土體層的影響,以基巖作為計(jì)算模型表面.由于采空區(qū)為地表50m以下,離水庫水平距離為30m,距水庫底垂直距離為30m左右,為了保證計(jì)算邊界大小,并充分考慮水庫滲流場的影響,本次建模取模型尺寸為800m×800m×320m,模型實(shí)地所在庫區(qū)范圍如圖2所示.將原地面取為模型頂面,采空區(qū)平均埋深為57.5m;水庫深度為20m,長寬取300m×140m,模型及水庫三維網(wǎng)格示意圖如圖3所示.

4.3 模型巖體參數(shù)確定

模擬計(jì)算采用的巖體物理力學(xué)參數(shù)見表1.

表1 巖體物理力學(xué)參數(shù)表

根據(jù)現(xiàn)場地質(zhì)調(diào)查和相關(guān)研究提供的巖石力學(xué)及滲流試驗(yàn)結(jié)果,考慮了尺寸效應(yīng)和地層構(gòu)造面的影響,對實(shí)驗(yàn)得出的巖石各項(xiàng)參數(shù)進(jìn)行相應(yīng)的調(diào)整和簡化.根據(jù)采空區(qū)現(xiàn)場觀測,巖體裂隙開口度小于0.001 m,每米平均裂隙密度為6～8條,裂隙平均間距為0.12m,采空區(qū)平均日涌水量為920m3,由公式(1)得到采空區(qū)等效滲透張量K,從而得到巖體彈性模型折減系數(shù)[13]為0.4,即數(shù)值模擬中巖體彈性模量取巖石彈性模量的0.4倍,采用巖體的破壞為塑形破壞(Drucker-Prager屈服準(zhǔn)則)和拉破壞.

4.4 模型邊界條件的確定

模型區(qū)域地應(yīng)力場最大主應(yīng)力方向?yàn)镾82°E,在計(jì)算域邊界上施加?xùn)|西向應(yīng)力,即在模型邊界上施加邊界法向水平構(gòu)造應(yīng)力,應(yīng)力邊界條件及其外加載模型如圖4所示.

圖4 計(jì)算模型應(yīng)力邊界條件

滲流場滿足水力邊界條件

式中,h(M,t)為時(shí)刻t點(diǎn)M的測壓管水頭值;φs(M,t)為邊界上給定的已知測壓管水頭值;M為邊界S1上的點(diǎn).

在水庫底部及側(cè)邊加載水力邊界條件,按水庫水體高度加載水頭,以水頭壓力直接加載.

5 下伏巖體滲流模擬結(jié)果分析

取采空區(qū)X-560剖面水頭進(jìn)行分析,從圖5可以看出,采空區(qū)滲流達(dá)到穩(wěn)定時(shí),在采空區(qū)區(qū)域整體上形成以采空區(qū)為中心,水頭由中心向四周擴(kuò)散逐漸變大的滲流區(qū),其中以采空區(qū)附近水頭最低,區(qū)域水頭最大值僅1.3m,水頭基本以采空區(qū)為中心呈梯度變化并趨向穩(wěn)定,距離采空區(qū)近,水頭變化較大,距離采空區(qū)越遠(yuǎn),水頭變化就越小,水力梯度也逐漸減小.

圖5 空區(qū)水頭壓力X-560剖面分析圖

取采空區(qū)X-550剖面孔隙水壓進(jìn)行分析,從圖6中可以看出,空區(qū)形成后,在采空區(qū)內(nèi)部整體上孔隙水壓力呈垂直梯度分布,最大值達(dá)到8567Pa,良好的水力通道在采空區(qū)上方形成,孔隙水壓在采空區(qū)上方附近呈直線下降趨勢.采空區(qū)下方孔隙水壓相對較大,為3191Pa,呈垂直梯度變化,水平層一定范圍內(nèi)孔隙水壓基本一致.

圖6 空區(qū)孔隙水壓力X-550剖面分析圖

通過對水庫下伏采空區(qū)滲流場模擬分析,可得出以下結(jié)論:(1)水庫下伏采空區(qū)區(qū)域,上層覆巖除了會產(chǎn)生由于開采引起覆巖移動(dòng)外,由于孔隙水壓力的減小,其產(chǎn)生的附加有效應(yīng)力同樣會導(dǎo)致覆巖產(chǎn)生一定量的位移,應(yīng)引起重視.(2)采空區(qū)上部破壞后的覆巖透水性較好,如果不進(jìn)行及時(shí)治理,水庫甚至地表水會成為滲流補(bǔ)給的來源,進(jìn)一步加大滲流的水力梯度和流速,水在流動(dòng)的過程中會帶走裂縫中的充填體,久而久之,會形成良好的水力通道,在采空區(qū)與水庫之間形成大量水力通道,一旦通道連通水庫底部,必然會產(chǎn)生嚴(yán)重后果.(3)建議采取措施:①水庫防滲加固處理:對水庫可能存在滲流裂隙帶進(jìn)行及時(shí)加固處理,并對水庫設(shè)立水位監(jiān)測點(diǎn)進(jìn)行觀測,及時(shí)發(fā)現(xiàn)水位陡減情況,排除險(xiǎn)情.②采空區(qū)填充處理:經(jīng)空區(qū)滲水量的監(jiān)測,雖然滲水量呈穩(wěn)定趨勢,但不排除在水流長期滲透沖刷下滲漏有進(jìn)一步惡化的趨勢,建議對該空區(qū)采用高水材料進(jìn)行填充以消除隱患.

6 結(jié) 論

擾動(dòng)后的巖體滲流特性研究主要表現(xiàn)在巖體應(yīng)力場與滲流場的相互作用研究,本文采用大型有限元分析軟件ANSYS,運(yùn)用物理環(huán)境耦合方法,在考慮擾動(dòng)應(yīng)力環(huán)境對巖體滲流參數(shù)的影響下,提出了滲透體積力的加載方程,即將滲透水壓力和體積力轉(zhuǎn)化為單元的結(jié)點(diǎn)外荷載而加以考慮,分析了滲流場影響下的巖體應(yīng)力場,并采用迭代有限逼近的方法,對擾動(dòng)巖體下伏采空區(qū)應(yīng)力場與滲流場進(jìn)行了耦合分析,得出了庫區(qū)下伏采空區(qū)滲流模擬剖面水頭壓力和孔隙水壓力,對實(shí)際工程處理有一定參考作用.

文中擾動(dòng)巖體應(yīng)力環(huán)境下的滲流特性的影響分析只考慮了巖體滲流場滲透系數(shù)k的變化,而影響巖體滲流特性的因素有很多,綜合考慮其他因素情況下,擾動(dòng)巖體滲流場與應(yīng)力場的耦合機(jī)理還有待進(jìn)一步的研究.

[1]任長吉,黃 濤.裂隙巖體滲流場與應(yīng)力場耦合數(shù)學(xué)模型的研究[J].武漢大學(xué)學(xué)報(bào):工學(xué)版,2004,37(2):8-12.

[2]M asanobu Oda.A n Equivalent Continuum M odel for Coupled Stress and Fluid Flow Analysis in Jointed Rock Masses[J].Water Resources Research,1986(13):98-103.

[3]Yang L Z,Huang T.Laboratory Study of Low Permeability Water Bearing Medium in Environmental Pollution Control[A].San Diego,Geotechnical Laboratory,Lehigh University,1996:263-271.

[4]毛昶熙.滲流計(jì)算分析與控制[M].北京:水利電力出版社,1988.

[5]許光祥,哈秋舲.裂隙巖體滲流與卸荷力學(xué)耦合作用及裂隙排水研究[M].重慶:重慶大學(xué)出版社,2003.

[6]楊立中,周 訓(xùn).深層地下水滲流的研究[M].成都:成都科技大學(xué)出版社,1995.

[7]周維桓.高等巖石力學(xué)[M].北京:水力水電出版社,1990.

[8]張金才,張玉卓.應(yīng)力對裂隙巖體滲流影響的研究[J].巖土工程學(xué)報(bào),1998,20(2):19-22.

[9]陳 平,張有天.裂隙巖體滲流與應(yīng)力耦合分析[J].巖石力學(xué)與工程學(xué)報(bào),1994,13(4):299-309.

[10]曹阿靜,許夢國.裂隙巖體滲流與應(yīng)力耦合有限元分析[J].礦業(yè)研究與開發(fā),2003,23(3):5-8.

[11]張社榮,楊璐玲,鐘登華.裂隙巖體滲流場分析及其三維有限元程序設(shè)計(jì)[J].巖土力學(xué),2005,26(1):46-50.

[12]蔣中明,Dashnor Hoxha.核廢料地質(zhì)貯存介質(zhì)黏土巖的三維各向異性熱-水-力耦合數(shù)值模擬[J].巖石力學(xué)與工程學(xué)報(bào),2007,26(3):493-501.

[13]Ouyang Zhihua,Elsworth D.Evaluation of Groundwater Flow into Mined Panels[J].Int J.Rock Mech.M in Sci,1993,30(2):71-79.