基于橫觀各向同性模型的片巖隧洞穩定性研究

孫 濤， 楊渝南， 高素芳， 劉 杰,2， 謝曉康， 韓紹康

(1.三峽大學 土木與建筑學院， 湖北 宜昌 443002； 2.三峽大學 三峽庫區地質災害教育部重點實驗室， 湖北 宜昌443002)

1 研究背景

片巖在深部高地應力區域大量存在，可近似為層厚極小的層狀巖體，其變形和強度具有明顯的各向異性特性，這一特點又會對工程巖體穩定性產生特殊影響[1]。在室內試驗方面，劉志勇等[2]利用MTS815巖石試驗系統研究石英云母片巖在加、卸荷應力路徑下的力學特性；周新等[3]探討片巖加速蠕變階段的軸向應變-應變率-偏應力關系特征；楊以榮等[4]系統研究高圍壓卸荷條件下石英云母片巖變形破壞特征及能量演化規律；何忠明等[5]通過室內試驗和數值模擬分析層狀巖體單軸壓縮情況下的應力應變響應以及強度特征；于德海等[6]對干燥和飽水狀態下的綠泥石片巖進行三軸壓縮試驗，獲取水對試樣強度和變形特性的影響規律；劉廷等[7]通過室內三軸試驗，研究不同含水率狀態下云母片巖的強度和變形特性；王章瓊等[8]研究片巖在片理面傾角分別為0°、45°、90°時的變形參數各向異性特征及破壞機制；鄭松林等[9]通過片巖在單軸應力作用下的巖石力學特性試驗研究,得出巖石試塊的強度特征、變形特征與作用應力的大小、作用應力與結構面之間夾角有密切聯系的一般規律；張曉平等[10]研究了與片巖裂紋擴展相關的應力門檻值。

針對強度理論和本構關系方面，肖明礫等[11]基于帶Kachanov 蠕變損傷因子的蠕變模型，研究丹巴水電站石英云母片巖的三軸蠕變特性及其各向異性特性；張春生等[12]以顆粒流方法為基礎,開展不同尺寸、不同結構組成的各向異性特征研究,揭示變形參數和強度參數的變化規律；文獻[13-17]針對層狀巖體彈性參數的各向異性進行了試驗與理論分析，并構建了相應的破壞準則，張學民等[18]研究了巖石材料的變形和強度各向異性性狀、各向異性破壞準則和彈塑性本構關系及各向異性特征對層狀圍巖隧道變形和穩定性的影響；張嘉威[19]根據單軸壓縮試驗及統計損傷理論,建立了針對單軸壓縮狀態下的巖石損傷本構模型。

在工程應用方面，已有許多結合大型工程的片巖研究成果，如Nasseri等[20]研究了片巖力學參數和變形破壞特征的各向異性特征；劉杰等[21]、許振華[22]及單治剛等[23]運用大型巖土分析軟件ADINA和單軸試驗闡明了丹巴石英云母片巖各向異性特征，且巖體的變形和破壞規律主要取決于地應力及巖體的結構特征；在基于隧道等含片巖巖體的大型工程方面，Hwang等[24]研究了各向異性的片理角度對隧道穩定性的影響；Shrestha等[25]指出當隧洞通過高度片狀化的巖層時，在塑性變形和剪切區域存在嚴重的隧道變形。

在上述各方面取得豐碩成果的同時，有以下幾方面需進一步深入研究：

(1)不同片理角度下巖體力學參數取值無試驗依據。很多學者在進行巖體參數的取值時往往按經驗進行參數折減或僅按數值反演來確定，人為主觀性太強，參數取值無試驗數據支持。

(2)現有研究未系統考慮巖樣在不同高度和不同片理角度下的巖樣抗拉強度尺寸效應的規律。

(3)現有研究多數采用的計算模擬范圍較大，在各邊界施加的應力邊界效應顯著，應力影響范圍較小，無法影響到內部局部精細區域，無法對局部區域進行精細的片理模擬和精細應力構造。

