熊子正 孟慶祥 馮 欣 滕志強
(1. 河海大學 巖土力學與堤壩工程教育部重點實驗室, 南京 210098; 2. 河海大學 江蘇省巖土工程技術工程研究中心, 南京 210098)
柱狀節理巖體是玄武巖的一種原生構造.它通常將玄武巖切割成六棱柱狀或其他不規則的棱柱狀.R.Mallet[1]在1875年開始對玄武巖柱狀節理特征和形成機制進行研究.近年來,國內外學者就玄武巖柱狀節理開展了專題研究.V.Moon和J.Jayawardane[2]對玄武巖風化特性展開了研究,國內黃國明等學者[3-4]對一般玄武巖的力學與變形特性進行了研究.根據巖體的空隙結構,變形機理和滲流特性,裂隙巖體滲流應力耦合模型分為連續介質模型、離散裂隙模型,雙重介質模型和多重裂隙模型.
連續介質模型耦合分析將裂隙巖體看作等效的多孔介質,基于比奧固結原理建立滲流應力的耦合關系[5],該方法適用于裂隙較密集的巖體.對于裂隙比較稀疏而巖塊透水性較大的巖體,可以用雙重介質模型或多重裂隙網絡模型進行耦合分析,根據巖體裂隙結構變化建立耦合關系.對于裂隙比較稀疏、巖塊滲透性很低,滲流主要沿裂隙網絡做定向移動的巖體,只能用離散裂隙網絡模型進行分析.目前,連續介質滲流應力耦合分析已相對比較成熟,而非連續介質方面的研究正處于發展階段.實際工程巖體中的相當部分應看作離散裂隙網絡巖體,利用非連續介質裂隙網絡模型解決工程實際問題是巖體水力學發展的重要方向[6-10].本文以模擬的節理巖體裂隙網絡為研究對象,基于離散單元法,分析柱狀節理巖體滲流應力耦合特性.
關于玄武巖柱狀節理形成機制有多種假說,目前為國內外學者所認可是“冷卻收縮說”.該學說認為,在巖流冷卻的過程中,熔巖冷凝面形成了規則且間隔排列的均勻收縮中心,產生了與收縮方向垂直的張力裂隙,體積收縮導致巖石向收縮中心聚集,巖石裂開,形成了多面柱狀節理.由于其地質成因特殊,所以地質特點與常見的巖體也有較大的差異,其變形和強度均具有較為顯著的各向異性.
根據柱狀節理的發育情況,可以將柱狀節理玄武巖分為三類.一類柱狀節理發育較好,柱體細長.二類柱狀節理發育不完全,尚未切割成完整的柱體.三類柱狀節理發育程度最差,對巖體工程特性影響較?。趯嶋H工程中,主要研究一類柱狀節理的特性,三類柱狀節理可按普通玄武巖對待.巖石的條件、巖體結構面條件、巖體的應力環境是影響巖體特性的三個主要因素.
室內試驗測得,玄武巖柱狀節理巖石塊體密度為2.90 g/cm3,顆粒密度為2.93 g/cm3,巖塊自然狀態下變形模量平均值為65.1 GPa,飽和條件下變形模量平均值為51.6 GPa[7].
柱狀節理玄武巖屬于堅硬巖,具有較高的強度與剛度.柱狀節理巖體結構面發育特征和巖體受力條件密切相關,同時二者也是柱狀節理玄武巖基本力學特征的重要影響要素.柱狀節理主要為不規則的多邊形(以四、五、六邊形為主)細長柱體.豎直向柱狀節理密集發育,與水平向節理相互切割,并常常存在交錯隱節理.因此柱狀節理巖體表現出比較規律的各向異性特征,且巖體完整性不好,是一種工程和力學特性較差的特殊節理巖體.
3DEC是ITASCA公司開發的基于離散單元法描述離散介質力學行為的三維計算分析程序.3DEC可以分析離散元介質在荷載,流體,溫度等作用下的靜態與動態響應問題,適合研究大型高邊坡穩定變形機理,地下深埋工程圍巖破壞,巖體結構裂隙滲流特性,還能解決多場耦合的問題.3DEC中滲流應力耦合的機理是:裂隙滲流產生的滲流力導致應力場重新分布,而應力場的變化又引起了裂隙隙寬的變化,使得裂隙巖體的滲透系數發生變化,從而導致巖體滲流場的重新分布.
本模型為隨機產生裂隙網絡的柱狀節理模型,模型長250 m,寬250 m,高100 m.巖體模型及裂隙網絡如圖1,圖2所示.根據地質資料及相關文獻,巖體及結構面采用的模型材料參數見表1,表2.

圖1 巖體模型 圖2 巖體初始裂隙網絡圖

表1 巖體模型的材料參數

表2 裂隙網絡的材料參數
1)初始條件:在巖體裂隙網絡上,從上邊界施加孔應力5×106MPa,均勻遞減至下邊界為0 MPa,以此來模擬水力梯度.
2)邊界條件:對模型外邊界采用位移約束條件.
對柱狀節理巖體離散元模型,分別采用耦合與非耦合算法,對比分析模型的主應力與位移分布.耦合算法與非耦合算法所得最大主應力分布圖如圖3~4所示;耦合與非耦合算法所得的位移云圖如圖5~6所示.
1)應力計算結果分析
耦合算法所得的最大主應力為 4.17 MPa,非耦合算法所得的最大主應力為 -2.01×105Pa.耦合算法所得最大主應力的絕對值要大于非耦合算法的.

