周　利，蘇　凱，周亞峰，文喜雨，伍鶴皋

（1.武漢大學 水資源與水電工程科學國家重點實驗室，湖北 武漢　430072；2.長江勘測規(guī)劃設計研究院，湖北 武漢　430010；3.長江巖土工程總公司，湖北 武漢　430010）

1　研究背景

隨著社會經濟與工程技術的發(fā)展，水工隧洞無論在數量上還是在規(guī)模上都達到了前所未有的高度，帶來了巨大的社會與經濟效益，同時也面臨諸多問題。例如：隧洞在內水壓力作用下，襯砌混凝土受拉開裂，并在施工缺陷以及溫度微裂紋等影響下，形成性態(tài)不一的宏觀裂縫，進而導致內水外滲［1-2］，此類襯砌結構均稱為透水襯砌。在高內壓下襯砌環(huán)向拉應力超過混凝土抗拉強度后，襯砌混凝土損傷破壞，表現為微裂紋萌生、擴展并貫通，直至產生宏觀裂縫，洞內高壓水沿裂縫通道外滲使得滲透水壓增大，在影響結構應力場分布的同時又加劇了損傷，如此相互影響，是一個復雜的多場耦合過程。卞康等［3］采用增量變彈性損傷剛度計算迭代的方法分析了襯砌水壓致裂過程中的耦合效應；Simanjuntak等［4］探討了襯砌開裂過程中水力耦合效應對隧洞滲流量的影響；Bian等［5］在分析實際工程時，提出了一種滲流-應力耦合方法；Fahimifar等［6］在平面軸對稱假定下，提出了一種解析方法以分析襯砌開裂過程中的耦合效應。然而從我國進行的6次較大規(guī)模的隧洞襯砌直接壓水試驗結果來看（詳見表1），當隧洞內壓增大到某一數值時，裂縫出現，此后無新的裂縫產生，只是縫寬隨內壓增大，裂縫呈少而寬分布［7］，而這一現象在現有的研究中卻鮮有體現?？紤]到水工隧洞襯砌裂縫分布的這一特點，劉秀珍［8］按照彈性地基曲梁結構推導了圓形有壓隧洞襯砌裂縫寬度的計算方法；沈威［9］采用室內水工隧洞結構模型試驗研究了襯砌開裂過程，同樣得出了裂縫少而寬的分布特征。

表1　隧洞壓水試驗參數及結果［7］

襯砌裂縫出現后，洞內高壓水勢必會沿著裂縫向外流動，使得圍巖變形增大，當襯砌因鋼筋的約束作用徑向變形與圍巖不一致時，兩者可能出現脫離現象，單獨承擔自身范圍內的水壓力［10-11］，體現了襯砌與圍巖典型的有條件聯(lián)合承載機理［12-13］。此外，工程監(jiān)測資料顯示，隧洞發(fā)電運行時，襯砌鋼筋應力普遍小于設計預期，并沒有達到既安全又經濟的效果，說明現有的襯砌設計方法仍存在一定的問題，以致高估了鋼筋應力值。例如：在廣蓄電站高壓水工隧洞充水加壓初期，鋼筋應力隨著內壓增大；襯砌開裂后，鋼筋應力值變化明顯放緩甚至出現下降的趨勢，最大值僅為49.5 MPa［14］。類似的情況也發(fā)生在周寧水電站水工隧洞充排水試驗中［15］。然而目前的襯砌設計理論、解析及數值方法尚不能很好地解釋這一現象［11］。

鑒于此，本文基于混凝土塑性損傷模型和透水襯砌理論，考慮在高內壓作用下襯砌損傷開裂、滲流與應力的耦合作用以及襯砌與圍巖有條件聯(lián)合承載特性，基于大型通用有限元軟件平臺ABAQUS進行二次開發(fā)，通過實用程序GETVRM獲取材料損傷，利用子程序USDFLD實現材料滲透系數隨損傷的動態(tài)更新，提出并實現了透水襯砌滲流-應力-損傷耦合算法，模擬分析了水工隧洞初次充水加壓過程，以此研究襯砌損傷開裂過程中的多場耦合作用，探討鋼筋應力的演化特征，為水工隧洞襯砌設計問題提供一定參考。