(4)現有研究未深入考慮圍巖片理與隧洞洞軸不同夾角下隧洞圍巖穩定性影響規律并給出最不利工況。

隨著我國西部大開發和可持續發展戰略的實施，西部區域的引水隧洞工程和交通隧洞工程正不斷涌現，針對地下工程中圍巖開挖后巖體各向異性顯著的研究越來越成為熱點。本文基于丹巴引水隧洞工程，通過對該工程中石英云母片巖的各項物理和力學性質進行試驗研究，并建立隧洞模型，從位移、塑性區、塑性應變值3個方面對隧洞開挖后圍巖穩定性進行評價，對于隧洞開挖加固設計有著重要的理論意義和工程參考價值。

2 工程概況

丹巴水電站位于四川省甘孜藏族自治州丹巴縣境內的大渡河干流上，引水隧洞位于大渡河左岸山體中，全長約16 710 m，垂直埋深基本為500～900 m，最大埋深可達約1 220 m，大部分洞室位于石英云母片巖軟巖區域，巖石堅硬程度初步分類屬較軟巖，巖體質量以Ⅳ類為主，成洞條件差，開挖后有可能出現較嚴重的大變形破壞現象，遇水也易發生軟化現象，如圖1所示。因此，丹巴水電站引水隧洞軟巖成洞問題突出，為樞紐開發的關鍵技術問題之一。

圖1 隧洞工程圍巖破壞圖

3 片理在不同空間角度下的參數測定

為獲得石英云母片巖在高地應力下三軸加卸載力學特性和變形破壞規律，本文開展了室內三軸壓縮試驗測定不同圍壓作用下巖樣的力學參數。

3.1 室內三軸試驗

3.1.1 三軸試驗方案 試驗以RMT-150C巖石力學試驗系統為平臺，三軸試驗采用3組不同片理傾角，4種圍壓循環加卸載。由于引水隧洞段實測最大主應力為36 MPa，故將圍壓分為0～10、0～15、0～20、 0～30、0～35、0～40 MPa共6級，軸向和圍壓同時加載，速率均為0.2 kN/s，每次循環2～4次，軸向應力線性增加直至巖樣破壞，不同傾角巖樣三軸加卸載試驗加載模式如圖2所示，部分巖樣三軸試驗破壞形態圖如下圖3所示。

圖2 不同片理傾角的三軸加卸荷試驗

圖3 部分巖樣三軸試驗破壞形態

3.1.2c、φ值各向異性的試驗分析 根據3組36個巖樣的三軸試驗成果，對石英云母片巖表現的特征繪制成圖，如圖4～5所示。

圖4 干燥巖樣彈性模量、黏聚力和內摩擦角與圍壓關系圖

圖5 飽和巖樣彈性模量、黏聚力和內摩擦角與圍壓關系圖

由圖4～5可知，隨著圍壓的增加，石英云母片巖的彈性模量、內摩擦角和黏聚力均與圍壓呈線性增大規律，但不同片理傾角巖樣的各力學參數對圍壓的敏感性有所差異。

(1)在干燥條件下，圍壓對巖樣力學參數影響分析：由圖4分析可知，①當軸向力與片理面平行時，隨著圍壓的增加巖樣參數彈性模量和內摩擦角線性增大速率最快，斜率分別為0.2245 和0.2873；②當軸向力與片理面斜交時巖樣黏聚力增大速率最快，斜率為0.0711；③根據斜率對比得出巖樣力學參數對圍壓敏感性排序為內摩擦角>彈性模量>黏聚力。

(2)在飽和條件下，圍壓對巖樣力學參數影響分析：由圖5分析可知，①當軸向力與片理面平行時，隨著圍壓的增加巖樣彈性模量、內摩擦角增大速率最快，斜率分別為0.135和0.2907；②當軸向力與片理面斜交時，巖樣黏聚力增大速率最快，斜率為0.0943；③根據斜率對比得出巖樣力學參數對圍壓敏感性排序為內摩擦角>彈性模量>黏聚力。

