庫(kù)水位變化和降雨作用下三門(mén)洞滑坡流固耦合分析

徐興倩 屈新 張新啟 竇思軍 陳肖肖

摘要：目前對(duì)三門(mén)洞滑坡的已有研究都重點(diǎn)關(guān)注穩(wěn)定系數(shù)和變形特征，缺乏對(duì)孔隙水壓力及應(yīng)力演化特征的系統(tǒng)性研究。采用Abaqus建立滑坡的三維計(jì)算模型，推導(dǎo)庫(kù)水位變化和降雨入滲的耦合邊界條件，提出基于莫爾一庫(kù)倫和最大拉應(yīng)力準(zhǔn)則相結(jié)合的復(fù)合準(zhǔn)則，模擬滑坡在庫(kù)水位變化和降雨共同作用下的變形演化過(guò)程。系統(tǒng)探討了滑坡的孔隙水壓力分布、位移場(chǎng)、應(yīng)力場(chǎng)和塑性區(qū)分布特征，并綜合數(shù)值模擬結(jié)果分析了庫(kù)水位變化和降雨入滲條件下的邊坡穩(wěn)定性問(wèn)題。結(jié)果表明：在強(qiáng)降雨和庫(kù)水位緩慢下降綜合作用下，三門(mén)洞滑坡在滑坡中后部和前緣分別出現(xiàn)局部的張拉和剪切破壞，滑體內(nèi)滲流場(chǎng)發(fā)生較大變化，導(dǎo)致滑體內(nèi)孔隙水壓力、應(yīng)力、位移和等效塑性應(yīng)變整體增大，但并未出現(xiàn)塑性破壞區(qū)，邊坡整體處于穩(wěn)定狀態(tài)。

關(guān)鍵詞：滑坡穩(wěn)定性;流固耦合;庫(kù)水位變化;降雨入滲;三門(mén)洞滑坡;三峽庫(kù)區(qū)

中圖法分類(lèi)號(hào)：P642.22 文獻(xiàn)標(biāo)志碼：A DOI：10.15974/j.cnki.slsdkb.2020.04.005

1研究背景

三峽庫(kù)區(qū)的地質(zhì)環(huán)境因素復(fù)雜，庫(kù)水位變化大，降雨頻繁，導(dǎo)致庫(kù)區(qū)邊坡極易出現(xiàn)不穩(wěn)定的情況，給當(dāng)?shù)厝嗣竦纳拓?cái)產(chǎn)安全造成了威脅。已有很多學(xué)者從庫(kù)水位變化或降雨角度對(duì)邊坡穩(wěn)定性開(kāi)展了大量研究工作。Ng和Shi對(duì)暴雨條件下土坡的暫態(tài)滲流場(chǎng)和穩(wěn)定性進(jìn)行了分析，探討了正負(fù)孔隙水壓力演化規(guī)律，并發(fā)現(xiàn)降雨入滲導(dǎo)致坡體內(nèi)部負(fù)孔隙水壓力減小，邊坡穩(wěn)定系數(shù)降低。劉才華等探討了庫(kù)水位上升邊坡失穩(wěn)的誘發(fā)機(jī)理。劉金龍等采用有限元模擬了庫(kù)水位變化對(duì)邊坡穩(wěn)定性的影響。岳慶河等分析了庫(kù)水位回落對(duì)土石壩邊坡穩(wěn)定性的影響。吳海燕等研究了降雨人滲導(dǎo)致的黃土邊坡破壞面形成過(guò)程及滑動(dòng)機(jī)理。齊云龍等討論了地震和降雨共同作用下邊坡穩(wěn)定性，針對(duì)性地提出了邊坡防治建議。婁一青等綜合有限元數(shù)值模擬分析結(jié)果得出了降雨條件下邊坡地下水的滲流特征。尚敏等采用Geo-studio模擬了白水河滑坡在降雨入滲作用下的滲流場(chǎng)，探討了地表變形幅度與降雨強(qiáng)度之間的關(guān)系。上述研究揭示了降雨或者庫(kù)水位變化對(duì)邊坡穩(wěn)定性的影響規(guī)律，但都是單一地從庫(kù)水位變化或降雨入滲的角度分析邊坡穩(wěn)定性。莫偉偉研究了庫(kù)水位變化和降雨條件下滑坡水巖流固耦合的作用機(jī)理，并基于BP神經(jīng)網(wǎng)絡(luò)技術(shù)和遺傳算法獲取了庫(kù)水位變化和降雨條件下巖土體的物理力學(xué)參數(shù)。帥紅巖采用有限元模擬了三峽庫(kù)區(qū)曬鹽壩滑坡在庫(kù)水位和降雨共同作用下的滑坡滲流場(chǎng)變化情況，并基于非飽和土抗剪強(qiáng)度理論和Fredlund的極限平衡法理論，分析了滑坡的穩(wěn)定性。秦洪斌提出了庫(kù)水和降雨共同作用下復(fù)活型滑坡的復(fù)活判據(jù)，并基于3個(gè)具體滑坡驗(yàn)證了該判據(jù)的有效性。上述研究探討了庫(kù)水位變化和降雨共同作用下的邊坡變形穩(wěn)定性問(wèn)題，為本文的研究工作提供了理論支撐。