2　理論基礎

水工隧洞水壓力由襯砌和圍巖承擔。本文采用彈塑性本構關系模擬圍巖，以Mohr-Coulomb準則作為巖體的破壞準則，由于研究重點在于透水襯砌的耦合分析，限于篇幅，在此不做贅述。對于襯砌，本文采用間接耦合的方法對隧洞初次充水期襯砌進行滲流-應力-損傷耦合分析，考慮滲流場通過等效節(jié)點荷載影響應力場分布，在結構自身所承擔的荷載以及等效節(jié)點荷載的作用下計算應力場后即可獲取結構應力、應變以及材料損傷，并建立損傷與滲流之間的耦合關系，根據損傷實時調整材料滲透系數并返回滲流場進行計算，直至計算收斂、加載結束，詳見圖1。

圖1　滲流-應力-損傷耦合關系

2.1滲流場與等效節(jié)點荷載的數值解法假設水和材料不可壓縮，水在材料孔隙間的滲流滿足微分方程［16-17］：

式中：θ為含水率；分別是x、y、z方向上的滲透系數；t為時間變量；水力勢H=p/γ+z，其中p為孔隙水壓力，z為位置水頭，γ為水容重。

計算空間域經過離散、插值和積分，應用Galerkin加權余量法及格林公式可得求解滲流場的矩陣方程：

式中：［Ks］為總滲透矩陣；{h}為節(jié)點水頭列向量；{A}為對滲流邊界積分得到的節(jié)點荷載。

滲流場計算完成后，水力梯度所產生的水荷載以滲透體積力的方式作用在結構內部，按式（3）換算成作用在模型節(jié)點上的等效荷載［3］：

式中：［N］為插值函數；Ω為節(jié)點積分域。

2.2應力場的計算將計算空間域離散化，得到滲流與應力的耦合方程：

式中：［Km］為結構剛度矩陣；{u}為節(jié)點位移列向量；{F}為結構的節(jié)點荷載。

2.3混凝土塑性損傷模型及鋼筋模型相較于預裂縫模型與分布裂縫模型，塑性損傷模型能夠較好地模擬混凝土材料軟化過程、卸載剛度退化及反復加載過程［18-19］，并且不必預先假定混凝土開裂位置，程序會根據結構的受力狀況判定并計算結構各個部位的損傷開裂情況，亦不會出現大面積開裂的情形，在工程分析中應用廣泛［20］。多軸應力狀態(tài)下，混凝土應力σ計算如下：

式中：d為單元損傷；為初始彈性剛度；ε為單元應變；εpl為塑性應變。

單元損傷d由等效塑性應變表示如下：

式中：t、c代表拉、壓狀態(tài)。

在ABAQUS中，混凝土損傷后，應力、損傷與塑性應變的關系可以按照列表形式輸入，本文襯砌材料采用C25混凝土，抗拉、抗壓強度取標準值1.78 MPa、16.7 MPa，結合《混凝土結構設計規(guī)范》（GB50010-2010）附錄C確定相關參數，具體應力、損傷與塑性應變之間的關系如圖2所示。

采用埋藏式模型模擬鋼筋，在迭代計算中鋼筋由線性桿單元模擬，選取理想彈塑性本構關系，應力σs與應變εs關系如下：

圖2　塑性損傷演化曲線

式中：Es為鋼筋彈性模量；fy為屈服應力。

2.4損傷開裂與滲流的耦合關系采用塑性損傷模型時，損傷將導致混凝土結構剛度矩陣的退化，以此模擬混凝土損傷開裂對結構受力狀態(tài)的影響，同時混凝土損傷開裂對結構的滲透特性亦會帶來顯著的影響。材料單元由損傷相與未損傷相兩部分組成，單元滲透系數計算如下［21］：

式中：km為未損傷相滲透系數；kd為損傷相滲透系數；損傷相的塑性體積應變?yōu)閱卧苄泽w積應變；若單元同時出現拉損傷與壓損傷，分別計算對應狀態(tài)下的單元滲透系數，k取較大值。

考慮到脆性材料出現宏觀裂縫后，其滲透特性較開裂前會出現突增現象，因此引入突跳因子ξ以計算損傷相滲透系數［22］：