(3)飽和石英云母片巖相對干燥狀態的彈性模量降低20%～40%，黏聚力降低25%～40%，內摩擦角降低20%-30%。表明水對片巖的彈性模量、內摩擦角和黏聚力有顯著影響，在實際工程設計應用中應當考慮水對巖體質量的劣化影響。

3.1.3 各向異性與圍壓關系分析 由于片理的存在，石英云母片巖具有顯著的各向異性特征，給出表征各向異性公式(1)和(2)。

(1)

(2)

式中：σ1為最大主應力，Pa；σ3為最小主應力(由于為單軸壓縮，故σ3=0)，Pa；k1、k2分別為石英云母片巖的單軸抗壓強度和變形模量各向異性參數；E∥，E⊥分別為平行、垂直節理的變形模量。k1和k2等于1時說明巖石各向同性，當其越偏離1時，各向異性越明顯。

根據表征各向異性的公式(1)和(2)計算得出各向異性參數k1和k2，根據計算結果得出各向異性程度與圍壓關系圖，如圖6所示，在不同圍壓下對不同片理傾角的石英云母片巖的各向異性進行評價。

由圖6可知，在低圍壓下k1>1且k2<1,各向異性顯著，隨著圍壓的增加，k2→1，而后逐漸減小；k1逐步降低，各向異性減小，而后逐步升高，各向異性增加。巖體強度受片理結構影響較大，在實際工程設計應用中應當充分考慮巖樣的各向異性特性。

3.2 抗拉強度各向異性劈裂試驗分析

3.2.1 試驗方案 該試驗基于筆者自主研發的單個巖樣多次劈裂測巖石抗拉強度的試驗方法來測定不同片理產狀的石英云母片巖抗拉強度。試驗分為5組，每組4個巖樣(巖樣高度分別為2.5、3.5、4.5和5.5 cm)，具體試驗方案如表1所示，通過對破壞后巖樣進行重粘接，旋轉巖樣調整加載力與片理交角實現通過單個巖樣測不同片理角度下巖石抗拉強度，兩次劈裂后巖樣如圖7～8所示：

表1 劈裂試驗方案表

3.2.2 不同高度巖樣抗拉強度隨片理角度的變化分析 由圖9可知，石英云母片巖的抗拉強度受控于其片理特征，相同高度下巖樣的抗拉強度隨片理角度線性增大；相同片理夾角下，巖樣抗拉強度隨巖樣高度遞增呈良好線性遞增關系。徐燕飛等[26]研究結果顯示，各種巖石的劈裂強度都存在尺寸效應，不同巖石劈裂強度尺寸效應主要受巖石的物質成分和結構所影響；硅質砂巖劈裂強度隨試樣厚度增加呈線性增大；本文研究試驗結果與該結論觀點相一致。

石英云母片巖的抗拉強度與其試驗巖樣高度和片理角度有關，為研究不同片理角度下巖樣抗拉強度尺寸效應的規律，定義巖樣抗拉強度尺寸效應(此處尺寸效應實際指巖樣的高徑比)系數為ksize=σ5.5/σ2.5，其中σ5.5是高度為5.5 cm巖樣的抗拉強度，σ2.5是高度為2.5 cm巖樣的抗拉強度，得出巖樣抗拉強度尺寸效應系數ksize隨片理角度的變化關系圖，如圖10所示。

圖6各向異性程度與圍壓關系圖 圖7第1次劈裂后巖樣 圖8第2次劈裂后巖樣

圖9不同高度巖樣抗拉強度隨片理角度的變化關系 圖10抗拉強度尺寸效應系數隨片理角度的變化關系

由圖10可知，不同片理角度下巖樣抗拉強度尺寸效應系數ksize均大于2，即相同角度下5.5 cm高度巖樣測得的抗拉強度一般為2.5 cm高度巖樣抗拉強度的2倍以上，可見巖樣的抗拉強度具有顯著的尺寸效應，且巖樣抗拉強度尺寸效應與片理角度呈線性遞減的關系。