圖3 耦合算法下的最大 圖4 非耦合算法下最大主應力分布 主應力分布
2)位移計算結果分析
耦合算法所得的最大位移為4.49×10-5m,非耦合算法所得的最大位移為4.43×10-5m.耦合算法所得的最大位移絕對值要大于非耦合算法的.

圖5 耦合算法下位移梯度分布 圖6 非耦合算法下位移梯度分布
綜上可知,耦合算法所得的最大應力與位移絕對值都要大于非耦合算法的,因此采用應力耦合算法是偏于安全且更符合工程實際的.
對已有的柱狀節理離散元模型,分別進行1 000時步,9 000時步,13 000時步17 000時步的耦合計算,對比分析所得孔應力及流量的變化.計算結果如圖7~14所示.
1)1 000時步計算結果分析
由圖7~8可知,1 000時步下孔應力與流速基本按照滲流方向梯度減小,滲流與應力耦合不明顯.

圖7 1 000時步下的孔應力分布 圖8 1 000時步下的流速分布
2)9 000時步計算結果分析
由圖9~10可知,9 000時步計算下,部分裂隙孔應力增大,部分孔應力減小,甚至為0.這說明滲流場與應力場已經產生耦合作用,部分裂隙已經閉合,此時裂隙中沒有水,因此孔應力與流速均為0;而部分裂隙的孔應力已經超過了初始孔應力,說明滲流場作用于應力場,使得空隙張開,孔應力增大,流量也增大,應力場又作用于滲流場,如此循環,動態地變化.

圖9 9 000時步下的孔應力分布 圖10 9 000時步下的流速分布
3)13 000時步計算結果分析
由圖11~12可知,相比于9 000時步的計算結果,孔應力與流速有明顯的減小,孔應力分布又呈現梯度遞減的趨勢,說明滲流應力耦合開始趨于穩定.

圖11 13 000時步下的孔應力分布 圖12 13 000時步下的流速分布
4)17 000時步計算結果
由圖13~14可知,與13 000時步計算結果相比,孔應力與流速的大小與分布均變化不大,最大孔應力等于初始孔應力5×106MPa,說明耦合已經基本穩定.

圖13 17 000時步下的孔應力分布 圖14 17 000時步下的流速分布
本文以柱狀節理巖體的三維裂隙網絡為研究對象,基于離散單元法以及裂隙網絡滲流原理,分析柱狀節理巖體在應力滲流耦合的作用下的變化特性.
通過對比耦合算法與非耦合算法的計算結果,發現耦合算法下的最大主應力與位移梯度均大于非耦合算法的計算結果.這說明考慮應力滲流耦合是偏于安全的,這也證明采用耦合算法分析柱狀節理巖體的滲透性的正確性與必要性.通過分析不同時步的孔應力與流速云圖,可分析應力滲流的耦合特性,柱狀節理中裂隙的滲流力,改變了應力場的分布,應力場的變化使得裂隙的隙寬發生變化,隙寬的變化又引起了滲透系數的變化,從而滲流場重新分布.裂隙隙寬,孔應力與滲流量相互影響,相互依賴,它們的變化是一個動態的過程.
本文基于離散單元法,構建三維離散元數值模型,使用3DEC的滲流模塊,研究了柱狀節理巖體滲流應力耦合作用特性.算例計算結果說明所用方法的正確性,同時也證明使用3DEC等軟件建立三維離散元模型來分析柱狀節理巖體滲流問題的可行性.
[1] Mallet R. On the Origin and Mechanism of Production of the Prismatic(or Columnar) Stmcture of Basalt[J]. Philosophical Magazine,1875, 4(50): 122-135, 201-226.
[2] Moon V, Jayawardane J. Geomechanical and Geochemical Changes During Early Stages of Weathering of Karamu Basalt, New Zealand[J]. Engineering Geology, 2004, 74(1/2): 57-72.
[3] 黃國明,黃潤秋.某壩址玄武巖巖體強度和變形特性[J].長江科學院院報,1998,15(6):19-22.
[4] 王 毅,楊建宏.玄武巖的巖體結構與力學性狀研究[J].巖石力學與工程學報,2002,21(9):1307-1310.
[5] 王 媛,速寶玉,徐志英.等效連續裂隙巖體滲流與應力全禍合分析[J].河海大學學報,1998,26(2):26-30.
[6] 張有天.巖石水力學與工程[M].北京:中國水利水電出版社,2005.
[7] 石安池,唐鳴發,周其?。鹕辰Q灘水電站柱狀節理玄武巖巖體變形特性研究[J].巖石力學與工程學報,2008,27(10):2079-2079.
[8] 張彥洪,柴軍瑞.巖體離散裂隙網絡滲流應力耦合分析[J].應用基礎與工程科學學報,2012,20(2):253-262.
[9] 仵彥卿.裂隙巖體應力與滲流關系研究[J].水文地質工程地質,1995(6):30-35.
[10] 盛金昌,速寶玉.裂隙巖體滲流應力耦合研究綜述[J].巖土力學,1998,19(2):92-98.