式中：對于壓剪損傷ξ=100，對于拉損傷

3　技術路線

在耦合計算時，首先根據式（1）—（3）計算滲流場與等效節(jié)點荷載，然后根據式（4）計算應力場，并依據塑性損傷演化曲線獲取材料損傷，根據式（11）和式（12）調整材料滲透系數后重新計算滲流場，如此循環(huán)直至計算收斂、加載結束。

為提高計算效率并保證計算精度［3］，本文內壓分3個階段進行施加以模擬隧洞初次充水運行工況：（1）通過試算確定混凝土即將出現損傷時所對應的內壓，該內壓采取一次性加載的方式進行施加；（2）襯砌損傷后，內壓水頭增量以較小值逐級緩慢施加；（3）貫穿縫出現后，內水外滲通道形成，適當增大內壓水頭增量。每施加一次內壓水頭增量，須保證在該級內壓作用下耦合迭代完全，滲流場達到穩(wěn)定，以保證耦合精度。

本文基于ABAQUS平臺，通過Fortran語言調用子程序USDFLD進行二次開發(fā)，使用USDFLD時，必須調用實用程序GETVRM以獲取材料點的數據，以此完成耦合迭代過程，具體步驟如下（詳見圖3）：

①施加上覆壓力以模擬初始應力場的影響，計算隧洞開挖后應力場，考慮隧洞開挖后，應力釋放完全后進行襯砌支護，即襯砌不承擔開挖釋放荷載，計算襯砌支護后應力場；

②襯砌內壁施加增量水頭計算滲流場，將各單元滲透體積力轉換為等效節(jié)點荷載；

③將等效節(jié)點荷載施加于模型節(jié)點計算得到應力、應變以及材料損傷，調用實用程序GETVRM獲取單元損傷值，對損傷單元根據式（11）和式（12）調用子程序USDFLD調整其滲透系數；

④根據新的滲透系數計算滲流場，檢查滲流場是否穩(wěn)定。當前后兩次滲流場計算結果滿足式（13）時，滲流場達到穩(wěn)定狀態(tài)；如不滿足，則繼續(xù)求解等效節(jié)點荷載施加于模型節(jié)點進行耦合迭代，直至滲流場穩(wěn)定；

式中：Pn+1、Pn為迭代n+1和n次后節(jié)點孔壓；δ為迭代控制閥值，為提高計算效率并保證計算精度，將δ取為 1%［3］。

⑤施加下一級增量水頭，在同一級內壓水頭下重復②—④步，直至加載完成、計算結束。

圖3　耦合迭代流程圖

4　耦合算例分析

4.1模型參數及邊界條件如圖4為一豎井段圓形隧洞軸對稱模型，X、Y軸平行于隧洞斷面，位于水平面內，豎直方向為Z軸，向上為正，襯砌內半徑r1為4.2 m，外半徑r2為5 m，圍巖沿徑向取30倍開挖洞徑［23］，為300 m，模型厚10 m，襯砌與圍巖均由八節(jié)點實體單元模擬，共計9865個節(jié)點和7888個單元，襯砌與圍巖之間共節(jié)點。在襯砌內外兩側鋪設直徑為28 mm的鋼筋，鋼筋到隧洞中心距離R1、R2分別為4.3 m、4.9 m。模型材料參數如表2所示。

圖4　有限元模型

隧洞埋深100 m，沿Z軸負向施加100 m的巖體上覆壓力于模型上表面以模擬初始應力場的影響。在初始時刻，襯砌及圍巖處于飽和狀態(tài)，結構內部為飽和滲流，襯砌及圍巖孔隙率設為1.5%、10%。采用增量水頭法計算滲流場及應力場［24］，忽略重力對增量水頭影響范圍即水力半徑的影響，針對增量水頭的影響范圍將圍巖外圈設為零水頭邊界［25］。模型頂面以及底面設定為不透水邊界，襯砌內壁設為定水頭邊界，計算過程中增大襯砌內壁的定水頭邊界得到增量水頭以模擬充水加壓過程。計算滲流場時，模型施加全位移約束；計算應力場時，模型底部和外圈節(jié)點施加法向位移約束。不考慮模型自身以及洞內水體的重力，只考慮襯砌及圍巖中水體的重力。內壓水頭從0.5 MPa逐級加載至1.1 MPa，每級增量水頭控制在0.025～0.1 MPa，保證每一級內壓水頭下耦合迭代完全。