3.3 橫觀各向同性本構模型研究

根據前期室內試驗成果和現場圍巖破壞特征的分析，石英云母片巖的變形與破壞均表現出顯著的橫觀各向同性特征。

根據橫觀各向同性理論可知，在局部坐標系下，增量型彈性應力-應變的關系為 ：

[Δσ′]=[K′][Δε′e]

(3)

式中：σ′為巖體彈性應力；ε′e為巖體線應變, [K′]為局部剛度矩陣，即:

(4)

其中：

式中：E1和μ12分別為各向同性平面的彈性模量和泊松比；E3和μ13分別為垂直與同性平面的彈性模量、泊松比。

通過對不同片理角度的巖樣進行室內巖石物理力學試驗(圖3、7、8)，得到平行片理和垂直片理方向巖樣的力學參數，如表2，將其作為片理模型的計算參數，為公式(4)提供參數支撐。

表2 巖體計算參數表

4 計算參數及模型

4.1 計算參數

采用橫觀各向同性本構模型，在ADINA軟件中通過空間兩點坐標來定位材料主軸，通過改變空間兩點坐標來調節材料軸的方位，如圖11，材料主軸z方向取垂直于片理參數，x和y軸均以平行片理參數選取，在測點地應力實測與計算值(表3)的基礎上得出空間任意角度的巖體片理模型的計算參數，實現各種片理類型及角度的模擬。

4.2 引水隧洞軟巖段地質體大區域三維模型

沿著引水隧洞方向建立三維模型，如圖12所示，模型順引水洞方向取10 358.215 m，垂直引水洞方向取2 066.367 m，根據地形，最大相對高度3 750 m，選取實測地應力測點位置為ZK2和ZK3，ZK2和ZK3測點剖面位置如圖13所示。

采用水壓致裂法測取測點ZK2和ZK3的地應力實測值，通過調節所建立大區域三維模型邊界應力的大小，反演分析得出與實測值相接近的反演計算值，如表3。

4.3 主應力構造精細應力場方法的提出

4.3.1 大區域模型反演分析存在的不足 采用以上大區域模型反演分析可較好地得出引水隧洞大區域的整體地應力場，但大區域模型反演分析存在以下不足：①隧洞開挖各工況下圍巖片理角度在不斷改變，無法對局部區域進行精細的片理模擬和精細應力構造；②大區域模型只是在各邊界面上施加應力，應力邊界效應顯著，應力影響范圍較小，無法影響到內部局部精細區域。因此本文對大模型局部區域建立小模型，進行精細計算模型應力場構造。

表3 測點地應力實測與計算值 MPa

建立小模型后，若直接采用測得各正應力、剪應力在各邊界上施加應力，各邊界面均有3個不同的應力，共18個應力，構造的應力場空間受力復雜，無法對局部區域進行精細應力場精準構造，空間單元體的受力示意圖如圖14所示。

4.3.2 精確主應力構造精細應力場的方法 本文基于此提出一種新的主應力構造精細應力場方法，提出的《一種精確主應力構造精細應力場的裝置及方法》(201811096241.5)發明專利已公示[27]，且該方法在丹巴工程中得以成功運用。采用實際精確的主應力方式進行構造，首先根據ZK2測點實測地應力的6個應力分量σx、σy、σz、σxy、τyz、τzx可得出3個不變量I1，I2，I3。

I1=σx+σy+σz

(5)

I2=σxσy+σyσz+σzσx+2τxy+2τyz+2τzx

(6)

I3=σxσyσz+2τxyτyzτzx-2σxτyz-

2σyτzx-2σzτxy

(7)

由公式(5)～(7)推導可得出主應力方程為：

σ3-I1σ2-I2σ1-I3=0

(8)