已有很多學(xué)者對(duì)三峽庫(kù)區(qū)的典型滑坡——三門(mén)洞滑坡進(jìn)行了研究。薛聰聰?shù)染C合考慮庫(kù)水位和降雨滲流的作用，分析了不同工況下三門(mén)洞滑坡的穩(wěn)定系數(shù)，并得出了定性結(jié)論。易慶林等n。基于三門(mén)洞滑坡的地質(zhì)勘探數(shù)據(jù)，綜合分析了該滑坡體的動(dòng)態(tài)變形特征和主要影響因素。上述研究均探討了三門(mén)洞滑坡的變形及穩(wěn)定系數(shù)，但很少對(duì)滑坡體的孔隙水壓力分布、位移場(chǎng)、應(yīng)力場(chǎng)和塑性區(qū)等進(jìn)行系統(tǒng)性研究。

基于此，本文擬以三門(mén)洞滑坡為研究對(duì)象，采用有限元軟件Abaqus模擬其在水庫(kù)水位變化和降雨共同作用下的變形演化過(guò)程，獲取該滑坡的孔隙水壓力分布、位移場(chǎng)、應(yīng)力場(chǎng)和塑性區(qū)，并綜合數(shù)值模擬結(jié)果分析庫(kù)水位變化和降雨對(duì)邊坡穩(wěn)定性的影響規(guī)律。

2三門(mén)洞滑坡地貌及有限元網(wǎng)格模型

三門(mén)洞滑坡位于長(zhǎng)江支流青干河右岸，距河口約8.5km。如圖1所示，滑坡前緣較緩，中部和后緣稍陡，西高東低?；鲁拾夹纹?，滑體前緣高程125m，后緣高程350m，呈圓弧狀，以基巖為界;左側(cè)以基巖山脊為界，右側(cè)以陡坎臨空面為界，平均坡度15°，均寬300m?；w長(zhǎng)830m，面積24.9萬(wàn)m²，平均厚度22m，體積548萬(wàn)m³，主滑方向60°。

根據(jù)滑坡的工程地質(zhì)勘察資料，分析了滑坡地質(zhì)結(jié)構(gòu)特征，建立了三維有限元計(jì)算網(wǎng)格模型，如圖2所示。該模型沿青干河水流方向?yàn)?080m，垂直青干河水流方向?yàn)?468m，底面高程為-142m，計(jì)算區(qū)域包含滑體、滑帶和基巖，對(duì)其劃分網(wǎng)格后，共生成94090個(gè)六面體單元，105644個(gè)節(jié)點(diǎn)。

在計(jì)算模型中考慮了滑坡體的3種典型土體：滑體、滑帶和滑床（見(jiàn)圖3），相應(yīng)的物理參數(shù)如表1所示。

3巖土體的屈服準(zhǔn)則

