蠕變對河流側蝕型黃土滑坡影響的數值模擬

張吉宏 段釗 唐皓

摘?要：為了解蠕變對黃土斜坡變形破壞的影響，給河流側蝕型黃土滑坡治理提供技術參考，對涇河下游舒唐王村一帶塬邊黃土進行蠕變試驗，建立適用于該區黃土的自定義蠕變模型并進行數值模擬，分析蠕變效應對河流側蝕型黃土滑坡變形破壞的影響，結果表明：蠕變對黃土滑坡變形破壞過程具有顯著的促滑效應，蠕變效應顯著加速了滑坡土體塑性區的發展及貫通，明顯加劇了滑坡土體剪應變增量的放大與集中、增大了滑坡的水平位移。

關鍵詞：黃土滑坡;蠕變效應;本構模型;河流側蝕;數值模擬

中圖分類號：S157.1;P642.22?文獻標志碼：A

doi：10.3969/j.issn.1000-1379.2020.02.030

Abstract：In order to understand the effect of creep on the deformation and failure of loess slope， and to provide technical reference for the treatment of loess landslide caused by river erosion， the creep test of loess in the loess plateau area of the downstream region of Jinghe River was carried out. A self-defined creep model suitable for the loess in this area was established and simulated with numerical method in this paper. It is concluded that the creep effect on the deformation and failure process of loess landslides is obvious， which accelerates the development and penetration of the plastic zone， and intensifies the enlargement and concentration of the shear strain increment in the soil， and the horizontal displacement of landslide is obviously increased.

Key words： loess landslide; creep effect; constitutive model; river erosion; numerical simulation

我國滑坡災害以西部地區最為嚴重[1]。近年來黃土滑坡災害的發生次數有逐年增加的趨勢，對社會生產及群眾生活構成了嚴重威脅，甚至引發重大災難性后果[2]，因此對黃土滑坡及其變形破壞規律的研究已成為一個熱點。河流側蝕型黃土滑坡具有厚度大、滑面深、潛伏期長、不易預警等特點，國內外眾多學者[2-7]對其誘發機理、破壞模式等進行了研究，但研究中對蠕變的影響因素考慮較少。許多學者[8-14]將蠕變理論用于研究滑坡產生原因、變形規律等，并提出不同土的蠕變經驗模型，但在黃土滑坡研究中的應用較少，黃土的蠕變特性對河流側蝕黃土滑坡變形破壞規律的影響研究尚不夠深入。陜西省涇河下游沿河分布有30多千米長的黃土臺塬，受地貌影響及涇河不斷側蝕，在涇河南岸誘發了大量滑坡，筆者以目前該區仍在侵蝕過程中的舒唐王村一帶塬邊典型滑坡為例，分析其基于蠕變效應的變形破壞規律，以期為該類黃土滑坡治理提供技術參考。

1?黃土蠕變試驗及特征分析

1.1?試驗設計

蠕變試驗設備為長春試驗機研究所研制的CSS-2901TS型土體三軸流變試驗機。試驗黃土取自涇陽南塬舒唐王村滑坡坡腳，為Q2原狀黃土，天然含水率為17%。三軸蠕變試驗試樣為圓柱形試樣，尺寸為39.1 mm×80 mm。分別在50、100、150、200、250 kPa圍壓下對天然含水率試樣及飽和試樣進行不排水蠕變試驗，蠕變試驗采用單試件分級加載方法，設計強度以三軸壓縮試驗確定的峰值強度為依據[3]，法向加載7級，單級加載時間為12 h。試驗數據采用Boltzmann疊加原理進行處理。

1.2?應變—時間曲線分析

以含水率為17%的試樣試驗結果為例，不同圍壓及偏應力條件下的應變—時間關系曲線見圖1（圖中ε為應變、t為時間、各曲線代表不同偏應力狀態），可以看出：試樣在開始加載階段的變形較快，在偏應力較小時以衰減蠕變為主，隨著偏應力的逐漸加大出現等速蠕變直至加速蠕變;蠕變變形量隨偏應力的增大而增加;圍壓越大，達到相同應變量所需要的偏應力更大，破壞應變也隨之增大。飽和狀態試樣的試驗結果與此類似，在此不再贅述。

1.3?應力—應變等時關系曲線分析

同樣以含水率為17%的試樣試驗結果為例，其在不同圍壓下的應力—應變等時關系曲線見圖2，可以看出試樣的應力—應變等時關系曲線由直線段與曲線段組成，兩者之間的轉折點所對應的應力即為試樣的屈服應力σs，當應力小于σs時試樣主要表現為線性黏彈性特征，當應力大于σs時試樣的非線性黏彈性特征表現得更為明顯。飽和狀態試樣的試驗結果與此類似，在此不再贅述。