表2　材料參數

4.2襯砌損傷開裂特征通過試算，內壓小于0.5 MPa時，襯砌不會出現損傷。內壓增至0.65 MPa時，損傷分布均勻，從內向外損傷值由大變小，最大值為0.108，詳見圖5（a）。在0.67 MPa下迭代4次達到穩(wěn)定，襯砌環(huán)向出現3條貫穿縫，等間距分布于0°（隧洞右腰）、120°和240°位置（逆時針），裂縫分布少而寬，詳見圖5（b）。此后直至加載完畢，不再出現新的裂縫，裂縫處損傷值隨內壓持續(xù)增大至0.999，且裂縫間混凝土中部出現拉損傷，原因在于裂縫出現后，襯砌結構類似于彈性地基（曲）梁，使得裂縫間襯砌外側出現反彎點，拉應力增大，出現明顯的拉損傷區(qū)［26］，詳見圖5（c）。通過耦合迭代計算得到的襯砌損傷開裂特征與隧洞壓水試驗結果（表1）吻合較好［7］，說明本文建立的耦合模型能夠較好地反映襯砌的開裂過程。

圖5　襯砌損傷開裂情況

4.3滲透孔壓演化特征為避開模型頂面以及底面邊界條件的影響，選取模型中部位置對滲透孔壓演化特征進行分析。內壓為0.5 MPa時，襯砌未進入塑性狀態(tài)，無損傷，襯砌承擔了大部分的水壓力，圍巖內孔壓相對較小，不超過0.115 MPa，詳見圖6（a）。在內壓達到0.67 MPa，襯砌開裂前，襯砌內側出現損傷，損傷區(qū)域孔壓較大，圍巖內孔壓有所增大，最大值為0.257 MPa，詳見圖6（b）；3條貫穿縫出現后，內水外滲，滲流場出現明顯的重分布現象，圍巖內孔壓升高，在襯砌環(huán)向0°、120°及240°裂縫位置滲流場變化明顯，孔壓梯度明顯降低，裂縫處圍巖承受的水壓接近內壓水頭，孔壓自裂縫位置向圍巖深部逐漸消散，滲流場分布與襯砌損傷演化情況相對應，基本滿足等間距分布，詳見圖 6（c）。

圖6　隧洞滲流場分布（單位：Pa）

隨著內壓的增大，襯砌損傷峰值與滲透系數峰值同步變化，表示材料損傷與滲透系數實現了同步更新。裂縫處混凝土損傷與滲透系數在襯砌開裂時出現陡增現象，此后隨內壓持續(xù)增大，加載完成時，滲透系數相較于初始1×10-9m/s提高了約3個量級，裂縫成為內水外滲的主要通道，詳見圖7。

4.4鋼筋應力演化特征充水加壓過程中，裂縫處鋼筋應力隨內壓增大，而未開裂位置鋼筋全程處于低應力狀態(tài)，應力值小于14 MPa。當內壓小于0.67 MPa時，鋼筋應力緩慢增長；貫穿縫形成后，裂縫位置混凝土應力釋放，大部分荷載轉由鋼筋承擔，裂縫處鋼筋應力突增至50.120 MPa；當內壓從0.67 MPa增至0.9 MPa時，裂縫處鋼筋應力緩慢增至62.770 MPa，當內壓達到1.0 MPa時，裂縫處鋼筋應力進一步增大，在內壓1.1 MPa迭代穩(wěn)定后裂縫處鋼筋應力達到228.800 MPa，詳見圖8。原因在于貫穿縫出現后，裂縫處的鋼筋所發(fā)揮的作用類似于錨栓，單獨承擔拉力，而未開裂位置鋼筋與襯砌混凝土共同承擔襯砌環(huán)向拉力，混凝土承擔大部分拉力，鋼筋所承擔的拉力較小。

