考慮滲流—應(yīng)力耦合的黏土心墻壩邊坡穩(wěn)定分析

黃 香，崔志瑞，郭瑩瑩，李 穎

(1.三峽金沙江云川水電開發(fā)有限公司，云南 昆明 650224；2.華北水利水電大學(xué)， 河南 鄭州 450046)

1 概述

對(duì)于黏土心墻壩而言，滲流與穩(wěn)定問題關(guān)系到大壩能否安全運(yùn)行[1]。而邊坡穩(wěn)定分析是黏土心墻壩設(shè)計(jì)的基礎(chǔ)，通過邊坡穩(wěn)定分析可反映出大壩設(shè)計(jì)方案的合理性，為大壩結(jié)構(gòu)設(shè)計(jì)提供科學(xué)依據(jù)[2]。目前關(guān)于黏土心墻壩邊坡穩(wěn)定的研究成果較多，如文獻(xiàn)[3-4]等均采用傳統(tǒng)的邊坡穩(wěn)定計(jì)算方法，對(duì)一些黏土心墻壩的邊坡穩(wěn)定情況進(jìn)行了計(jì)算分析。許多研究成果[5-6]表明，在天然狀態(tài)下，滲流與土體孔隙空間相互關(guān)聯(lián)，滲流場(chǎng)與應(yīng)力場(chǎng)是一種耦合關(guān)系。對(duì)土石壩(包括黏土心墻壩)進(jìn)行邊坡穩(wěn)定分析時(shí)，只考慮單獨(dú)的滲流場(chǎng)或應(yīng)力場(chǎng)進(jìn)行計(jì)算，得到的結(jié)果會(huì)與實(shí)際情況存在差異，采用傳統(tǒng)的邊坡穩(wěn)定計(jì)算方法，會(huì)存在一些不足。另一些文獻(xiàn)[7-8]考慮了滲流—應(yīng)力雙場(chǎng)耦合的作用，并對(duì)黏土心墻壩的應(yīng)力變形分布規(guī)律和工作性態(tài)進(jìn)行了計(jì)算分析，但對(duì)雙場(chǎng)耦合作用對(duì)邊坡穩(wěn)定的影響并沒有進(jìn)行更深入的研究。綜合分析文獻(xiàn)可知，目前考慮滲流—應(yīng)力雙場(chǎng)耦合的土石壩邊坡穩(wěn)定分析，較多研究偏重于尾礦壩或均質(zhì)土壩[9-10]，而考慮雙場(chǎng)耦合的黏土心墻壩邊坡穩(wěn)定的研究則比較少。

另外，從研究方法層面，目前大多數(shù)研究采用有限元強(qiáng)度折減法進(jìn)行土石壩邊坡穩(wěn)定計(jì)算分析。通過此法進(jìn)行邊坡穩(wěn)定計(jì)算時(shí)，邊坡失穩(wěn)判據(jù)大致有三種：判據(jù)一，數(shù)據(jù)計(jì)算是否收斂；判據(jù)二，塑性區(qū)是否貫通；判據(jù)三，特征點(diǎn)位移突變[11]。根據(jù)判據(jù)進(jìn)行土石壩邊坡穩(wěn)定計(jì)算分析時(shí)，大多數(shù)研究采用Drucker-Prager模型等彈塑性模型進(jìn)行壩體的位移和應(yīng)力計(jì)算[9-11]。盡管上述理想彈塑性模型能很好描述巖土體的非線性特性，但此類模型如屈服面存在尖角導(dǎo)致計(jì)算繁瑣，且收斂速度緩慢，甚至不收斂的缺陷[9]。因此，當(dāng)采用此類模型進(jìn)行邊坡穩(wěn)定計(jì)算分析，并結(jié)合判據(jù)一和判據(jù)二進(jìn)行邊坡失穩(wěn)破壞判斷時(shí)，容易造成誤判。而鄧肯—張模型做為一種增量型彈性模型，可克服上述彈塑性模型計(jì)算時(shí)不易收斂的缺陷，而且可以在一定程度上反映土體的彈塑性變形，目前鄧肯—張模型在巖土體的應(yīng)力變形計(jì)算分析中得到了廣泛的應(yīng)用，但考慮滲流—應(yīng)力雙場(chǎng)相互耦合，相互影響，并基于鄧肯—張模型的邊坡穩(wěn)定有限元強(qiáng)度折減法的研究則相對(duì)尚較少。