運用Matlab軟件，可求出三元一次方程組(8)，得到各主應力值大小，根據主應力值的大小，由公式(9)～(10)聯立求解后按公式(11)計算主應力傾角αi和方位角θi，測點ZK2各主應力、方位角和傾角如表4。

(9)

(10)

公式(9)中li、mi、ni(i=1、2、3)分別為主應力σi對于X、Y、Z軸的方向余弦。

(11)

選取局部小區域，根據ZK2測點各主應力傾角和方位角建立隧洞模型，在模型邊界上施加與之垂直的相應面力，其大小為3個主應力值，即為精細計算模型應力場構造，如圖15所示。該方法用實際精確的3個主應力進行精細應力場精準構造，大大優化了空間多個應力構造應力場的復雜性并提高了精準度。

表4 測點ZK2各主應力、方位角和傾角計算值

根據表4，建立尺寸為240 m×240 m×240 m的立方體隧洞模型，隧洞斷面為圓形，直徑為16 m，片理層面產狀取N65°、WNE∠80°，巖體強度參數和變形參數參照表2進行取值，對開挖后隧洞圍巖穩定性進行評價。

5 隧洞開挖后圍巖穩定性評價

5.1 空間多角度片理與隧洞軸向夾角關系說明

由于施工現場開挖各工況下圍巖片理角度在不斷改變，使得圍巖片理與隧洞軸向夾角也在不斷改變，施工現場實地考察結果顯示，片理與洞軸方向水平面夾角的度數分布集中在45°～90°之間，局部地段也可能出現其他類型及各種角度穿切情況。由此設定了3類片理走向類型及0°、15°、45°、75°、90°片理進行模擬計算，如圖16所示，假定洞軸方向為x方向，定義第1類45°片理面為與平面ACHF平行，第2類45°片理面為與平面ABHG平行，第3類45°片理面為與平面BCGF平行，其他角度亦按照此進行分類，以此對空間多角度片理與隧洞軸向夾角的關系進行說明。

5.2 空間多角度片理產狀下隧洞圍巖穩定性分析

本文以洞長36 m，洞徑8 m的試驗洞為例，研究洞室圍巖的穩定性隨片理產狀的變化規律，從而預測最不利片理隧洞工況和隧洞失穩破壞可能性最大的位置。通過調整材料主軸來實現3類片理走向及0°、15°、45°、75°、90°片理進行模擬計算，下面主要從應力、位移、塑性區3個方面對開挖后隧洞圍巖穩定性進行分析，計算云圖中位移和應力數據均為以開挖隧洞為中心點的徑向位移和徑向應力，其方向以指向隧洞圓心為正，反之為負。

圖11空間任意角度巖體片 圖12引水隧洞軟巖段大區域 圖13ZK2和ZK3所在y向

理模型示意圖 三維模型圖x=1671m處的位置圖

圖14空間單元體的受力示意圖 圖15主應力構造精細應力場構造圖 圖16片理的假設走向

5.2.1 應力分析 隧洞開挖后洞室圍巖處于受壓狀態，各工況下洞室周圍應力分布基本一致，主要集中于洞室的左下方拱腳和右上方拱肩，水平方向應力均較小，沒有明顯的變化規律。

由計算結果得出不同片理工況下應力變化趨勢圖，如圖17所示。

由圖17分析可知，不同片理角度工況下圍巖應力隨著片理角度總體趨勢是先遞增而后遞減。第2類和第3類不同片理產狀的隧洞圍巖的應力變化規律較為一致，均在片理角度為15°時，隧洞圍巖應力值達到最大，分別為25.81和23.07 MPa，第1類片理變形規律略有差異，在片理角度為45°時，隧洞圍巖應力值達到最大，為24.14 MPa，表明工程中3類片理走向類型在片理角度分別為15°、15°、45°時為最不利工況，圍巖應力對隧洞穩定性影響較大，應重點加固并監測該工況下隧洞圍巖的變形情況。

5.2.2 位移分析 根據計算結果統計得出3類片理角度位移云圖對比表，如表5所示。