圖7　襯砌損傷峰值及滲透系數峰值隨內壓的演化

圖8　鋼筋應力演化特征

5　襯砌與圍巖有條件聯(lián)合承載機理研究

襯砌開裂前，透水性遠弱于圍巖，內部水力梯度較大，承擔大部分的水荷載，而圍巖內部水力梯度較小，承擔的水荷載較小，因此襯砌沿徑向朝外變形的趨勢更大，襯砌與圍巖處于壓緊狀態(tài)；襯砌開裂后，透水性顯著增強，內部水力梯度大幅降低，所承擔的水荷載明顯減小，此時圍巖的透水性相對較弱，內部水力梯度相對較大，承擔大部分的水荷載，因此圍巖沿徑向朝外變形的趨勢更為顯著，圍巖拉著襯砌向外變形，同時襯砌受到內置鋼筋的約束，在襯砌與圍巖交界面上產生拉應力，當拉應力超過交界面粘結強度時，襯砌與圍巖將出現脫離現象，表現出有條件聯(lián)合承載特性［10，27-29］。

在原有模型的基礎上，對計算模型進行局部修改，考慮襯砌與圍巖之間不共節(jié)點，在襯砌與圍巖之間建立面-面接觸單元以模擬襯砌與圍巖有條件聯(lián)合承載特性，接觸單元切向摩擦系數取為0.5，局部襯砌-接觸單元-圍巖模型如圖9所示。此時襯砌與圍巖接觸面只存在法向力與切向摩擦力，當法向力為壓力時，兩者處于壓緊狀態(tài)，一旦法向受拉，兩者脫離，不再傳遞法向應力，也不再滿足位移連續(xù)條件，即默認襯砌與圍巖交界面粘結強度為零，接觸壓力與間隙的關系如圖10所示。在應力場計算中，接觸面的作用主要是傳遞法向應力；在滲流場計算中，接觸可以保證孔隙流體在不同材料之間的滲透作用。

5.1襯砌與圍巖的接觸狀態(tài)內壓增加到0.625 MPa前，襯砌未開裂，襯砌滲透特性遠弱于圍巖，襯砌承擔了大部分水荷載，在0.625 MPa內壓下迭代前，10 m厚度范圍內襯砌內壁內水外滲總流量僅為1.109 L/min，此時襯砌與圍巖不存在脫離現象，詳見圖11（a）；在0.625 MPa內壓下迭代完成后，襯砌環(huán)向在290°位置（襯砌右腰為0°位置，逆時針旋轉）形成一條徑向貫穿裂縫，同時襯砌除貫穿縫外的其他部位的損傷值在迭代過程中有所增大，滲透性增強，兩部分因素的作用下使得襯砌內壁內水外滲總流量陡增至6.850 L/min。此后隨著內壓增大，不再出現新的裂縫，襯砌損傷開裂特征與隧洞壓水試驗一致（表1）。由于裂縫處水力梯度的降低以及襯砌向外變形受到內置鋼筋的約束，襯砌與圍巖首先在290°開裂位置出現脫離現象，產生的空隙由滲透水充滿，脫離范圍在迭代過程中逐步擴大，迭代11次后達到穩(wěn)定狀態(tài)，脫離范圍達75%，詳見圖11（b）。隨著內壓增至0.635 MPa，迭代4次后達到穩(wěn)定狀態(tài)，襯砌與圍巖已經完全脫離，詳見圖11（c）。

圖9　襯砌-接觸單元-圍巖模型

圖10　接觸壓力與間隙的關系

圖11　襯砌與圍巖接觸狀態(tài)（深色為脫離部分）