2?黃土蠕變本構模型

2.1?本構模型的建立

由試驗結果可以看出，研究區黃土蠕變變形呈現衰減穩定的特征，即由彈性變形和速率衰減的黏彈性變形共同組成蠕變變形，而Kelvin模型可以較好地描述這種變形。前已述及，研究區黃土蠕變變形在開始階段增長較快，這與單一Kelvin模型蠕變曲線有一定的差異，因此為了加快蠕變初期變形增長速度，選用五元件廣義Kelvin模型來描述，其由一個Hooke體和兩個Kelvin體串聯而成，如圖3所示（圖中：EH為Hooke體元件彈性模量，EK1、EK2為兩個Kelvin體元件彈性模量，ηK1、ηK2為兩個Kelvin體元件黏滯系數，五元件模型可表示為“H-（H/N）-（H/N）”，其中“-”表示串聯、“/”表示并聯、“H”為Hooke體、“H/N”為Kelvin體），應滿足的應力應變條件見式（1）～式（5）。

式中：σ與ε分別為總應力與總應變;σ1與ε1分別為Hooke體應力與應變;σ2、σ3、ε2、ε3分別為兩個Kelvin體的應力與應變。

當t=0（t為時間）時，ε=σ/EH，ε2=ε3=0，對式（1）～式（5）進行整理可得研究區黃土蠕變本構模型為

在式（6）的基礎上，推廣到三維狀態下的蠕變本構模型差分形式為

式中：SNij、SOij分別為一個時間增量步內的新、老應力偏量;Δeij為應變增量;ΔePij為偏應變增量;eK1，Oij、eK2，Oij分別為兩個Kelvin體在一個時間增量步內的老應變偏量;Ai=1+GKi2ηKiΔt，Bi=1-GKi2ηKiΔt，GKi、ηKi分別為Kelvin體的剪切模量與黏滯系數;a=12Ge+Δt4（1A1ηK1+1A2ηK2），b=12Ge-Δt4（1A1ηK1+1A2ηK2），Ge為Hooke體剪切模量，Δt為時間增量。

球應力的差分形式為

式中：σN0、σO0分別為一個時間增量步內的新、老球應力;K為體積模量;Δevol為體應變增量。

2.2?本構模型參數

得到蠕變本構模型后，選用ExpAssoc函數進行參數擬合，其函數形式為

將蠕變本構模型即式（6）展開，可得：

對比式（9）和式（10），可得：EH=σ/y0，EK1=σ/A1，EK2=σ/A2，ηK1=EK1t1=σt1/A1，ηK2=EK2t2=σt2/A2。

選用含水率w=17%的黃土試樣、σ3=150 kPa的試驗結果進行擬合，得到不同偏應力條件下本構模型參數值，見表1，擬合曲線見圖4。由圖4可以看出，試驗數據與擬合曲線吻合度較高，說明該自定義本構模型能很好地描述研究區黃土的蠕變特征。

3?河流側蝕型黃土滑坡變形破壞數值模擬

3.1?計算模型及參數

基于研究區滑坡的地層分布規律，將地層簡化為飽和黃土、含水率為17%的非飽和黃土及飽和砂卵石3類，且僅考慮黃土地層的蠕變。為對比分析蠕變效應對斜坡變形破壞的影響規律，屈服準則分別采用Mohr-Coulomb模型及自定義蠕變模型。采用有限差分法軟件FLAC3D進行數值模擬。

滑坡計算模型為左、右和底邊界約束，地下水位位于斜坡頂面下9 m。設定側蝕速率為每次1/3塊即5 m，把3次連續向坡腳水平側蝕15 m、綜合計算不收斂、剪應變增量區或塑性區貫通作為滑坡變形破壞判別準則?？紤]蠕變的數值模型中蠕變時間設定為1 a及5 a。計算模型見圖5，模型參數見表2、表3。

3.2?蠕變效應對斜坡變形破壞的影響

（1）塑性區的影響。Mohr-Coulomb模型和自定義蠕變模型對斜坡塑性區數值模擬結果見圖6，可以看出：在蠕變效應的影響下，土體塑性區的發展及貫通顯著加速，在同一個側蝕階段，蠕變效應顯著增大了土體塑性區范圍，在斜坡塑性剪切帶基本貫通、滑動即將發生時，蠕變效應的影響表現得更為明顯。