在庫(kù)水位變化和降雨共同作用下，滑坡前緣和中部主要發(fā)生剪切破壞，后緣出現(xiàn)張拉裂隙，發(fā)生張拉破壞。巖土體發(fā)生剪切破壞的屈服準(zhǔn)則常采用廣義米賽斯準(zhǔn)則（D-P準(zhǔn)則）和莫爾一庫(kù)侖準(zhǔn)則（M-C準(zhǔn)則）。D-P準(zhǔn)則在主應(yīng)力空間上的屈服面為一圓錐面，在π平面上為圓形，不存在尖頂處的數(shù)值計(jì)算問(wèn)題。M-C準(zhǔn)則的屈服面為不規(guī)則的六角形截面的角椎體表面，在π平面上為不等角六邊形，存在尖頂和菱角。運(yùn)用現(xiàn)有技術(shù)，已能較好處理尖頂，同時(shí)M-C準(zhǔn)則在邊坡工程中得到廣泛應(yīng)用。因此，本文選用M-C準(zhǔn)則來(lái)判定巖土體是否發(fā)生剪切破壞。巖土體發(fā)生張拉破壞常采用的屈服準(zhǔn)則為最大拉應(yīng)力準(zhǔn)則。

復(fù)合準(zhǔn)則在（σ1，σ3）平面上的描述如圖4所示。圖4中，A點(diǎn)到B點(diǎn)為M-C屈服準(zhǔn)則f^s=0，其中f^s可以表示為：

4邊界條件

采用飽和非恒定滲流與應(yīng)力耦合理論進(jìn)行三門(mén)洞滑坡流固耦合計(jì)算分析，滲流邊界條件和應(yīng)力邊界條件的確定方法如下。

4.1滲流邊界條件

（1）滲流初始邊界條件。基于工程勘察資料確定坡體初始地下水位，然后采用Abaqus確定坡體內(nèi)初始滲透場(chǎng)的孔隙水壓力分布。如表1所示，采用不同巖土體滲透系數(shù)對(duì)應(yīng)的滲流速度來(lái)模擬計(jì)算模型四周邊界的透水條件。

（2）非恒定滲流邊界條件。降雨引起坡面入滲及產(chǎn)流條件由降雨強(qiáng)度和坡面巖土體的人滲率決定。當(dāng)降雨強(qiáng)度小于坡面巖土體的入滲率時(shí)，降雨產(chǎn)生的水流全部入滲至坡體;當(dāng)降雨強(qiáng)度大于坡面巖土體的入滲率時(shí)，坡面上將產(chǎn)生徑流?；诮涤暝谄旅嫒藵B產(chǎn)流的主要特點(diǎn)，坡面入滲產(chǎn)流的邊界條件可通過(guò)如下方法確定：

式中，v為雨水在坡面的滲流流速（由坡面巖土體的滲透系數(shù)決定）;q為降雨強(qiáng)度。

本文采用第一類(lèi)邊界條件來(lái)模擬庫(kù)水位變化時(shí)坡體表面的水頭變化情況。數(shù)值計(jì)算模型底邊為不透水邊界，四周及坡面不接觸降雨或庫(kù)水的邊界為自由透水面。

4.2應(yīng)力邊界條件

數(shù)值計(jì)算模型的底部邊界為法向約束，前后兩側(cè)及左右兩側(cè)均采用法向約束，坡體表面為自由邊界。

最終建立的滲流邊界和應(yīng)力邊界如圖5所示。

5滑坡變形現(xiàn)狀

三門(mén)洞滑坡坡體上共布設(shè)ZG360-ZG365等6個(gè)GPS監(jiān)測(cè)點(diǎn)，監(jiān)測(cè)時(shí)間為2006年9月至2009年12月，監(jiān)測(cè)結(jié)果如圖6-7所示。

從圖6可看出，三峽庫(kù)區(qū)庫(kù)水位每年均會(huì)經(jīng)歷不同程度的漲落，基本與各監(jiān)測(cè)點(diǎn)的位移變形變化相互對(duì)應(yīng);庫(kù)水位上漲時(shí)，滑坡變形變化不明顯，而庫(kù)水位下降時(shí)，滑坡變形變化大，且每次庫(kù)水位下降均會(huì)導(dǎo)致滑坡累計(jì)位移曲線上揚(yáng)，說(shuō)明其對(duì)滑坡穩(wěn)定性會(huì)產(chǎn)生很大影響。