5.2鋼筋應力演化特征引入接觸單元后，襯砌鋼筋應力表現出明顯的差異。襯砌開裂前鋼筋應力較小，內壓達到0.625 MPa時，襯砌開裂，鋼筋應力呈整體上升趨勢，且存在較大波動，迭代過程中裂縫處鋼筋應力最大值達到327.500 MPa，此后隨迭代過程鋼筋應力逐漸趨于穩(wěn)定，并在水壓力緩慢增加的過程中，鋼筋應力基本穩(wěn)定并呈緩慢降低的趨勢。在1.1 MPa內壓作用下，鋼筋應力最大值為56.560 MPa，整體量值較小，詳見圖12。出現此現象的原因在于襯砌開裂后，襯砌與圍巖逐步脫離，圍巖成為水荷載的主要承擔對象，而襯砌承擔的水荷載逐漸減小，導致鋼筋應力量值較小。同時與完全聯(lián)合承載情況不同，完全聯(lián)合承載時襯砌在開裂之后仍然會隨著圍巖向外擴張變形，導致鋼筋應力持續(xù)增大。相比之下襯砌與圍巖相互脫離，即有條件聯(lián)合承載極有可能是隧洞發(fā)電運行過程中鋼筋應力普遍小于設計預期的原因，并且鋼筋應力的這一演化特征和量值很好地解釋了相關工程的監(jiān)測結果［14-15，30］。

圖12　有條件聯(lián)合承載時鋼筋應力演化特征

6　結語

本文基于混凝土塑性損傷模型與透水襯砌理論，通過等效節(jié)點荷載實現了滲流場對應力場的作用，建立了損傷開裂與滲流之間的耦合關系，以此提出了透水襯砌滲流-應力-損傷耦合算法。基于ABAQUS平臺進行二次開發(fā)，通過Fortran語言調用實用程序GETVRM獲取材料點的損傷值，利用子程序USDFLD實現了材料滲透系數隨損傷的動態(tài)更新，完成了隧洞初次充水期間襯砌損傷開裂的耦合迭代過程。研究了襯砌損傷開裂、滲透孔壓以及鋼筋應力的演化特征，并采用面-面接觸單元探討了襯砌與圍巖有條件聯(lián)合承載機理，可以發(fā)現：（1）本文提出的耦合算法所得到的襯砌損傷開裂特征與襯砌壓水試驗結果吻合，表現為襯砌在一定內壓下開裂，此后不再出現新的裂縫，只在裂縫間襯砌外側出現拉損傷，較好地反映了水工隧洞初次充水加壓工況中襯砌損傷開裂的過程。（2）襯砌開裂前后，滲流場表現出明顯的差異。開裂后，裂縫處孔壓梯度明顯降低，滲透系數呈現數量級的增大，裂縫成為內水外滲的主要通道，這與工程普遍認知是相符的，說明本文所提出的耦合算法是可行的，具有一定的推廣性。（3）完全聯(lián)合承載時，襯砌與圍巖之間徑向變形協(xié)調，裂縫處鋼筋應力持續(xù)增大，未開裂區(qū)鋼筋全程維持在低應力水平；考慮有條件聯(lián)合承載特性后，襯砌開裂之后與圍巖逐步脫離，圍巖成為承擔水荷載的主體，襯砌鋼筋應力出現應力回縮現象，在加載完畢時應力量值維持在較小水平，揭示了有條件聯(lián)合承載是隧洞發(fā)電運行過程中鋼筋應力普遍小于設計預期的原因，解釋了相關工程的監(jiān)測結果，說明本文所提出的耦合算法可以為隧洞襯砌設計問題提供一定參考。

參考文獻：

［1］楊林德，丁文其.滲水高壓引水隧洞襯砌的設計研究［J］.巖石力學與工程學報，1997，16（2）：17-22.

［2］杜小凱，任青文，陳偉 .有壓引水隧洞內水外滲作用研究［J］.四川大學學報（工程科學版），2008，40（5）：63-68.

［3］卞康，肖明.水工隧洞襯砌水壓致裂過程的滲流-損傷-應力耦合分析［J］.巖石力學與工程學報，2010，29（S2）：3769-3776.

［4］SIMANJUNTAK T D Y F，MARENCE M，MYNETT A E，et al.Mechanical-hydraulic interaction in the lining cracking process of pressure tunnels［J］.Int.J.Hydropow.Dams，2013，20（5）：98-105.

［5］BIAN K，XIAO M，CHEN J.Study on coupled seepage and stress fields in the concrete lining of the underground pipe with high water pressure［J］.Tunnelling and Underground Space Technology，2009，24：287-295.

［6］FAHIMIFAR A，MOHAMMAD R Z.A theoretical solution for analysis of tunnels below groundwater considering the hydraulic-mechanical coupling［J］.Tunnelling and Underground Space Technology，2009，24：634-646.