綜上所述，本文利用ANSYS軟件的APDL命令流語言并結(jié)合Fortran語言，編制模塊化的滲流—應(yīng)力耦合計(jì)算程序，并基于鄧肯—張模型的邊坡穩(wěn)定有限元強(qiáng)度折減法，對(duì)某黏土心墻壩的邊坡穩(wěn)定進(jìn)行計(jì)算分析，以期為黏土心墻壩的邊坡穩(wěn)定分析提供借鑒和參考。

2 計(jì)算原理

2.1 滲流—應(yīng)力耦合計(jì)算原理

目前，滲流—應(yīng)力雙場(chǎng)耦合計(jì)算方法主要為直接耦合法和間接耦合法，考慮直接耦合計(jì)算時(shí)需建立復(fù)雜的數(shù)學(xué)模型且求解困難，因此，本文采用間接耦合方法進(jìn)行雙場(chǎng)的耦合計(jì)算。

1)滲流場(chǎng)對(duì)應(yīng)力場(chǎng)的影響

滲流場(chǎng)對(duì)應(yīng)力場(chǎng)的影響有動(dòng)水壓力(滲透體積力)和靜水壓力(滲透壓力)兩種形式。靜水壓力(滲透壓力)為作用與接觸面上的力，計(jì)算時(shí)可按面力作用形式計(jì)入節(jié)點(diǎn)荷載，即：

(1)

式中NT(x)為形函數(shù)；Γ為滲流邊界；p為滲透壓力，其可表示為：

p=γw(H-y)

(2)

式中γw為水的容重；H為水頭；y為位置y向坐標(biāo)。

動(dòng)水壓力(滲透體積力)可視為與自重類似，計(jì)算時(shí)按體力計(jì)入節(jié)點(diǎn)荷載，即：

fs=?ΩNT(x)fpdΩ

(3)

式中NT(x)為形函數(shù)；Ω為滲流計(jì)算區(qū)域；fp為動(dòng)水壓力，其可表示為：

(4)

考慮滲流場(chǎng)影響的應(yīng)力場(chǎng)控制方程，則可表示為：

(5)

2)應(yīng)力場(chǎng)對(duì)滲流場(chǎng)的影響

在荷載作用下，應(yīng)力場(chǎng)的變化，會(huì)導(dǎo)致巖土體的體積及孔隙發(fā)生變化，從而導(dǎo)致其滲透系數(shù)發(fā)生改變。應(yīng)力與滲透系數(shù)的關(guān)系可用文獻(xiàn)[6]推薦的關(guān)系式表達(dá)，即：

K=K0exp[-β(σ-p)]

(6)

式中σ為土體應(yīng)力；K0為初始滲透系數(shù)；β為經(jīng)驗(yàn)系數(shù)，取為3e-7；p為靜水壓力。

考慮應(yīng)力場(chǎng)影響的滲流場(chǎng)控制方程，可表示為：

(7)

式中H為水頭分布；K(σij)為滲透系數(shù)，可由式(6)表示；Ω為滲流計(jì)算區(qū)域；Γ1為已知流量邊界；Γ2為自由面邊界；Γ3為已知水頭邊界；n2和n3為對(duì)應(yīng)邊界的法線方向。

2.2 基于鄧肯—張模型的強(qiáng)度折減法原理

利用有限元強(qiáng)度折減法時(shí)，黏聚力c、內(nèi)摩擦角φ按下式進(jìn)行計(jì)算：

cf=c/Fs,φf=arctan(tanφ/Fs)

(8)

式中Fs為折減系數(shù)(即安全系數(shù))；cf、φf分別為折減后的材料內(nèi)聚力和內(nèi)摩擦角。

對(duì)于鄧肯—張模型，其切線模量Et可由下式表示：

(9)

式中Ei為初始彈性模量；Rf為破壞比；兩者的具體表達(dá)式可參見文獻(xiàn)[5]。而上式中(σ1-σ3)f表示破壞時(shí)的主應(yīng)力差，可由Mohr-Coulomb定律推導(dǎo)得到，具體為：

(10)

將式(10)及Ei的表達(dá)式代入式(9)，可得切線模量的表達(dá)式為：

(11)

從式(11)可以看出，只要對(duì)式中的黏聚力c和內(nèi)摩擦角φ進(jìn)行折減，即可實(shí)現(xiàn)與Drucker-Prager模型或Mohr-Coulomb模型相似的材料強(qiáng)度折減效果，另外本文在進(jìn)行強(qiáng)度折減時(shí)，還考慮將初始切線模量k進(jìn)行折減，以充分滿足土體的變形性質(zhì)，即：

