拉西瓦拱壩壩基滲流應(yīng)力耦合分析

水利水電工程專業(yè)畢業(yè)論文 摘 要 巖體的應(yīng)力場導(dǎo)致裂隙的幾何特性發(fā)生變化,從而影響巖體的滲透性,引起滲 流場的改變。 反過來， 滲流場的存在將改變滲透體積力的分布,必將對應(yīng)力場產(chǎn)生 影響。因此，對于水工結(jié)構(gòu)，考慮巖體滲流場與應(yīng)力場的耦合作用是十分必要的。 本文在總結(jié)國內(nèi)外學(xué)者研究成果的基礎(chǔ)上，結(jié)合拉西瓦工程，分析了壩基滲流應(yīng) 力耦合特性，主要內(nèi)容和研究成果歸納如下: （1） 從試驗(yàn)研究、 理論模型研究和耦合模型計算方法研究等方面較為深入地 總結(jié)了巖體滲流應(yīng)力耦合問題的研究進(jìn)展。 （2）介紹了應(yīng)力場與滲流場耦合分析的基本原理和方程。 （3）編制了 Matlab 與 COMSOL Multiphysics 的鏈接程序，通過導(dǎo)入控制坐 標(biāo)，建立了拉西瓦拱壩和壩基的三維有限元數(shù)值模型。 （4） 計算分析了不考慮耦合情況下拉西瓦工程在水庫正常運(yùn)行期的滲流場與 應(yīng)力場；確定滲流應(yīng)力耦合模型的相關(guān)參數(shù)和邊界條件，計算分析了考慮耦合情 況下壩體和壩基在水庫正常運(yùn)行期的耦合滲流場和應(yīng)力場，總結(jié)了其變化規(guī)律。 結(jié)果表明上游面滲透水由河道往兩岸的逐漸減小，下游面滲透水頭由河道往兩岸 的逐漸增大， 總體上由上游到下游呈減小趨勢； 壩體順河向位移量值隨高程增加， 同一高程從拱冠至左、右拱端逐漸減小，壩體左、右半拱順河向位移具有良好的 對稱性。 （5） 分別對比分析了考慮耦合和不考慮耦合情況下滲流場和應(yīng)力場的計算成 果，結(jié)果表明在考慮耦合情況下，滲流場等勢線偏向下游，下游的滲透坡降明顯 增大；各位移分量都有不同程度的增加；壩體上游巖體的拉應(yīng)力增加，壩肩和壩 基下游巖體的壓應(yīng)力增加，對壩體的穩(wěn)定不利。 關(guān)鍵詞：拉西瓦拱壩；滲流；應(yīng)力；耦合；COMSOL Multiphysics 水利水電工程專業(yè)畢業(yè)論文 Abstract The change of the rock stress field will lead to the crack geometry characteristic changing, which affects the permeability and causes the seepage field change. In turn, will change the distribution of permeability volume force and influence the stress field. Therefore, It is imperative for Hydraulic Structure to be designed based on the coupling of seepage field and stress field. On the basis of research achievements of domestic and foreign scholars, this paper mainly combines with the Lashiwa engineering to analyze the coupling characteristics further. The main content and research is summarized as follows: (1) Research progresses of the coupling of seepage field and stress field in rocks are comprehensively summarized from several aspects, such as the experimental research , the theoretical model research and the calculation method research of coupling model. (2) The basic principles and equations to analysis coupled seepage and stress fields are obtained. (3) Work out a program which connects the Matlab software with the COMSOL Multiphysics software, through importing coordinate control point, to establish the 3-D finite element numerical model of the Laxiwa arch dam and its foundation. (4) Considering uncoupled cases, calculate and analyze the seepage field and stress field of the Laxiwa engineering in the reservoir normal operation periods to determine the related parameters and boundary conditions of the coupling model of seepage and stress. Then, analyze the seepage field and stress field in the reservoir normal operation periods under the coupled situations, with which the distributive regularity of seepage fields and stress fields are summarized. The results show that in upstream faces hydraulic pressure head reduces from the river to both sides gradually, but in downstream faces that the distribution is the contrary. In general, the hydraulic pressure head reduces from upstream to downstream. The displacement of dam along river direction increases with elevation and decreases gradually from the crown to arch abutment at the same elevation, and it has a good symmetry perpendicular to river direction. (5) The calculation results of seepage field and stress field between uncoupled and coupled situations are compared and analyzed respectively. It shows that the isopotential lines of seepage field are to downstream, the seepage slope of downstream increases obviously, each displacement component increases at different rates, the tensile stress of the upstream rock increases and the compressive stress of the downstream rock also increases in the coupling cases, which is bad for the stability of the dam. Keywords: Laxiwa arch dam; seepage; stress; coupling; COMSOL Multiphysics 水利水電工程專業(yè)畢業(yè)論文 目 錄 第一章 緒論1 1.1 問題的提出 1 1.1.1 工程破壞實(shí)例 . 1 1.1.2 討論. 2 1.2 國內(nèi)外研究進(jìn)展 4 1.2.1 試驗(yàn)研究成果 . 4 1.2.2 理論模型研究成果 8 1.2.3 耦合模型計算方法研究 9 1.2.4 工程應(yīng)用前景展望 10 1.3 研究內(nèi)容和方法 .11 1.4 技術(shù)路線 12 第二章 滲流應(yīng)力耦合分析數(shù)學(xué)模型13 2.1 基本假設(shè) 13 2.2 滲流場與應(yīng)力場耦合分析數(shù)學(xué)模型 13 2.2.1 應(yīng)力場的基本方程 13 2.2.2 滲流場的基本方程 16 2.2.3 耦合場的基本方程 17 2.3 COMSOL Multiphysics 軟件 18 2.3.1 COMSOL Multiphysics 軟件簡介 . 18 2.3.2 算例驗(yàn)證. 21 2.4 本章小結(jié) 23 第三章 拉西瓦拱壩壩基滲流應(yīng)力耦合分析24 3.1 拉西瓦基本資料 24 3.1.1 工程概況. 24 3.1.2 地形地貌. 25 3.1.3 地層巖性. 25 3.1.4 壩址區(qū)水文地質(zhì)條件 26 3.1.5 巖體滲透性. 27 3.2 計算模型及計算參數(shù). 27 3.2.1 計算模型. 27 拉西瓦拱壩壩基滲流應(yīng)力耦合分析 3.2.2 計算參數(shù). 29 3.3 計算模型有限元網(wǎng)格. 30 3.4 計算工況 31 3.5 計算成果及成果分析. 31 3.5.1 滲流場計算成果及成果分析 31 3.5.2 應(yīng)力場計算成果及成果分析 38 3.6 本章小結(jié) 45 第四章 總結(jié)與展望.46 4.1 總結(jié) 46 4.2 研究展望 46 參考文獻(xiàn).48 致 謝51 水利水電工程專業(yè)畢業(yè)論文 1 第一章 緒論 1.1 問題的提出 隨著國民經(jīng)濟(jì)建設(shè)事業(yè)的快速發(fā)展，越來越多的水利水電工程、交通、能源、 工民建和防護(hù)工程等被建筑在巖石地區(qū)。天然巖體中存在大量的不連續(xù)的地質(zhì)結(jié) 構(gòu)面，如節(jié)理、斷層、裂隙、巖層層面、片理面等，這些結(jié)構(gòu)面不但大大改變了 巖體的力學(xué)性質(zhì),也嚴(yán)重影響了巖體的滲透特性1。 國內(nèi)外發(fā)生的巖土工程事故中，90%以上的巖體邊坡破壞和地下水滲透力有 關(guān)，60%的礦井事故與地下水作用有關(guān)，30%40%的水電工程大壩失事是由滲透 作用引起的2。此外，隧道開挖、地下核廢料存儲、油氣開發(fā)等都涉及到巖體應(yīng) 力、工程干擾力和地下水滲透力的相互作用及其耦合問題。因此，對裂隙巖體的 滲流應(yīng)力耦合分析事關(guān)工程建設(shè)的成敗與安全，對整個工程的可行性、安全性和 經(jīng)濟(jì)性起到重要的制約作用，并在很大程度上影響著工程的投資及使用效益。目 前，我國壩工建設(shè)正在向高壩發(fā)展，對裂隙巖體的滲流應(yīng)力耦合特性的研究極其 重要。 1.1.1 工程破壞實(shí)例 （1）法國 Malpasset 拱壩潰壩 Malpasset 雙曲薄拱壩位于法國南部 Rayran 河上，壩高 66m，于 1954 年末建 成并蓄水。12 月 2 日，一場大雨加快了水位的上升速度，21 點(diǎn) 20 分，大壩突然 潰決，當(dāng)時庫水位為 100.12m。根據(jù)壩下游 1.5km 對這一災(zāi)難少數(shù)幸存者描述， 他們首先感到一陣劇烈顫動， 隨之聽到 類似動物吼叫的突發(fā)巨響， 然后感到強(qiáng) 烈的空氣波。 最終他們看到巨大的水墻 順河谷奔騰， 洪水出峽谷后流速仍高達(dá) 20km/h 左右，同一時間電力供應(yīng)中斷。 次日清晨發(fā)現(xiàn)大壩已被沖走， 僅右岸靠 基礎(chǔ)部分有殘留拱壩， 一些壩塊被沖到 下游 1.5km 處， 左岸壩基巖體被沖出深 槽 （如圖 1.1 所示） 。 下游 12km處 Frejus 城鎮(zhèn)部分被毀，死亡 421 人，財產(chǎn)損失 圖 1.1 左岸壩肩， 壩基巖體被沖出深槽 拉西瓦拱壩壩基滲流應(yīng)力耦合分析 2 達(dá) 300 億法郎3。 Malpasset 拱壩潰壩引起了工程界的廣泛重視，因?yàn)槭鹿拾l(fā)生在壩工建設(shè)方 面，尤其是在拱壩建設(shè)方面為世界最先進(jìn)的國家，并且是由最負(fù)盛名的設(shè)計大師 Andce Coyne 設(shè)計的；它是當(dāng)時潰壩記錄中最高的壩，而在此之前尚沒有拱壩潰 壩的先例；潰壩幾乎毀滅了 Frejus 市，使最富的地中海地區(qū)遭受了重大損失。 （2）意大利 Vajont 拱壩近壩庫岸巖體大滑坡 Vajont 雙曲拱壩位于意大利 Piave 河支流 Vajont 河上，壩高 265m，水庫設(shè)計 正常高水位為 722.5m，于 1959 年建成，是當(dāng)時世界上最高的拱壩。1963 年 10 月 9 日夜，Vajont 水庫水位達(dá) 700m 高程，大壩上游近壩庫左岸約 2.5 億 m3巨大 巖體突然發(fā)生高速滑坡，以 25m/s 的速度沖入水庫，使 5500 萬 m3的庫水產(chǎn)生巨 大涌浪，約有 3000 萬 m3的水翻越 壩頂泄入底寬僅 20m 的狹窄河谷。 翻壩的水流在右岸超出壩頂高程 達(dá) 250m，左岸達(dá) 150m。水流以巨 大的流速滾向下游， 經(jīng) Vajont 河沖 入 Piave 河時，經(jīng)接近直角的彎轉(zhuǎn) 后席卷了 Longarone 小鎮(zhèn)及幾個附 近村莊，造成生命財產(chǎn)的巨大損 失，約有 2500 人死亡，大滑坡的 石渣掩埋了水庫，使大壩和水庫完 全報廢，如圖 1.2 所示。 1.1.2 討論 Malpasset 拱壩失事至今已 50 多年，對其失事的原因至今尚未取得完全一致 的認(rèn)識，比較權(quán)威的有以下兩種認(rèn)識：Bellier(1976)認(rèn)為，Malpasset 大壩壩肩巖 體的結(jié)構(gòu)面走向與拱壩推力方向平行，在壩肩巖體中形成高的壓應(yīng)力區(qū)，引起壩 肩巖體結(jié)構(gòu)面閉合，巖體的滲透系數(shù)減小到約為初始值的 1/100 甚至更小，導(dǎo)致 巖體滲流受阻，而產(chǎn)生了等于水庫全水頭的壓力作用于該巖體之上，使壩肩巖塊 沿下游斷層滑移而失穩(wěn)；Wikkte 和 Leonards(1985)認(rèn)為，拱壩弧部受拉力，使大 壩壩踵附近巖體受拉，傾向下游的巖體結(jié)構(gòu)面張開，裂縫使帷幕短路。庫區(qū)蓄水 后，庫水沿張開裂縫滲透，由于下游斷層封堵了滲透通道，致使張開裂縫中產(chǎn)生 等于水庫全水頭的靜水壓力，使壩肩巖塊失穩(wěn)4。盡管這兩種分析出發(fā)點(diǎn)不同， 但卻得出同一結(jié)論，即巖體的滲透壓力作用引起的巖體變形導(dǎo)致了大壩的破壞。 圖 1.2 Vajont 拱壩右岸滑坡 水利水電工程專業(yè)畢業(yè)論文 3 意大利瓦依昂大壩上游庫區(qū)大滑坡，主要是由于庫水位上升，引起巖體中地 下水位抬升，增加了巖體的浮托力而引起的一種滲透變形和破壞。由于該水庫兩 岸巖體滲透性良好，庫水位上升迅速引起地下水位上升，滯后時間短。因而，當(dāng) 庫水位上升到 700m 時，在庫岸邊坡坡腳處產(chǎn)生極大的揚(yáng)壓力，正應(yīng)力減小，抵 抗滑坡的摩擦阻力減小而引起滑坡。 這兩次重大工程事故均發(fā)生在 20 世紀(jì)六七十年代， 當(dāng)時是水利水電工程建設(shè) 的高潮時期，但是工程經(jīng)驗(yàn)相對不足，特別是對巖石的水力特性了解很少，對巖 體中的滲流問題還是完全照搬土體滲流學(xué)，即孔隙介質(zhì)滲流學(xué)的方法及經(jīng)驗(yàn)來解 決3。然而巖體與作為孔隙介質(zhì)的土體有很大差別：土體是相對較均一的多孔介 質(zhì)，土顆粒及孔隙尺寸都很小，其樣本體積與所研究滲流域相比甚小，因此可以 把土體抽象為連續(xù)滲流介質(zhì)；而巖體則是裂隙一塊體離散體介質(zhì)，其中存在著各 種各樣的裂隙結(jié)構(gòu)面及斷層，這些結(jié)構(gòu)面不僅把完整的巖石分割開，而且它們既 是巖體中的軟弱面，又是巖體中的主要透水通道，而巖塊本身的透水性在工程上 往往可以忽略不計。 根據(jù)國內(nèi)外大壩失事實(shí)例的報道，在比較多的情況下大壩破壞是由于巖土介 質(zhì)中的滲流場和應(yīng)力場的相互作用引起的；而通常在壩設(shè)計時沒有考慮滲流與應(yīng) 力的相互作用。由此可見，進(jìn)行大壩工程設(shè)計時，必須動態(tài)地考慮大壩及其周圍 巖土中滲流與應(yīng)力的相互作用，而對滲流場與應(yīng)力場進(jìn)行耦合分析則是對其相互 作用進(jìn)行定量研究的必由之路。 目前，對大壩進(jìn)行滲流場的分析與應(yīng)力場的分析大多是分開各自進(jìn)行的，即 在進(jìn)行滲流場的分析時不考慮壩體以及壩基巖土體的應(yīng)力-應(yīng)變關(guān)系和應(yīng)力狀態(tài) 圖 1.3 巖體滲流應(yīng)力耦合作用示意圖 滲流場 流體 水頭 H、水量 Q 應(yīng)力場 巖體 應(yīng)力、應(yīng)變 力學(xué)參數(shù) 水力參數(shù) 巖體中孔隙體積 或裂隙開度變化 拉西瓦拱壩壩基滲流應(yīng)力耦合分析 4 及其改變對滲流場的影響；在進(jìn)行滲流場的計算時，根據(jù)滲流計算結(jié)果，賦予不 同區(qū)域巖土體不同容重來進(jìn)行穩(wěn)定及應(yīng)力分析。盡管這種方法簡單易行且已積累 了一定的工程實(shí)踐經(jīng)驗(yàn)，但是沒有真實(shí)客觀地反映滲流場與應(yīng)力場之間的相互作 用、相互影響。事實(shí)上，滲流場與應(yīng)力場是相互影響、不可分離的，如圖1.3所示: 一方面，滲流場的改變會導(dǎo)致滲流體積力和滲透壓力的改變，使作用在壩體上的 外荷載發(fā)生變化，從而引起壩體應(yīng)力場的重新分布；另一方面，應(yīng)力場的改變， 會引起壩體體積應(yīng)變的改變，使壩體內(nèi)各部分的孔隙率發(fā)生改變，滲透系數(shù)也隨 之變化，從而改變壩體滲流場的分布；兩者的這種相互作用結(jié)果，會使壩體達(dá)到 一種雙場的耦合平衡狀態(tài)，形成滲流場影響下的穩(wěn)定應(yīng)力場和應(yīng)力場影響下的穩(wěn) 定滲流場。研究巖體滲流應(yīng)力耦合問題不僅有重要的理論意義，而且在工程中有 廣泛的實(shí)際應(yīng)用價值。 1.2 國內(nèi)外研究進(jìn)展 隨著近年來技術(shù)水平和標(biāo)準(zhǔn)要求的不斷提高該，國內(nèi)外的學(xué)者在滲流、應(yīng)力 及兩者之間耦合關(guān)系做了大量的研究工作，主要包括試驗(yàn)研究、理論模型研究和 耦合模型計算方法研究三大類。 1.2.1 試驗(yàn)研究成果 1.2.1.1 單裂隙巖體試驗(yàn)研究 立方定律是描述平行板裂隙面水流運(yùn)動規(guī)律的一個著名定律，它將裂隙抽象 為兩片光滑、平直、無限長的平行板，粘性不可壓縮流體在平板內(nèi)做定常運(yùn)動， 根據(jù)流體力學(xué)原理推導(dǎo)出平行板裂隙面水流公式：J ge q 12 3 （其中，q為單寬流 量；e為裂隙寬度；g為重力加速度； 為水流運(yùn)動粘滯系數(shù)；J 為水力梯度） ， 此式表明裂隙面的單寬流量與隙寬的 3 次方成正比。 Lomize、Louis 等進(jìn)行了單裂隙的水流試驗(yàn)研究，證明了在層流時立方定理 的有效性。Romm 通過對微裂隙和極微裂隙的研究，提出只要隙寬大于 0.2m， 立方定理總是成立的5。 但在天然情況下的裂隙面大多是粗糙不平的，很難滿足平行板裂隙的假定， 許多學(xué)者進(jìn)行了仿天然裂隙的試驗(yàn)研究，對立方定理提出了各種各樣的修正，具 體如表 1.1 所示： 水利水電工程專業(yè)畢業(yè)論文 5 表 1.1 各學(xué)者對立方定理提出的修正611 學(xué)者 公式 符號含義 說明 Lomize 5 . 1 )(61 e C Louis 5 . 1 ) 2 (8 . 81 e C 速寶玉 75. 0 )(2 . 11 e C Amadei C J eg q 1 12 3 2 . 1 )(6 . 01 e C e 仿天然裂隙的試驗(yàn)得出， 其中：e為平均隙寬，C 1 為立方定理的修正系數(shù)， 與裂隙面的粗糙度及隙 寬情況有關(guān)，為裂隙粗 糙度； e 為隙寬均方差 Barton J ge JRC q m 12 1 6 5 . 7 h e為等效水力隙寬； m e 為力學(xué)隙寬 通 過 大 量 試 驗(yàn) ， 提 出 JRC(節(jié)理粗糙度系數(shù))修 正法 耿克勤 n m AeQ 對于小開度裂隙層流， 1.7 n3.0；對于中開度， 0.8n1.4；對于大開 度，0.3n0.48 人工、 天然光滑和粗糙裂 隙的試驗(yàn) Nolte n m CeQQ 0 n隨隙寬的變化范圍為 7.69.8 經(jīng)驗(yàn)公式 Iwai 1 1 0 Q Q Q為相應(yīng)于接觸面積率 的流量， 0 Q為0時 的流量，為經(jīng)驗(yàn)系數(shù) Walsh 和周創(chuàng)兵分別通過 模擬熱傳導(dǎo)理論和數(shù)學(xué) 推導(dǎo)得出相似的結(jié)論。 1.2.1.2 應(yīng)力對裂隙滲流影響的試驗(yàn)研究 立方定律和修正的立方定律都沒有考慮應(yīng)力對裂隙滲流的影響，而建立滲流 與應(yīng)力的關(guān)系是分析裂隙巖體滲流場與應(yīng)力場耦合的基礎(chǔ)和關(guān)鍵環(huán)節(jié)。因此，許 多學(xué)者在此方面做了很大的努力，具體見表 1.2 所示： 拉西瓦拱壩壩基滲流應(yīng)力耦合分析 6 表 1.2 應(yīng)力對裂隙滲流影響的研究成果1217 學(xué)者 公式 符號意義 說明 Snow (1968) )( )2( 0 2 N f0f S bK KK f K為水平裂隙的滲透系數(shù)； f0 K 為初始滲透系數(shù)； N K為裂隙的 法向剛度；S為裂隙的隙間距； 0 為初始應(yīng)力 多組平行裂隙 滲透系數(shù)與法 向應(yīng)力關(guān)系 Louis (1974) ekkf 0 0 k為初始滲透系數(shù)；為有效正 應(yīng)力PH ,H為研究點(diǎn) 上的上覆巖層的重量，P為孔隙 水壓力；為經(jīng)驗(yàn)系數(shù)，取決于 巖石中的裂隙狀態(tài) 根據(jù)均質(zhì)裂隙 巖體不同深度 的鉆孔壓水試 驗(yàn)成果建立的 巖體滲透系數(shù) f k與法向應(yīng)力 的經(jīng)驗(yàn)關(guān)系 式 Jones (1975) 3 0 ) p lg( k p KK p為有效應(yīng)力； k p為 f K時愈合 有效應(yīng)力 碳酸鈣巖石裂 隙滲透系數(shù)的 經(jīng)驗(yàn)方程 Kranz 等 n PPQAK )-( c0 A為過水面積； c P為總壓力；P 為內(nèi)部孔隙水壓力；n為常數(shù) Barre 花崗巖的 裂隙滲透系數(shù) 的經(jīng)驗(yàn)公式 Nelson -n f pkBA P 為有效壓力；A，B，n 均為常 數(shù) 砂巖裂隙滲透 系數(shù)經(jīng)驗(yàn)公式 Galet f T 12/gb3 h f T；、為參數(shù)； 為法向應(yīng)力 由花崗巖、大理 巖、玄武巖三種 巖體裂隙的室 內(nèi)試驗(yàn)得出 張 玉 卓 等 4 1 yx BAQ x 、 y 為雙向主應(yīng)力；A，B 均為常數(shù) 裂隙滲流與應(yīng) 力耦合雙向等 壓和雙向不等 壓試驗(yàn) 水利水電工程專業(yè)畢業(yè)論文 7 鄭 少 河 等 a f pkk )( 3120 f k為天然單裂隙滲透系數(shù)； 1 、 2 、 3 為三向主應(yīng)力，其中 2 垂直裂隙面， 1 、 o z y x 平行于f裂 隙面；為泊松比；p為裂隙水 壓力；a為系數(shù)，取決于裂隙面 的粗糙度 通過大量天然 裂隙滲流實(shí)驗(yàn)， 揭示了三維應(yīng) 力裂隙滲流規(guī) 律，提出了裂隙 閉合量與巖體 三維應(yīng)力，巖體 滲透系數(shù)與裂 隙閉合量的關(guān) 系式 劉繼山 3 0 0 1 1 i f m i R kk （zyxi,） m R為巖體分類指數(shù)； 0 f 為初始 孔隙度； ii k,分別為對應(yīng)方向 的滲透系數(shù)及應(yīng)變； 0 k為常數(shù) 單裂隙和兩正 交裂隙受正應(yīng) 力時的滲流公 式 忤 彥 卿 (1995) f 0 D KK K 為應(yīng)力等于 s 時的滲透系數(shù)； 0 K為應(yīng)力等于 0 時的滲透系 數(shù) ； f D為 裂 隙 的 分 維 數(shù) （20 f D） 當(dāng)0 f D時， 用 Louise 負(fù)指 數(shù)關(guān)系式；當(dāng) 21 f D時， 說明巖體中裂 隙比較發(fā)育 劉 亞 晨 等 )exp( 3 TAkf 3 為法向壓應(yīng)力； T 為溫度； A， ，均為系數(shù)，可通過試驗(yàn)獲 得 通過單、正交裂 隙花崗巖的試 驗(yàn)，研究了高 溫、高壓下的裂 隙巖體滲透特 性 以上試驗(yàn)均沒有考慮剪應(yīng)力的影響，裂隙巖體滲流與剪應(yīng)力耦合作用所做的 研究較少。剪應(yīng)力作用下裂隙面幾何特征的變化比較復(fù)雜，而且試驗(yàn)難度也大。 Esaki T. 等采用一種新的實(shí)驗(yàn)技術(shù)對花崗巖人工的裂隙試樣進(jìn)行了滲流與剪 應(yīng)力耦合研究。Chen Z 也同樣得出滲透系數(shù)隨剪切位移增加的曲線。劉才華等18 通過對規(guī)則、均勻、粗糙裂隙的滲流剪切實(shí)驗(yàn)，研究了巖體裂隙在剪切荷載作用 下的滲流特性。劉才華等19通過室內(nèi)模擬實(shí)驗(yàn)，探討了低應(yīng)力、低水頭條件下裂 拉西瓦拱壩壩基滲流應(yīng)力耦合分析 8 隙巖體受剪應(yīng)力作用的滲流特性，建立了剪應(yīng)力作用下巖石裂隙的滲流模型。 Takahashi 在砂巖和花崗巖的剪切和滲流耦合分析方面進(jìn)行了較深入的研究。 耿克 勤20對剪切變形與滲流耦合進(jìn)行試驗(yàn)，結(jié)果顯示在不同壓應(yīng)力作用下裂隙面由于 剪切變形會使?jié)B流量發(fā)生先減少（相應(yīng)于剪縮） 、后來又增大（相應(yīng)于剪脹）的變 化。裂隙面剪切變形對滲流系數(shù)的影響取決于裂隙面的幾何形態(tài)、剪切位移和法 向應(yīng)力的大小。 1.2.2 理論模型研究成果 從 20 世紀(jì)六七十年代開始,國內(nèi)外很多學(xué)者對裂隙巖體的滲流模型作了許多 研究,建立了等效連續(xù)介質(zhì)、裂隙網(wǎng)絡(luò)介質(zhì)、雙重介質(zhì)等模型。 1.2.2.1 等效連續(xù)介質(zhì)模型 等效連續(xù)介質(zhì)模型把孔隙裂隙系統(tǒng)等效成連續(xù)介質(zhì)，用連續(xù)介質(zhì)理論描述 滲流方程。模型將裂隙水流等效平均到整個巖體，把裂隙巖體模擬為具有對稱滲 透張量的各向異性連續(xù)體，利用經(jīng)典的連續(xù)介質(zhì)理論進(jìn)行分析21。其優(yōu)點(diǎn)是可以 直接應(yīng)用較成熟的孔隙介質(zhì)飽和 滲流分析方法，可以沿用各向異性 連續(xù)介質(zhì)理論來描述，可操作性 好、理論基礎(chǔ)雄厚。不足之處主要 在于把裂隙網(wǎng)絡(luò)等效為連續(xù)介質(zhì)， 不能很好地刻畫裂隙的特殊導(dǎo)水 作用，由于不是所有巖體均可等效 為連續(xù)介質(zhì)； 另外， 典型單元體(定 義如圖1.4所示)的大小和等效水力 參數(shù)也較難確定，模型的適用范圍 經(jīng)常受限制。 肖裕行22等從力的運(yùn)動效果和對應(yīng)應(yīng)力狀態(tài)兩方面探討了滲透作用力，證明 了水力等效連續(xù)介質(zhì)中滲透作用力運(yùn)動效果等效的必然性和對應(yīng)應(yīng)力狀態(tài)性質(zhì)的 一致性，繼而建立裂隙巖體水力等效連續(xù)介質(zhì)是否存在的二維判別方法。Long， 仵彥卿(1996)23，王媛等2425(1998,2000)均在這方面做了大量的研究。 1.2.2.2 裂隙網(wǎng)絡(luò)介質(zhì)模型 裂隙網(wǎng)絡(luò)介質(zhì)模型把裂隙介質(zhì)視為由不同大小和方向的裂隙在空間交叉構(gòu)成 的網(wǎng)絡(luò)狀系統(tǒng)，水只能在其中運(yùn)動。該模型是在查清每條裂隙的空間方位、隙寬 圖 1.4 裂隙巖石典型單元體的定義 水利水電工程專業(yè)畢業(yè)論文 9 等幾何參數(shù)的前提下，以單個裂隙水流基本公式為前提，根據(jù)流入和流出各裂隙 交叉點(diǎn)的流量相等來求其水頭值， 進(jìn)而將水壓力作用在每條裂隙上進(jìn)行受力分析。 它較好地描述了裂隙巖體的非均勻各向異性，故當(dāng)巖體中巖塊的滲透性和裂隙的 滲透性相比非常小時， 可以認(rèn)為水流只在裂隙網(wǎng)絡(luò)中定向流動。 該模型擬真性好、 精度高，但當(dāng)裂隙較多時，其工作量相當(dāng)大，特別是三維問題，甚至是不可能實(shí) 現(xiàn)。另外，由于裂隙分布具有隨機(jī)性，要建立離散的真實(shí)裂隙系統(tǒng)也是十分困難。 Wilson 和 Witherspoon 分別以三角形單元或線單元模擬巖體中的裂隙，提出 模擬二維裂隙網(wǎng)絡(luò)水流的兩種有限元技術(shù)，并以算例表明裂隙交叉點(diǎn)的水流干擾 是可以忽略的， 從而闡明了采用線單元的優(yōu)越性和可行性26。 Long 在假定巖體中 裂隙發(fā)育呈圓盤狀的基礎(chǔ)上提出三維圓盤裂隙網(wǎng)絡(luò)模型，并采用混合解析一數(shù)值 方法對此進(jìn)行求解。Dershowitz 假定裂隙在巖體中互相切割成多邊形，建立多邊 形裂隙網(wǎng)絡(luò)模型。國內(nèi)毛昶熙27，王恩志28，仵彥卿29，高海鷹30等對裂隙網(wǎng)絡(luò) 介質(zhì)模型的研究做了大量的貢獻(xiàn)。 1.2.2.3 雙重介質(zhì)模型 雙重介質(zhì)模型把巖體看作由孔隙和裂隙組成的雙重介質(zhì)空隙結(jié)構(gòu)，孔隙介質(zhì) 和裂隙介質(zhì)均布在滲流區(qū)域內(nèi)，形成連續(xù)介質(zhì)系統(tǒng)。在該系統(tǒng)內(nèi)，孔隙介質(zhì)儲水， 裂隙介質(zhì)導(dǎo)水，每一點(diǎn)處存在兩種液體壓力，研究兩者之間的水交換過程及由此 帶來的兩種介質(zhì)的變形，可細(xì)分為雙孔隙度模型、雙滲透性模型以及雙孔隙度 雙滲透性模型。 雙重介質(zhì)模型在一定程度上刻畫了優(yōu)先流的現(xiàn)象，而且考慮了巖塊與裂隙間 客觀存在的能量交換，具有較好的擬真性；但數(shù)值分析工作量更大。王恩志31等 把控制滲流總體分布且起主導(dǎo)滲透作用的大裂隙定義為裂隙巖體的主干裂隙網(wǎng) 絡(luò)，將主干裂隙網(wǎng)絡(luò)間的巖塊定義為裂隙巖塊，再建立各自的滲流方程，并通過 兩類系統(tǒng)間的水量交換方程構(gòu)建雙重裂隙系統(tǒng)滲流模型。 1.2.3 耦合模型計算方法研究 國內(nèi)外許多學(xué)者通過對滲流應(yīng)力耦合分析的廣泛研究，提出了不同的耦合計 算方法，主要包括迭代耦合法、直接耦合法、等效連續(xù)法和離散法。 1.2.3.1 迭代耦合法 迭代法通常包含兩個或多個分析模塊,不同的模塊模擬不同物理環(huán)境中的變 量演化過程,各個模塊的計算結(jié)果互為外加荷載或邊界條件32。這種方法通過兩 場分開計算,然后進(jìn)行兩場間交叉迭代計算,最后達(dá)到求解的目的。陳平、張有天 拉西瓦拱壩壩基滲流應(yīng)力耦合分析 10 33以裂隙滲流理論和變形本構(gòu)關(guān)系為基礎(chǔ)進(jìn)行了滲流耦合分析，裂隙的滲透系數(shù) 采用經(jīng)驗(yàn)公式。 此法在計算方法上容易實(shí)現(xiàn),可以利用已成熟的固體變形場和流體 滲流場的求解程序,但未能將滲流場和應(yīng)力場作為同一場考慮,從方法上說不能稱 為嚴(yán)格的耦合算法。 1.2.3.2 直接耦合法 直接耦合法，也稱為多自由度法，是根據(jù)耦合模型涉及的所有物理變量，建 立與各變量相應(yīng)的剛度矩陣、影響矩陣并形成總耦合矩陣，以及相應(yīng)的耦合荷載 向量34。直接法的物理意義更為清晰直觀,計算過程中的參數(shù)易于控制和調(diào)整。 耦合模型的物理環(huán)境及相應(yīng)的偏微分控制方程一般不少于兩個，而且對于同一個 單元結(jié)果， 不同的物理變量往往具有不同的屬性和自由度數(shù),由此給相關(guān)的數(shù)值計 算尤其是組剛過程中的編碼帶來一定的難度。 此外在處理多相流和非線性介質(zhì)時, 很難把握方程的收斂性。對于這種方法的研究和運(yùn)用，王媛35提出一種考慮巖體 滲透性隨應(yīng)力變化的、且適用于裂隙巖體滲流場和應(yīng)力場耦合分析的四自由度全 耦合法，楊志錫36等進(jìn)行了各向異性飽和土體的滲流耦合分析和數(shù)值模擬。 1.2.3.3 等效連續(xù)法 Noorishad 等將土體固結(jié)分析方法用于建立巖體滲流與應(yīng)力耦合關(guān)系， 進(jìn)而用 有限元方法研究了非連續(xù)介質(zhì)中的固液兩相介質(zhì)的耦合問題。這一方法忽略了滲 透性隨應(yīng)力的變化，同時用節(jié)理單元處理如同斷層一樣的幾個主要不連續(xù)面，確 實(shí)是一種有效手段，但是用節(jié)理單元解決所有不連續(xù)面是不切實(shí)際的。小的節(jié)理 不能單獨(dú)在滲流中起重要作用，而它們連通后卻也成為主要滲流路徑，不考慮這 樣一些裂隙將導(dǎo)致滲流分析中的錯誤的結(jié)果，因此目前這一方法在實(shí)際工程中還 難以應(yīng)用。 1.2.3.4 離散法 采用離散單元法時，由熱平衡方程、塊體和不連續(xù)體的力學(xué)動力方程以及不 連續(xù)體的水力方程組成控制方程。離散法在進(jìn)行裂隙巖體滲流應(yīng)力耦合分析時， 大多采用將滲流場和應(yīng)力場分開計算再交叉迭代求解來達(dá)到耦合分析的目的。其 求解較為繁瑣，而且隨著裂隙數(shù)量的增加，該方法甚至不可行。 1.2.4 工程應(yīng)用前景展望 裂隙巖體滲流應(yīng)力耦合研究在處理實(shí)際工程中，尤其是在大壩壩基穩(wěn)定性方 面發(fā)揮了重要作用。它的應(yīng)用前景相當(dāng)廣泛,主要體現(xiàn)在以下幾個方面： 水利水電工程專業(yè)畢業(yè)論文 11 （1）水利水電工程方面：南水北調(diào)、長江三峽等大型水利工程、拉西瓦拱壩 等地質(zhì)條件復(fù)雜地區(qū)都必須考慮裂隙巖體滲流與應(yīng)力的相互影響； （2）石油開采方面：對裂縫性砂巖油藏的開采具有重要意義； （3） 邊坡工程方面： 對于高陡邊坡,為了避免出現(xiàn)滑坡現(xiàn)象,必須考慮地表水 入侵和巖體裂隙滲流作用下的邊坡穩(wěn)定性問題； （4）城市地下水的抽取方面：可以解決抽水量和地面沉降等相關(guān)環(huán)境問題; （5）隧道及地下工程方面：由于基巖中存在裂隙水,需要考慮地下水作用下 工程的穩(wěn)定性問題,例如近年來西部大開發(fā)中一些鐵路、 公路受地下水影響非常明 顯，隧道工程中的滲流應(yīng)力耦合問題也越來越突出； （6）采礦工程方面：工作面和巷道突水問題發(fā)生頻繁且后果嚴(yán)重，隨著采礦 工程向深部開展,需要考慮深部采動裂隙巖體的應(yīng)力場、 溫度場、 滲流場三場耦合 問題。 1.3 研究內(nèi)容和方法 本文結(jié)合拉西瓦工程，在總結(jié)國內(nèi)外學(xué)者研究成果的基礎(chǔ)上，進(jìn)一步分析了 滲流應(yīng)力耦合特性。建立滲流場與應(yīng)力場耦合的等效連續(xù)介質(zhì)模型，進(jìn)行計算分 析，同時與不耦合情況分別進(jìn)行對比，總結(jié)得到了巖體中滲流場與應(yīng)力場的相互 作用關(guān)系。具體研究內(nèi)容和方法如下： （1） 從試驗(yàn)、 理論模型和耦合模型計算方法三方面簡要闡述國內(nèi)外學(xué)者對滲 流、應(yīng)力及兩者之間的耦合關(guān)系的研究成果。 （2）介紹滲流場與應(yīng)力場耦合分析的基本原理和方程。 （3） 編制 Matlab 與 COMSOL Multiphysics 的鏈接程序， 通過導(dǎo)入控制坐標(biāo)， 建立拉西瓦拱壩和壩基的三維有限元數(shù)值模型。 （4）利用多物理場全耦合分析軟件 COMSOL Multiphysics 對拉西瓦拱壩在 不考慮滲流應(yīng)力耦合的情況下的滲流場和應(yīng)力場進(jìn)行計算分析。 （5） 確定滲流應(yīng)力耦合模型的相關(guān)參數(shù)和邊界條件， 并對正常蓄水情況下壩 體和壩基的耦合滲流場和應(yīng)力場進(jìn)行計算分析，研究耦合條件下滲流場和應(yīng)力場 的變化規(guī)律。 （6） 比較考慮滲流應(yīng)力耦合和不考慮滲流應(yīng)力耦合的計算分析成果， 得出滲 流場與應(yīng)力場的相互作用關(guān)系。 拉西瓦拱壩壩基滲流應(yīng)力耦合分析 12 1.4 技術(shù)路線 圖 1.5 技術(shù)路線 利用 COMSOL 建立數(shù)值模型 國內(nèi)外研究進(jìn)展 試驗(yàn)研究成果 理論模型研究成果 計算方法研究成果 工程應(yīng)用前景 應(yīng)力場 基本方程 數(shù)學(xué)模型的建立 對比分析 得出結(jié)論 拉西瓦工程 基本資料 計算考慮耦合情況的 滲流場和應(yīng)力場 計算不考慮耦合情況 的滲流場和應(yīng)力場 收集相關(guān)文獻(xiàn)資料、 工程實(shí)例及分析 滲流場 基本方程 耦合場 基本方程 水利水電工程專業(yè)畢業(yè)論文 13 第二章 滲流應(yīng)力耦合分析數(shù)學(xué)模型 壩基巖體是被水飽和的。一方面，大壩基坑抽水、壩基開挖、水庫蓄水等必 然在巖體中產(chǎn)生非穩(wěn)定滲流場，該滲流場的作用將影響巖體的受力和變形；另一 方面，巖體應(yīng)力場的變化將引起節(jié)理裂隙寬度的變化，從而引起巖體滲透系數(shù)的 變化，并影響滲流場。因此本章根據(jù)拉西瓦工程的實(shí)際情況和地質(zhì)構(gòu)造的特點(diǎn)， 介紹滲流場分析和滲流場與應(yīng)力場耦合分析的數(shù)學(xué)模型。 2.1 基本假設(shè) 在進(jìn)行滲流應(yīng)力耦合分析時，需要對巖體作如下簡化假設(shè)： （1）巖體是分區(qū)均質(zhì)的，各向同性的連續(xù)體； （2）巖體是飽和的，其中水的滲流運(yùn)動符合達(dá)西定律； （3）巖體骨架變形是微小的； （4）巖體礦物不可壓縮，即巖體變形主要來源于隙變形； （5）慣性力影響忽略不計。 2.2 滲流場與應(yīng)力場耦合分析數(shù)學(xué)模型 2.2.1 應(yīng)力場的基本方程 在巖體中截取一微分單元，其受力如圖 2.1 所示。 圖 2.1 微分單元體受力平衡 假如體積力僅考慮重力，并規(guī)定應(yīng)力符號以拉應(yīng)力為正，則由微元體的平衡 可得平衡方程如下： dz z z z dz z zy zy dx x xz xz dx x xy xy dx x x x r x zy z zx dz z zx zx xz xy o z y x 拉西瓦拱壩壩基滲流應(yīng)力耦合分析 14 0 zyx xz xy x （2-1a） 0 zyx yzyxy （2-1b） 0 r z yz xz zyx （2-1c） 式中： ij 為總應(yīng)力的張量，i、j表示x、y、z； r 為巖體的飽和容重。 根據(jù)太沙基（Terzaghi）有效應(yīng)力原理，有 p （2-2） 式中： 為總應(yīng)力， T zxyzxyzyx ； 為有效應(yīng)力， T zxyzxyzyx ； p為空隙水壓力， T pppp000。 巖體的壓力應(yīng)變關(guān)系表示為 0 D （2-3） 式中： 為應(yīng)變， T zxyzxyzyx ； 為考慮非線性時巖體的初應(yīng)變； D為彈性矩陣或彈塑性矩陣。 由幾何方程，應(yīng)變與位移的關(guān)系為 uB （2-4） 式中： u為位移， T zyx uuuu ； B為幾何矩陣，且 水利水電工程專業(yè)畢業(yè)論文 15 T xyz zxy zyx B 000 000 000 （2-5） 為便于推導(dǎo)，下面暫取彈性本夠關(guān)系。將式(2-2)、(2-3) 、(2-4)代入式（2-1） ， 便可得以位移分量和空隙水應(yīng)力表示的平衡微分方程，即 0 21 02 X x p x G uG v x （2-6a） 0 21 02 Y y p y G uG v y （2-6b） r v z Z z p z G uG 02 21 （2-6c） 式中：G為剪切模量， )1 (2 E G； 2 為拉普拉斯算子， 2 2 2 2 2 2 2 zyx ； 為體積應(yīng)變， z u y u x u z y x ； 000 ,ZYX為由初應(yīng)變 0 引起的等價體積應(yīng)力。 一般情況下，用總水頭（即水力勢）表示空隙水應(yīng)力比用空隙水應(yīng)力表示更 具普遍意義，因此，將式（2-6）中的孔隙水應(yīng)力改用總水頭表示。總水頭h與空 隙水應(yīng)力p之間的關(guān)系為 )(zhp w （2-7） 式中： w 為水的容重； z為位置高程。 總水頭h用h w 表示，并仍記為h，則 拉西瓦拱壩壩基滲流應(yīng)力耦合分析 16 w z h z p y h y p x h x p （2-8） 將式（2-8）代入式（2-6） ，便可得用位移分量和總水頭表示的平衡方程，即 0 21 02 X x h x G uG v x （2-9a） 0 21 02 Y y h y G uG v y （2-9b） wr v z Z z h z G uG 02 21 （2-9c） 2.2.2 滲流場的基本方程 根據(jù)質(zhì)量守衡定律，水在飽和微元體中的增減速率等于水進(jìn)出該微元體流量 速率之差，如圖 2.2 所示，即 Vn t dxdydzv z v y v x wzwywxw （2-10） 式中： w 為水的密度； x v、 y v、 z v為滲流速度； n為巖體的孔隙率； V為微元體的體積，dxdydzV 。 考慮水的壓縮性。設(shè)水的壓縮系數(shù)為，它是水的彈性模量 w E 的倒數(shù)，則式 （2-10）右端可表示為 dxdydz t h n t Vn t w V ww （2-11） 由假設(shè)(2),水的的滲流運(yùn)動服從達(dá)西（Darcy）定律,即 w hgradkv （2-12） 式中： k為滲透張量，其主值為 x k， y k， z k。 水利水電工程專業(yè)畢業(yè)論文 17 圖 2.2 微分單元體流量平衡 將式（2-11）和（2-12）代入（2-10） ，即可得流體的連續(xù)性方程。為簡便計， 滲透主向量坐標(biāo)系表示，仍記為x，y，z，則得 t h n tz h k zy h k yx h k x w v wzyx （2-13） 其中，總水頭h的因次是帕(Pa)，即與應(yīng)力的因次相同。 2.2.3 耦合場的基本方程 巖體內(nèi)的應(yīng)力狀態(tài) ij 決定了巖體的滲透張量 ij k，即)( ijij fk；反過來，滲 透張量 ij k又通過滲透動水壓力 d F 影響著巖體的應(yīng)力狀態(tài)， 也即)( ijij kg。 因此， 巖體應(yīng)力場與滲流場是相互影響的，或兩者是耦合的，需要反復(fù)迭代才能使耦合 應(yīng)力場達(dá)到穩(wěn)定。 式（2-9）和式（2-13）聯(lián)立，即組成應(yīng)力場和滲流場耦合分析的基本方程， 即著名的比奧方程。其邊界條件包括： （1）位移邊界條件 0 uu (2-14a) （2）應(yīng)力邊界條件 0klklm 3 , 2 , 1, lk (2-14b) （3）水頭勢邊界條件 o z y x Vz Vx Vy 拉西瓦拱壩壩基滲流應(yīng)力耦合分析 18 0 hh (2-14c) （4）流量邊界條件 0 1 q n h kn w (2-14d) 式中： l m 為邊界外法線的方向余弦。 在適當(dāng)?shù)倪吔鐥l件下求解方程，便可得到滲流場和應(yīng)力場。 關(guān)于滲透系數(shù)的變化，一般假設(shè)為應(yīng)力狀態(tài)的指數(shù)函數(shù)。假設(shè)滲流主向與應(yīng) 力主向重合，則主滲透系數(shù)取以下形式 )exp(00 0)exp(0 00)exp( 3 2 1 0 kk (2-15) 式中： 0 k 為巖體表面滲透系數(shù)； 為影響系數(shù)，由實(shí)驗(yàn)確定； i 為有效主應(yīng)力，3 , 2 , 1i。 設(shè)主應(yīng)力與整體坐標(biāo)之間的轉(zhuǎn)換矩陣為 T，則在整體坐標(biāo)下的滲透張量為 TkTk T （2-16） 對于滲流對介質(zhì)應(yīng)力的影響，考慮滲流體積力F作用于介質(zhì)上，即 )(gradpF （2-17） 由式（2-7） ，并注意到總水頭h的因次，可得 )(grad)grad(zhF w （2-18） 其中，)(grad zF ws 是靜水壓力， 其作用效果可用有效應(yīng)力考慮， 而)(grad hFd 是動水壓力，計算時作為外荷載考慮。 2.3 COMSOL Multiphysics 軟件 2.3.1 COMSOL Multiphysics 軟件簡介 COMSOL Multiphysics 是一款基于有限元計算方法的大型的高級數(shù)值仿真軟 件，由瑞典的 COMSOL 公司開發(fā)，被當(dāng)今世界科學(xué)家稱為“最專業(yè)的多物理場 水利水電工程專業(yè)畢業(yè)論文 19 全耦合分析軟件”37。COMSOL Multiphysics 專注于解決多物理場耦合問題，以 高度的開放性、杰出的易用性和優(yōu)秀的計算性能實(shí)現(xiàn)了任意多物理場高精確的數(shù) 值仿真，廣泛