［7］張有天.水工隧洞建設的經驗和教訓（上）［J］.貴州水力發(fā)電，2001，15（4）：76-84.

［8］劉秀珍.關于圓形有壓水工隧洞襯砌裂縫計算的探討［J］.水力發(fā)電，1983，12：39-46.

［9］沈威.圓形水工壓力隧洞襯砌變形特性與限裂設計研究［D］.大連：大連理工大學，2012.

［10］蘇凱，伍鶴皋，周創(chuàng)兵 .內水壓力下水工隧洞襯砌與圍巖承載特性研究［J］.巖土力學，2010，31（8）：2407-2412.

［11］ZHOU Y，SU K，WU H.Hydro-mechanical interaction analysis of high pressure hydraulic tunnel［J］.Tunnelling and Underground Space Technology，2015，47：28-34.

［12］GRAZIANI A，BOLDINI D.Influence of hydro-mechanical coupling on tunnel response in clays［J］.Journal of Geotechnical&Geoenvironmental Engineering，2012，138（3）：415-418.

［13］YOO C.Interaction between tunneling and groundwater-numerical investigation using three dimensional stress-pore pressure coupled analysis［J］.Journal of Geotechnical&Geoenvironmental Engineering，2005，131（2）：240-250.

［14］葉冀升.廣蓄電站水工高壓隧洞設計施工的若干問題［J］.水力發(fā)電學報，1998（2）：39-50.

［15］李煊明.周寧水電站透水襯砌高壓鋼筋混凝土岔管研究［J］.水電站設計，2007，23（2）：15-17.

［16］吳夢喜.飽和—非飽和土中滲流Richards方程有限元算法［J］.水利學報，2009，40（10）：1274-1279.

［17］紀佑軍，劉建軍，程林松 .考慮流-固耦合的隧道開挖數值模擬［J］.巖土力學，2011，32（4）：1229-1233.

［18］LUBLINER J，OLIVER J，OLLER S，et al.A plastic-damage model for concrete［J］.International Journal of Solids&Structures，1989，25（3）：299-326.

［19］LEE J，FENVES G L.Plastic-damage model for cyclic loading of concrete structures［J］.Journal of Engineering Mechanics，1998，124（8）：892-900.

［20］蘇凱，張偉，伍鶴皋，等.考慮摩擦接觸特性的鋼襯鋼筋混凝土管道承載機理研究［J］.水利學報，2016，47（8）：1070-1078.

［21］賈善坡，陳衛(wèi)忠，于洪丹，等.泥巖隧道施工過程中滲流場與應力場全耦合損傷模型研究［J］.巖土力學，2009，30（1）：19-26.

［22］周亞峰.水工隧洞滲流-應力-損傷/開裂耦合理論及方法研究［D］.武漢：武漢大學，2016.

［23］周亞峰，蘇凱，伍鶴皋.水工隧洞鋼筋混凝土襯砌外水壓力取值方法研究［J］.巖土力學，2014，35（S2）：198-203.

［24］張有天，張武功.隧洞水荷載的靜力計算［J］.水利學報，1980（3）：52-62.

［25］潘家錚，唐家驥，張有天，等.隧洞水荷載靜力計算（學術討論）［J］.水利學報，1981（5）：75-81.

［26］李寧，劉波 .城門洞型壓力隧洞襯砌限裂配筋設計方法研究［J］.巖石力學與工程學報，2013（12）：2409-2414.

［27］張有天，張武功，王鐳.再論隧洞水荷載的靜力計算［J］.水利學報，1985（3）：21-31.

［28］BOBET A，NAM S W.Stresses around pressure tunnels with semi-permeable liners［J］.Rock Mechanics and Rock Engineering，2007，40（3）：287-315.

［29］FERNANDEZ G.Behavior of pressure tunnels and guidelines for liner design［J］.Journal of Geotechnical Engineering，1994，120（10）：1768-1791.

［30］劉素琴，張秀麗，劉啟釗，等.高水頭鋼筋混凝土岔管結構計算與觀測對比分析［J］.水利水電技術，2006，37（8）：24-28.