從圖7可看出，ZG360點(diǎn)累積水平位移為382.2mm，ZG361點(diǎn)累積水平位移為718.7mm，ZG362點(diǎn)累積水平位移為1205.3mm，ZG363點(diǎn)累積水平位移為433.7mm，ZG364點(diǎn)累積水平位移為413.9mm，ZG365點(diǎn)累積水平位移為1557.6mm。

基于監(jiān)測(cè)數(shù)據(jù)資料，查明三峽庫(kù)區(qū)水位調(diào)度過(guò)程（見(jiàn)圖8）及當(dāng)?shù)?0a內(nèi)每月最大的連續(xù)3d降雨量（見(jiàn)圖9）。

6計(jì)算工況及結(jié)果分析

6.1計(jì)算工況

為充分表現(xiàn)庫(kù)水位變化和降雨的綜合作用效應(yīng)，將對(duì)比分析普通工況（工況1）和極端工況（工況2）條件下滑坡孔隙水壓力分布、應(yīng)力場(chǎng)、位移場(chǎng)和塑性區(qū)分布特征。其中，工況1為水庫(kù)蓄水至175m水位;工況2為庫(kù)水位從175.0m降至145.0m，并遇上50a一遇暴雨，荷載組合為自重和地表荷載。

6.2結(jié)果分析

（1）孔隙水壓力計(jì)算結(jié)果及分析。水庫(kù)蓄水至175m后，孔隙水壓力如圖10所示：滑坡前緣高程約140m處出現(xiàn)最大孔壓，數(shù)值約為0.6MPa;滑坡后緣底部孔隙水壓力主要受地下水位控制，最大孔壓值為2.60MPa（見(jiàn)圖10）。

隨著庫(kù)水位由175m降至145m，地下水位也相應(yīng)下降但存在滯后，這是由于滑坡前緣孔隙水壓力差消散需要一段時(shí)間;滑坡后緣地下水基本不受庫(kù)水位變化影響，這是由于滑坡后緣離庫(kù)水位變化帶的距離相對(duì)較遠(yuǎn)，但受強(qiáng)降雨影響，滑坡后緣底部最大孔壓值由2.60MPa增至2.83MPa（見(jiàn)圖11）。

（2）應(yīng)力計(jì)算結(jié)果及分析。水庫(kù)蓄水至175m水位后大部分區(qū)域拉、壓應(yīng)力分布變化較小。但是，由于庫(kù)水作用滑體中下部產(chǎn)生較大變化，最大主應(yīng)力的最大值為-0.20MPa，最小主應(yīng)力的最小值為-0.71MPa（圖12-13）。

當(dāng)庫(kù)水位由175m降至145m，并遇上50a一遇暴雨時(shí)，滑坡整體應(yīng)力增加，最大、最小主應(yīng)力場(chǎng)整體變化較大。滑體中最大主應(yīng)力的最大值為-1.51MPa，最小主應(yīng)力的最小值為-2.60MPa;滑體中上部最大主應(yīng)力的最大值為-1.10MPa，最小主應(yīng)力的最小值為-3.14MPa。該工況應(yīng)力場(chǎng)變化較大，主要原因是庫(kù)水下降和降雨的共同作用引起滲流場(chǎng)的改變，而滲流場(chǎng)與應(yīng)力場(chǎng)相互耦合作用，導(dǎo)致應(yīng)力場(chǎng)發(fā)生改變（見(jiàn)圖14～15）。

（3）位移計(jì)算結(jié)果及分析。水庫(kù)蓄水至175m后，滑體下部發(fā)生較大水平位移，最大值為0.015m;滑體中部水平位移相對(duì)較小，最大值為0.004m;滑體上部水平位移較小，最大值為0.001m。這主要是由于庫(kù)水對(duì)滑體產(chǎn)生浮托作用，使滑體抗滑力降低，導(dǎo)致滑體產(chǎn)生水平位移。在滑體上緣與滑床交界的部位出現(xiàn)一些突變位移，這是因?yàn)樵诖藚^(qū)域滑體與滑床巖土材料屬性和地形均產(chǎn)生了不均勻特性（見(jiàn)圖16）。此時(shí)，滑體上部發(fā)生較大向下的垂直位移，最大值為0.040m;滑體中部垂直位移相對(duì)較小，最大值為0.029m;滑體前緣產(chǎn)生向上的垂直位移，最大值為0.006m。這主要是由于庫(kù)水對(duì)滑體前緣的浮托作用，使滑體產(chǎn)生向上的垂直位移（見(jiàn)圖17）。