（2）剪應變增量的影響。Mohr-Coulomb模型及自定義蠕變模型對斜坡剪應變增量的數值模擬結果見圖7?？梢钥闯?，在蠕變作用下斜坡土體剪應變增量顯著提高，并且在滑動面附近出現峰值聚集，尤其在斜坡坡腳處剪應變增量集中現象尤為明顯，剪應變增量帶也隨之加寬。

（3）水平位移的影響。Mohr-Coulomb模型及自定義蠕變模型對斜坡水平位移數值模擬結果見圖8?？梢钥闯觯紤]了蠕變效應的斜坡水平位移明顯增大，在側蝕初期（側蝕1/3塊），蠕變模型最大水平位移已經達到11 m，而未考慮蠕變的邊坡位移僅3.6 m;隨著蠕變作用的不斷加強，到側蝕末期（側蝕3/3塊）蠕變模型最大水平位移可達到112.8 m，而未考慮蠕變的邊坡位移僅35.6 m。

綜上所述，蠕變對黃土滑坡變形破壞過程具有顯著的促滑效應，基于蠕變效應的斜坡穩定性明顯降低。側蝕發生之前斜坡土體應力平衡完成時在自重作用下已存在一定的蠕變變形，但此時土體應力水平較低，蠕變變形量相對較小，而側蝕作用改變了斜坡的初始應力狀態，斜坡需要重新分布應力場以形成新的應力平衡，在此過程中產生了新的蠕變變形，降低了斜坡的穩定性，可導致滑坡的發生。

4?結?論

對涇河下游舒唐王村一帶塬邊黃土進行蠕變試驗及數值模擬，分析了蠕變對黃土斜坡變形破壞的影響，得到以下主要結論：根據蠕變試驗的結果，自定義的五元件廣義Kelvin蠕變模型能較好地描述研究區黃土的蠕變特性;蠕變效應顯著加速了斜坡土體塑性區的發展及貫通，在同一個側蝕階段，蠕變效應顯著增大了土體塑性區范圍，這一特性在斜坡塑性剪切帶基本貫通、滑動即將發生時表現得更為明顯;考慮蠕變效應的斜坡土體剪應變增量的集中與放大明顯加劇、水平位移明顯增大，在側蝕量相同的情況下，隨著蠕變時間的延長，剪應變增量增大、剪應變增量帶變寬、斜坡水平位移增大，這一點在側蝕量增大的情況下表現得更為明顯。

參考文獻：

[1]?黃潤秋.20世紀以來中國的大型滑坡及其發生機制[J].巖石力學與工程學報，2007，26（3）：433-454.

[2]?張茂省，李同錄.黃土滑坡誘發因素及其形成機理研究[J].工程地質學報，2011，19（4）：530-540.

[3]?段釗.黃土滑坡觸發機理研究：以涇河下游南岸黃土塬區滑坡為例[D].西安：長安大學，2013：144-149.

[4]?STEFANO Dapporto， MASSIMO Rinaldi， NICOLA Casagli. Failure Mechanisms and Pore Water Pressure Conditions： Analysis of a Riverbank Along the Arno River （Central Italy） [J]. Engineering Geology， 2001， 61： 221-242.

[5]?袁湘秦，趙法鎖，段釗.河流作用誘發黃土滑坡機理[J].煤田地質與勘探，2018，46（4）：154-160.

[6]?張常亮，王阿丹，邢鮮麗，等.侵蝕作用誘發黃土滑坡的機制研究[J].巖土力學，2012，33（5）：1585-1592.

[7]?石菊松，李濱，吳樹仁，等.寶雞渭河北岸黃土塬邊大型滑坡成因機制研究[J].工程地質學報，2013，21（6）：938-949.

[8]?張先偉.結構性軟土蠕變特性及擾動狀態模型[D].長春：吉林大學，2010：205-207.

[9]?李曼，孫強.斜坡穩定性的流變模型分析[J].人民黃河，2009，31（9）：112-113.

[10]?胡貴壽，李西寧，楊世新，等.黃南Ⅷ號滑坡變形破壞機制[J].人民黃河，2009，31（3）：70-71.

[11]?葛苗苗，李寧，鄭建國，等.基于一維固結試驗的壓實黃土蠕變模型[J].巖土力學，2015，36（11）：3164-3170.

[12]?牛強.基于FLAC3D二次開發的黃土蠕變模型及滑坡失穩時間預測[D].蘭州：蘭州理工大學，2017：63-64.

[13]?張朝鵬.黃土的非線性流變本構模型及其在逆作法工程中的應用[D].西安：西安理工大學，2000：71-72.

[14]?顏斌，倪萬魁，劉海松.公路黃土地基流變模型研究[J].中國公路學報，2009，22（4）：15-19.

【責任編輯?張智民】