kf=k/Fs

(12)

將式(8)和式(12)代入式(11)，即可得到考慮有限元強(qiáng)度折減法的鄧肯—張模型切線模量計(jì)算表達(dá)式，如式(13)所示：

(13)

3 工程實(shí)例及結(jié)果分析

3.1 工程實(shí)例

某黏土心墻壩，工程總庫容為110.87萬m3，工程規(guī)模為小(1)型，工程等別為Ⅳ等，永久性水工建筑物級(jí)別為4級(jí)，臨時(shí)性水工建筑物級(jí)別為5級(jí)?？绾硬贾每傞L(zhǎng)為185.00 m，壩頂高程為1 804.30 m，壩頂設(shè)高1.2m混凝土防浪墻，墻頂高程為1 805.50 m，壩頂寬度為5.00 m。大壩位于河谷處，最大壩高斷面及材料分區(qū)如圖1所示。限于篇幅，本文中只給出正常蓄水工況的計(jì)算結(jié)果，大壩正常蓄水位為1 800.08 m。根據(jù)實(shí)驗(yàn)資料[12]，本次計(jì)算分析時(shí)，壩體材料采用鄧肯—張E-v模型，具體材料參數(shù)見表1，壩基及帷幕各材料假定為線彈性材料(見表2)。

圖1 黏土心墻壩計(jì)算斷面及材料分區(qū)示意(單位：m)

表1 壩體鄧肯—張E-v模型材料參數(shù)

表2 壩基彈性模型材料參數(shù)

3.2 計(jì)算結(jié)果分析

進(jìn)行滲流—應(yīng)力耦合計(jì)算時(shí)，計(jì)算斷面的有限元網(wǎng)格如圖2所示。有限元模型的頂邊界截至壩頂高程為1 840.30 m，底邊界截取至壩基帷幕最深處以下，至高程為1 680.00 m處。本次有限元模擬計(jì)算，滲流計(jì)算邊界條件取為：壩基底部為不透水邊界，左側(cè)為上游水頭邊界，右側(cè)無水，因此不加水頭，但下游側(cè)斜邊視為可能溢出邊界；應(yīng)力計(jì)算部分邊界條件為：壩基底部為x、y向位移約束，壩基上下游側(cè)x向位移約束。模型的坐標(biāo)系選用笛卡爾直角坐標(biāo)系，x軸指向下游為正，y軸向上為正。

圖2 計(jì)算斷面有限元網(wǎng)格示意

1)滲流計(jì)算結(jié)果分析

從圖3可以看出:考慮和不考慮耦合作用時(shí)，壩體的水頭等值線分布規(guī)律大致相同。但從中也可看出考慮耦合作用時(shí)，水頭為90.737 m和92.884 m的等值線更靠近下游，這是由于考慮耦合作用后，受應(yīng)力場(chǎng)的影響使材料的滲透系數(shù)有所增大，心墻削減水頭的能力減弱。

a 考慮耦合作用

b 不考慮耦合作用

2)應(yīng)力變形計(jì)算結(jié)果分析

從圖4可以看出，考慮耦合后壩體的y向位移與不考慮耦合時(shí)相比變大，這是由于耦合作用滲透水壓力的豎向分力使壩體的土體有下沉趨勢(shì)，導(dǎo)致考慮耦合效應(yīng)后，壩體的y向位移增大。

從圖5可以看出，考慮耦合和不考慮耦合情況下，壩體的y向應(yīng)力均符合一般規(guī)律。但不考慮耦合情況下，心墻y向應(yīng)力小于兩側(cè)壩殼料的應(yīng)力，具有拱效應(yīng)，這是因?yàn)樾膲Φ淖冃文Ａ枯^兩側(cè)壩殼料小，心墻位移大，兩側(cè)壩殼料位移小，心墻部分自重傳遞到壩殼，反映了心墻內(nèi)部應(yīng)力小于兩側(cè)壩殼料的應(yīng)力?？紤]耦合情況下，心墻部位也產(chǎn)生了明顯的拱效應(yīng)，但考慮耦合情況下，心墻上下游兩側(cè)應(yīng)力分布不成對(duì)稱分布，且心墻上下游兩側(cè)同一高程處，上游側(cè)y向應(yīng)力大于下游側(cè)，這是由于耦合情況下，考慮了滲透水壓力的影響，上游水頭較高，滲透水壓力較大，導(dǎo)致上游土體壓力大于下游土體的壓力。