當(dāng)庫(kù)水位由175m降至145m，并遇50a一遇暴雨時(shí)，滑體下部發(fā)生較大水平位移，最大位移值為0.049m;滑體中部水平位移相對(duì)較小，最大位移為0.026m;滑體上部水平位移相對(duì)最小，最大位移為0.0013m。該工況水平位移較大，主要原因是庫(kù)水下降使滑體內(nèi)產(chǎn)生向外的滲透力，同時(shí)降雨作用使向外的滲透力增大，使滑體下滑力增大，導(dǎo)致滑體下部水平位移增大（見(jiàn)圖18）。此時(shí)，滑體上部發(fā)生較大向下的垂直位移，最大值為0.102m;滑體中部垂直位移相對(duì)較小，最大值為0.063m;滑體下部垂直位移相對(duì)最小，最大值為0.052m。該工況的垂直位移較大，主要原因是庫(kù)水下降引起滑體向下滑移，同時(shí)在自重和降雨作用下發(fā)生沉降，使滑體產(chǎn)生向下的垂直位移增大?；w前緣產(chǎn)生向上的垂直位移，最大值為0.021m，主要原因是：庫(kù)水對(duì)滑體前緣的浮托作用，使滑體向上的垂直位移增大（見(jiàn)圖19）。

（4）塑性區(qū)計(jì)算結(jié)果及分析。工況1和2計(jì)算過(guò)程中滑坡未出現(xiàn)塑性破壞區(qū)（見(jiàn)圖20～21），說(shuō)明滑坡整體處于較穩(wěn)定狀態(tài)。

（5）邊坡整體穩(wěn)定性分析。對(duì)比分析兩種工況的數(shù)值計(jì)算模擬結(jié)果和監(jiān)測(cè)數(shù)據(jù)后發(fā)現(xiàn)，數(shù)值模擬結(jié)果與監(jiān)測(cè)數(shù)據(jù)匹配度較高，證明了本文模擬結(jié)果的有效性。在強(qiáng)降雨和庫(kù)水位緩慢下降綜合作用下，三門(mén)洞滑坡在滑坡中后部和前緣分別出現(xiàn)局部的張拉和剪切破壞，滑體內(nèi)滲流場(chǎng)發(fā)生很大變化，導(dǎo)致滑體內(nèi)孔隙水壓力、應(yīng)力、位移和等效塑性應(yīng)變整體增大，但并未出現(xiàn)塑性破壞區(qū)，邊坡整體處于穩(wěn)定狀態(tài)。

7結(jié)語(yǔ)

本文以三峽庫(kù)區(qū)三門(mén)洞滑坡為研究對(duì)象，建立了基于Abaqus的滑坡三維計(jì)算模型、庫(kù)水位變化和降雨人滲的耦合邊界條件、莫爾一庫(kù)倫和最大拉應(yīng)力準(zhǔn)則相結(jié)合的復(fù)合準(zhǔn)則，分析了滑坡在水庫(kù)蓄水、庫(kù)水位緩慢降落及強(qiáng)降雨過(guò)程中的流固耦合及三維穩(wěn)定性問(wèn)題。

三門(mén)洞滑坡的變形破壞是在三峽水庫(kù)運(yùn)行期水位波動(dòng)條件下，耦合降雨過(guò)程長(zhǎng)期發(fā)展的結(jié)果。由于篇幅有限，本文只考慮了水庫(kù)蓄水、庫(kù)水位緩慢回落和強(qiáng)降雨等工況。若要得到邊坡變形破壞發(fā)展的真實(shí)趨勢(shì)，須考慮庫(kù)水位上升、下降、再上升、下降等反復(fù)循環(huán)的整個(gè)過(guò)程。