由表5分析可知，第1類和第2類不同片理產狀的隧洞圍巖的變形規律較為一致，其最大變形多發生在左側拱肩和右側拱腳，第3類片理不同片理傾角下，變形規律略有差異，45°時最大變形位于拱頂和拱底。工程中應針對變形較大的區域重點加固監測。

表5 3類片理角度位移云圖對比表

由計算結果得出不同片理工況下位移值隨片理角度變化圖，如圖18所示。

由圖18分析可知，3類片理走向類型均在片理角度為90°時，隧洞圍巖位移值達到最大，分別為10.88、10.88、9.72 cm；表明工程中3類片理走向類型均在片理角度為90°時為最不利工況，隧洞圍巖會出現較大變形或失穩，應著重考慮隧洞開挖的安全性。

5.2.3 塑性區分析 根據計算結果統計得出3類片理角度塑性區云圖對比圖，如表6所示。

對表6各工況進行分析可知，開挖卸荷后3類各角度片理均出現較大塑性區，塑性區開展部位均類似，其主要集中于隧洞左側拱底和右側拱肩處。多數片理角度工況下，隧洞開挖后左側拱底塑性區開展深度較右側拱肩塑性區略大，僅在第1類45°和第3類45°工況下，右側拱肩塑性區深度大于左側拱底。

由計算結果得出不同片理工況下圍巖塑性區深度和面積變化趨勢圖，如圖19和20所示。

表6 3類片理角度塑性區云圖對比表

圖17不同片理工況下應力變化趨勢圖 圖18不同片理工況下位移變化趨勢圖

圖19不同片理工況下圍巖塑性區深度變化趨勢圖 圖20不同片理工況下圍巖塑性區面積變化趨勢圖

由圖19～20分析可知，不同片理角度工況下的圍巖塑性區深度及面積差異很大，但總體趨勢是隨著片理角度的增加圍巖塑性區深度和面積先遞增而后遞減。

(1)圍巖塑性區深度及面積分析：①第1類片理角度隧洞開挖后，45 °時塑性區開展深度最大，為8.2 m，15°時塑性區面積最大，為139.7 m2；②第2類片理角度隧洞開挖后，45 °時塑性區開展深度最大，為8.3 m,且45°時塑性區面積也最大， 為200.1 m2；③第3類片理角度隧洞開挖后，45 °時塑性區開展深度最大，為9.0 m，75°時塑性區面積最大,為188.5 m2。

(2)從圍巖塑性區深度和面積分析可得出：①3類片理均在45°時塑性區深度最大，為最不利片理角度；②不同類片理塑性區面積無明顯規律性，第1類在15°時為最不利工況，第2類在45°時為最不利工況，第3類則在75°時為最不利工況。

6 結 論

(1)隨著圍壓的增加，石英云母片巖的彈性模量、內摩擦角和黏聚力均與圍壓呈線性增大的關系，但不同片理傾角巖樣的各力學參數對圍壓的敏感性有所差異。在干燥和飽和兩種條件下，巖樣力學參數對圍壓敏感性排序均為內摩擦角>彈性模量>黏聚力。

(2)飽和石英云母片巖相對干燥狀態的彈性模量降低20%～40%，黏聚力降低25%～40%，內摩擦角降低20%～30%。表明水對片巖的彈性模量、內摩擦角和黏聚力有顯著影響，在實際工程設計應用中應當考慮水對巖體質量的劣化影響。

(3)隨著圍壓的增加，石英云母片巖的各向異性程度呈線性遞減趨勢；且石英云母片巖的抗拉強度具有顯著的尺寸效應。

(4)本文提出一種新的主應力構造精細應力場方法，避免了大計算模型的邊界效應，可對局部區域進行精細的片理模擬和精細應力構造。

(5)根據片理走向及不同傾角選定15個不同片理空間展布模型進行計算分析，分別得出應力、位移、塑性區隨片理角度的變化規律，并指出不同片理展布條件下的最不利工況。