a 考慮耦合作用

b 不考慮耦合作用

a 考慮耦合作用

b 不考慮耦合作用

3)邊坡穩(wěn)定計(jì)算結(jié)果分析

考慮邊坡失穩(wěn)判據(jù)1和判據(jù)2的特點(diǎn)及局限性以及失穩(wěn)判據(jù)與鄧肯—張模型的相容性，本文采用判據(jù)3對(duì)邊坡的失穩(wěn)破壞進(jìn)行判定。采用判據(jù)3時(shí)，需確定特征點(diǎn)后將水平或豎向位移作為位移突變的判據(jù)。因?yàn)轲ね列膲ψ冃伪葍蓚?cè)壩殼料大，更易導(dǎo)致邊坡的失穩(wěn)破壞，故本文選取心墻頂?shù)? 672號(hào)節(jié)點(diǎn)作為特征點(diǎn)(圖2中A點(diǎn))。在開始計(jì)算時(shí)，折減系數(shù)以0.1為間隔進(jìn)行變化，以找出安全系數(shù)的大致范圍，再以0.025或更小的間隔進(jìn)行折減，直到找出更精確的安全系數(shù)。根據(jù)上述方法得到的特征點(diǎn)豎向位移與安全系數(shù)的關(guān)系曲線如圖6所示。

圖6 特征點(diǎn)豎向位移與安全系數(shù)的關(guān)系曲線示意

從圖6可以看出，當(dāng)考慮耦合時(shí)，折減系數(shù)由1.2增大為1.3時(shí)，特征點(diǎn)的位移有明顯的突變。對(duì)折減系數(shù)進(jìn)行加密計(jì)算，當(dāng)安全系數(shù)在1.250～1.275之間時(shí)，有明顯的的位移突變，因此，考慮耦合計(jì)算時(shí)，本文所研究黏土心墻壩的邊坡穩(wěn)定安全系數(shù)為1.250。同樣可以看出，當(dāng)不考慮耦合計(jì)算時(shí)，當(dāng)折減系數(shù)由1.3增大為1.4時(shí)，特征點(diǎn)的位移有明顯的突變。同樣，對(duì)折減系數(shù)進(jìn)行加密，可以看出安全系數(shù)在1.300～1.400之間時(shí)，有明顯的的位移突變，因此，不考慮耦合計(jì)算時(shí)，可確定本文所研究黏土心墻壩的邊坡穩(wěn)定安全系數(shù)為1.375。從圖6還可以看出，考慮耦合后計(jì)算所得壩體邊坡的穩(wěn)定安全系數(shù)比不考慮滲流—應(yīng)力耦合時(shí)小，究其原因是，一方面考慮耦合作用時(shí)，滲透水壓力豎向分力使壩體的土體有下沉趨勢(shì)，滲流力加大了壩坡的滑移力，從而增大了壩體邊坡土體的下滑趨勢(shì)，另一方面考慮耦合作用時(shí)，受滲流場(chǎng)影響，土體的變形模量變小，土體軟化，更易產(chǎn)生變形，對(duì)壩坡抗滑穩(wěn)定不利。

4 結(jié)語

本文利用ANSYS軟件，建立了某黏土心墻壩的有限元模型，計(jì)算結(jié)果表明：

1)考慮耦合作用后，受應(yīng)力的影響，壩體土料的滲透系數(shù)增大后會(huì)使心墻削減水頭的能力有所減弱。

2)考慮耦合作用后，滲透水壓力豎向分力使壩體的土體有明顯下沉趨勢(shì)，豎向位移較不考慮耦合時(shí)明顯增大，而且滲透水壓力的作用使得壩體上游側(cè)土體的應(yīng)力明顯比下游大。

3)考慮耦合作用進(jìn)行邊坡穩(wěn)定計(jì)算分析時(shí)，滲流場(chǎng)作用不僅會(huì)影響壩坡的滑移力，而且會(huì)影響土體的變形模量，對(duì)壩坡的抗滑穩(wěn)定不利，從而使計(jì)算得到邊坡穩(wěn)定安全系數(shù)明顯比不考慮耦合作用的小。為偏安全，在進(jìn)行黏土心墻壩有限元計(jì)算分析及邊坡穩(wěn)定計(jì)算時(shí)，建議考慮滲流場(chǎng)—應(yīng)力場(chǎng)的雙場(chǎng)耦合效應(yīng)。