Deformation characteristics and reactivation mechanism of giant ancient landslide in Wumeng Mountain area: case study of Daguan ancient landslide
-
摘要: 巨型古滑坡的影响范围大、隐蔽性强，严重威胁到峡谷区城镇居民生命财产的安全。为了掌握巨型古滑坡的目前状态和未来发展趋势，以乌蒙山区云南大关古滑坡为例，采用高精度遥感解译、无人机航测、现场精细调勘查及数值模拟计算，详细分析了该滑坡的基本特征、变形过程和复活机理。调查结果显示：大关滑坡平面面积约385×104 m2，体积约2.1×108 m3。发育于3.5万年前左右。大关滑坡属于巨型古滑坡。大关古滑坡地处峡谷地带，按地形可分为三级平缓斜坡，按变形程度可分为4个变形区。大关古滑坡体及周边发育次级滑坡40处，以推移式滑坡为主。受到降雨、地震、地质环境条件和人类工程活动等多重因素影响，目前滑坡体多处发生不同程度蠕滑变形。数值模拟结果表明：在100 a一遇降雨条件下，古滑坡体上多处发生滑动，整体稳定性系数为0.98，存在整体滑动的可能；在Ⅷ度强震条件下，古滑坡体前部和中后部区域可能出现了深层滑动，整体稳定性系数为0.93。研究成果可为峡谷区此类巨型古滑坡复活研究及防灾减灾提供借鉴意义。Abstract: The giant ancient landslide has a large influence area and strong concealment, which seriously threatens the safety of lives and properties of the urban residents in the canyon area. In order to understand the current state of the giant ancient landslide in Daguan and predict its future development trend, the basic characteristics, deformation process and resurrection mechanism of the ancient landslide are analyzed by using the high-precision remote sensing interpretation, UAV survey, field fine-adjustment survey, indoor rock and soil mass tests and numerical simulation. The landslide include three gentle slopes according to its topography, which are accumulated in multiple periods during the evolution of complex slopes in geological history. The 14C dating of organic matter in the deep slip belt reveals that the landslide was developed about 35000 years ago. The plane area of the ancient landslide is about 385×104 m2, and the volume is about 2.1×108 m3. According to the deformation degree, it can be divided into four deformation zones. Under the influences of multiple factors such as rainfall, earthquake, geological environment and human engineering activities, creep deformation occurs in many parts of the slope at present. The numerical simulation results show that under the condition of one rainfall in 100 years, several secondary landslides slide on the ancient landslide, and the overall stability coefficient is 0.98, with the possibility of overall sliding. Under the strong earthquake condition, there may be deep sliding in the front and middle and rear areas of the ancient landslide, and the overall stability coefficient is 0.93. The research results may provide reference for the studies on the revival of such giant ancient landslide in the canyon area and disaster prevention and mitigation.
-
0. 引言
降雨入渗是造成滑坡最为重要的诱因之一[1]，降雨入渗会降低土体抗剪强度，增加土体重度，进而导致边坡下滑力增大，诱发边坡失稳。发展一种有效的降雨入渗模型进而在最大程度上降低降雨对边坡的危害，是降雨作用下边坡稳定性分析的关键。目前常用的边坡降雨入渗分析方法主要有[2]：①基于Richards方程的数值分析方法；②基于物理入渗模型的解析计算方法，如Green-Ampt（GA）入渗模型、Smith模型和Philip模型等。数值分析方法需要不断迭代求解降雨入渗控制方程，计算过程较为繁琐，对于复杂问题将十分耗时，并且若是网格尺寸和时间步长选择不当，将会导致计算不收敛。相较之下，以GA入渗模型为代表的解析计算方法，模型参数具有明确的物理意义，计算过程便捷，在降雨入渗下边坡稳定性研究中得到了广泛应用[3-14]。如Chen等[3]基于GA入渗模型研究了边坡降雨入渗过程，推导了边坡降雨入渗率计算表达式。张杰等[4]考虑气阻作用建立了分层假定的雨水入渗计算分析的GA模型。Zhang等[5]基于GA模型建立了可以计算时变破坏概率和失效时间的可靠度分析框架，并开发了基于电子表格的降雨型边坡破坏概率预测方法。Dou等[6]考虑土体饱和渗透系数变异性，基于GA入渗模型建立了雨水入渗-重分布分析模型。潘永亮等[7]建立了可综合考虑地下水位、初始含水率及不同降雨工况的改进GA模型。蒋水华等[8]提出了一种可考虑土体参数空间变异性的修正GA模型，并探讨了降雨入渗与多参数空间变异性作用下边坡失稳机理。
虽然上述研究大大地拓展了GA入渗模型的应用范围，但是忽略了土壤中客观存在的湿润区非饱和过渡层。大量试验发现湿润锋至土体入渗表面之间的土体并非处于完全饱和状态，在饱和区和天然区之间至少存在一个过渡层，并且过渡层内不同深度处土体体积含水量也不尽相同[9-10]。为此，张杰等[11]基于土体分层假定，推导了GA模型降雨历时与湿润锋深度之间的函数关系。胡海军等[12]在分层GA模型基础上，合理考虑湿润锋处吸力作用和水分剖面形状，分析了积水和降雨作用下非饱和重塑黄土的降雨入渗过程。Yao等[13]将饱和层与非饱和层厚度设定为入渗区深度的50%，分别对无限长边坡和马家沟滑坡进行了稳定性分析，并发现分层GA模型对研究强降雨事件下边坡稳定性问题十分有效。温馨等[14]针对黄土入渗中存在的非饱和区，提出了一种修正的GA模型来表征土体实际渗透系数随埋深的变化规律。
另一方面，受沉积、后沉积等地质作用的影响，岩土体参数存在固有的空间变异性。现有研究[8, 15-16]表明土体水力参数的空间变异性是影响降雨入渗下边坡渗流计算的关键因素。为客观评价降雨入渗下边坡失稳机理及稳定性，需要合理表征岩土体参数空间变异性的影响。然而，目前缺乏适用于考虑土体渗透系数空间变异性多层非均质边坡的降雨入渗模型。为此，提出了一种改进的Green-Ampt模型，通过4种常见的边坡工况验证该模型的有效性，进而探讨降雨入渗下多层非均质边坡含水率分布及稳定性问题。
1. 改进Green-Ampt入渗模型
1.1 Green-Ampt入渗模型
Green-Ampt（GA）入渗模型是1911年由Green等[17]提出的，最初用于解决初始干燥土体在薄层积水时的入渗问题。该模型假定雨水在入渗时存在可将湿润区和非湿润区明确分开的水平湿润锋面（图 1）。GA模型假设土体体积含水量呈矩形分布，湿润区为土体饱和层，其体积含水量为饱和体积含水量θs，非湿润区为天然层，其体积含水量为初始体积含水量θi。
![]()
图 1 经典GA入渗模型的含水率分布Figure 1. Distribution of water content of classic Green-Ampt model据此，土体入渗率f可表示为
f=−ksdhdz=−ks(hi−zf)−h0zf−0=ks(1+h0−hizf)=ks(1+h0+Sfzf)。 (1) 式中：ks为土体饱和渗透系数；h0为地面处积水水头；zf为湿润锋深度；hi为湿润锋下方的水头；Sf为湿润锋处土体概化基质吸力水头，是一个关于含水量的函数，可通过相关试验或下面的推导方程进行求解[18]：
Sf=∫1Kr(θi)ψ(Kr)dKr1−Kr(θi)=3λ+23λ+1ψb[1−Kr(θi)3λ+13λ+2]1−Kr(θi)≈3λ+23λ+1ψb。 (2) 式中：ψb为土体进气值；Kr为湿润锋处相对饱和渗透系数；λ为表征土体孔径分布特征的拟合系数。进而，可获得累积降雨入渗量I为
I = (θs−θi)zf。 (3) 可见，该模型计算简单，且模型参数具有明确的物理意义。然而，针对边坡降雨入渗问题，该模型存在如下不足：①在模拟降雨入渗时，该模型没有考虑到降雨强度R对入渗率f的影响，将会高估入渗量I，导致计算的湿润锋深度偏大。②该模型没有考虑客观存在的过渡层，大量试验研究[19-21]表明，湿润锋至土体入渗表面之间的土体并非完全处于饱和状态，即在饱和区和天然区之间应至少存在一个过渡层。如不考虑客观存在的过渡层，将会高估湿润锋以上的土重，低估边坡安全系数。③当只考虑参数（如饱和渗透系数）变异性时，传统的GA模型可以采用临界渗透系数划分降雨入渗边界条件，但是一旦考虑土体参数空间变异性，整个土层剖面的渗透系数呈空间非均匀分布，不同位置处数值大小不同，此时采用临界渗透系数划分降雨入渗边界条件不再适用[8]。
1.2 改进Green-Ampt入渗模型
为克服上述不足，十分必要对经典GA入渗模型进行改进，使其能够处理考虑土体参数空间变异性的多层非均质边坡降雨入渗问题。实际工程边坡或天然边坡通常存在一定的坡度，因而积水易形成坡面径流，其中的关键难题是如何确定不同时刻下边坡降雨入渗边界条件。鉴于采用临界渗透系数划分降雨入渗边界条件不适用考虑参数空间变异性的多层非均质边坡，为此蒋水华等[8]提出采用开始积水时间tp划分入渗边界，本文将该GA模型称为“修正GA模型”。由GA入渗模型[11, 13]可知，土体入渗率f随着累计入渗量I的增加而逐渐减小。对于倾角为α的边坡而言，当土体的累计入渗量达到某一累积入渗量Ip值时，边坡坡面开始积水[5, 7]，此时f = Rcosα，Ip的计算表达式为
Ip=(θs−θi)Sfcosα(R/ks−1)。 (4) 式中：Ip为积水时刻tp所对应的累积入渗量；α为边坡倾角；R为降雨强度。据此可得边坡开始积水的时间tp为
tp=Ip/(Rcosα)。 (5) 根据开始积水时间tp可划分入渗边界条件如下：
（1）当降雨历时t < tp时，入渗区土体处于非饱和状态，雨水将全部入渗至土体中，此时的入渗边界由流量控制，可得到入渗控制方程为
f=Rcosα=k(θ)+ks[ψr(θ)−ψr(θi)]zf。 (6) 式中：θ为土体体积含水量；ψr(θ)为湿润锋处土体相对吸力水头[8]，
ψr(θ)=∫−∞ψk(θ)ksdψ=∫−∞ψ(ψbψ)3λ + 2dψ=(ψbψ)3λ + 1ψb3λ + 1=Se3+1/λψb3λ + 1。 (7) 式中：θr为土体残余体积含水量；Se为有效饱和度，Se=(θ−θr)/(θs−θr)。
随后，依据GA模型的矩形入渗原理可推导相应的湿润锋深度为
zf=Iθ−θi = Rtcosαθ−θi。 (8) 将式（8）代入式（6）中，可得到非线性方程：
f=Rcosα=k(θ)+ks[ψr(θ)−ψr(θi)](θ−θi)Rtcosα。 (9) 求解式（9）所示的非线性方程，可得到不同降雨历时下土体体积含水量θ分布，进而结合式（8）可以计算得到对应的湿润锋深度。
（2）当降雨历时t ≥ tp时，此时边坡表面开始积水，因降雨强度超过土体允许入渗率，故土体表面入渗区处于饱和状态，此时的入渗边界由水头控制[2, 6, 8]。降雨入渗过程可分为自由入渗和积水入渗两个阶段。其入渗率f和累计入渗量I可为[11, 13]
f(t)={Rcosα(t⩽tp)ks(cosα+Sf/zf)(t>tp), (10) I=∫t0f(t)dt。 (11) 按照式（11）进行积分得到累积入渗量，进而得到入渗饱和区湿润锋深度zf为
zf=Iθs−θi。 (12) 因此，一旦确定了土体饱和渗透系数ks和相应的降雨历时t，便可针对上述两种降雨入渗边界，分别求解对应的土体湿润锋深度和含水量分布。
另一方面，大量相关试验研究[9-10, 19]表明：湿润锋至土体入渗表面之间的土体真实体积含水量并非完全处于饱和状态，即在饱和层和天然层之间应至少存在一个过渡层。根据文献[9~13]，将GA入渗模型进行分层处理，增加客观存在的过渡层，如图 2（b）所示。从图 2中可以看出，改进后的土体体积含水量分布包含饱和层、过渡层和天然层3部分。其中饱和层深度为zs，过渡层深度为zt。
![]()
图 2 改进的GA入渗模型的含水量分布Figure 2. Distribution of water content of improved Green-Ampt model文献[9，10]的相关试验数据表明，采用椭圆曲线分布来表征过渡区含水量分布具有较高的计算精度，故本文也采用椭圆曲线来描述过渡层含水量分布，由此可得修正的土体含水量分布为
θ(z)={θs (饱和层, 0≤z≤zs)θi+(θs−θi)√1−((z−zs)2zt2) (过渡层, zs≤z≤zh)θi (天然层, z>zh)。 (13) 式中：zh为整个入渗区深度，
zh=zs+zt。 (14) 随着降雨入渗的推进，饱和层、过渡层和天然层的深度都会发生变化。彭振阳等[10]通过试验发现过渡层深度与总入渗区深度zh存在一个比例关系，并且该比例与整个入渗区深度之间具有良好的线性关系：
zt=ηzh, η=azh+b。} (15) 式中：η为过渡层占整个入渗区深度的比例；a和b为拟合系数。
在过渡层修正过程中，应保持累计入渗量不变，即I=IM。I可通过式（8）和（11）计算得到，IM为增加过渡层后的累计入渗量，其值应等于饱和层入渗量Is和过渡层入渗量It之和。其中对于饱和层，入渗量Is =(θs−θi)zs；对于过渡层，由于其为椭圆曲线分布，故其入渗量等于1/4的椭圆面积It = 0.25π(θs−θi)zt，据此可得总的累计入渗量为
I=IM=Is+It=(θs−θi)zs+0.25π(θs−θi)zt。 (16) 联立式（14）~（16），便可计算得到饱和层深度zs和过渡层深度zt，进而可求得整个边坡土层的入渗区深度zh及体积含水量分布。本文将所提模型称为“改进GA模型”
2. 多层非均质边坡稳定性分析
降雨诱发滑坡通常是沿边坡浅层面发生，在剪切力作用下趋于破坏，滑动面深度远小于滑动面长度，故常采用图 3所示的无限长边坡模型分析浅层滑坡稳定性。根据莫尔-库仑破坏准则，基于极限平衡分析的无限长边坡模型安全系数计算公式为[15, 22]
FS = τfτm = c′+[(σn−ua)−σs]tanφ′Wsinαcosα。 (17) 式中：τf为土的抗剪强度；τm为沿潜在破坏面任意一点的剪应力；W为单位土体切片的重量；c′，φ′分别为土体有效黏聚力和有效内摩擦角；(σn−ua)为单位土条底部净法向力，ua为孔隙气压力。根据文献[2，6，8，15]，取ua = 0，σn=γtzcos2α，其中γt为土体有效重度，z为滑动面深度，α为无限长边坡倾角；σs为吸应力，
σs=−Seψ(θ)。 (18) 式中：ψ(θ)为基质吸力，是关于土体体积含水量θ的函数。
为了考虑降雨入渗过程中因体积含水量变化而引起的土体单位重量变化，对滑动面以上土体部分进行积分得到土体重量为
W=∫z0γtdz。 (19) 式中：土体有效重度γt[15]为
γt=γd+θ(z)γw。 (20) 式中：γd和γw分别为土体干重度和水重度。
综上，基于本文提出的改进GA模型进行降雨入渗下非均质边坡稳定性分析的步骤如下：
（1）收集并确定土体参数统计特征（包括均值、标准差、自相关函数、波动范围及概率分布等），沿垂直方向将边坡划分为H/d个等厚的均质土层（H为坡高，d为随机场单元尺寸）。
（2）基于土体参数统计特征，采用蒙特卡洛模拟（MCS）方法产生N组独立标准正态随机样本，随后基于Karhunen-Loève（KL）级数展开方法[15-16, 23]生成N组土体参数随机场实现值，并将其一一赋值给对应的土层。
（3）采用式（5）计算开始积水时间tp，并与降雨历时t进行比较，划分降雨入渗边界条件，再根据式（8），（12），（13）~（16），求解考虑客观存在过渡层的降雨入渗下边坡体积含水量分布。
（4）采用式（17）计算边坡每个土层及湿润锋处安全系数，从中可获得（H/d+1）个安全系数，并将其中最小的安全系数视作该边坡安全系数，对应的滑动面则为该边坡最危险滑动面。
3. 算例分析
下面以图 3所示的非饱和无限长边坡为例[6, 24]，验证提出的改进GA模型的有效性。边坡坡高H = 3 m，倾角α = 50°，下部为不透水基岩，降雨强度R恒定为5 mm/h。采用Brook等[25]提出的土水特征曲线模型描述土体渗透系数和基质吸力与体积含水量之间的函数关系为
k(θ)=ksSe3+2/λ, (21) ψ(θ)=ψbSe−1/λ。 (22) 式中：k(θ)为土体渗透系数。根据文献[6，24]，得到土体物理力学参数取值如表 1所示。另外根据文献[10]的试验结果，取式（15）中的系数a和b分别为-0.003，0.8712。
表 1 土体物理力学参数取值Table 1. Values of physical parameters of soil计算参数 取值 计算参数 取值 饱和渗透系数ks/(mm·h-1) 3 有效黏聚力c′/kPa 5 有效内摩擦角φ′/(°) 28 进气值ψb/kPa 2.752 土体干重度γ /(kN·m-3) 16.217 湿润锋处土体概化基质吸力Sf /mm 424.3 饱和含水量θs 0.335 残余含水量θr 0.068 初始含水量θi 0.148 土体孔隙分布特征参数λ 0.319 通过以下4种常见的边坡工况，验证采用本文提出的改进GA模型进行降雨入渗下边坡体积含水量计算及稳定性分析的有效性：①均质边坡；②两层边坡（上层渗透性大、下层渗透性小）；③两层边坡（上层渗透性小、下层渗透性大）；④考虑任意一次参数随机场典型实现的非均质边坡。同时，选用Richards方程数值解和修正GA入渗模型[8]计算结果进行对比分析。其中，Richards方程数值解采用一维流动模型来模拟边坡垂直入渗过程，描述一维水流在边坡上的入渗过程控制方程为[26]
∂θ∂t=∂∂z[k(θ)(∂h∂z−cosα)]。 (23) 式中：θ∈[θr,θs]；h为压力水头。
本文借助Hydrus-1D软件[27]对式（23）所示Richards方程进行数值求解，可得到不同降雨历时t下不同深度z处压力水头h(z, t)和土体体积含水量θ(z, t)分布。本文将Richards方程数值解视作精确解进行判断对比分析。
3.1 工况1：均质边坡
基于表 1中的土体参数均值，对图 3所示均质边坡模型进行降雨入渗分析。图 4比较了采用修正GA模型和提出的改进GA模型计算的体积含水量分布和安全系数及Richards方程数值解。
![]()
图 4 不同降雨历时下边坡体积含水量和安全系数分布的比较（工况1）Figure 4. Comparison of distribution of water content and factors of safety of slope under different rainfall durations (Scenario Ⅰ)由图 4可知，采用修正GA模型计算的体积含水量分布与提出的改进GA模型计算结果和Richards方程数值解均存在较大差别，这是由于修正GA模型假定雨水呈矩形入渗，没有考虑客观存在的过渡层。在此基础上进行边坡稳定性分析，其边坡FS沿坡体深度的变化如图 4（b）所示。由图 4可知，随着降雨的进行，雨水入渗区也不断下移，饱和层内土体基质吸力基本丧失，导致入渗区内边坡FS随着埋深的增大而逐渐减小。到达过渡层后，由于土体并非完全处于饱和状态，基质吸力的作用未完全丧失，边坡FS先下降而随后又呈增加趋势。接着，虽然土体天然层内土体含水量和基质吸力均保持初始状态，但由于土体深度不断增加，导致土条自重逐渐增加，坡体滑动力增大，此时的边坡FS表现为随着埋深的增大而逐渐降低，并在埋深z = 3.0 m的不透水层处达到最小值。
表 2比较了不同方法计算的降雨历时20，36，60 h下边坡安全系数和最危险滑动面深度。表中，入渗区是指雨水入渗至坡体内的区域，包括饱和层和过渡层（例如降雨历时60 h，图 4（a）中入渗区可视为距离坡面0~1.25 m区域）；整个边坡则包含入渗区和天然区。将其分开进行处理，可分析不同降雨历时下边坡最危险滑动面位置和确定边坡失稳是否因降雨所诱发。若边坡最危险滑动面位于入渗区，可认为边坡失稳因降雨造成。由表 2可知，边坡安全系数随着降雨历时的增加而减小。相较于修正GA模型，提出的改进GA模型与基于Richards方程数值方法计算的安全系数吻合。如降雨历时60 h，提出的改进GA模型计算的边坡安全系数为1.22，与基于Richards方程数值方法计算的1.19差别很小。表 2中采用修正GA模型计算的边坡安全系数偏小的原因是，该模型假定雨水呈矩形入渗，导致计算的体积含水量θ偏大，进而高估了土条重量W，造成边坡安全系数偏小。相较之下，提出的改进GA模型因考虑了客观存在的过渡层，计算的边坡体积含水量、安全系数以及最危险滑动面位置均与Richards方程数值解吻合，这表明提出的改进GA模型能够有效分析均质边坡降雨入渗问题。此外，降雨历时20，36 h时，3种方法均识别出边坡最危险滑动面位于z = 3.0 m不透水层处。降雨历时60 h，相较于修正GA模型计算结果（1.03 m），本文提出的改进GA模型结果（0.99 m）与精确解（0.99 m）更吻合，说明提出方法能够准确识别边坡最危险滑面位置。
表 2 边坡安全系数和最危险滑动面深度的比较（工况1）Table 2. Comparison of factors of safety and depths of critical slip surface of slope (Scenario Ⅰ)类别 降雨历时/h Richards方程数值解 修正GA模型 提出的改进GA模型 入渗区安全系数 20 2.79 2.35 2.78 36 1.72 1.53 1.74 60 1.19 1.10 1.22 整个边坡安全系数 20 1.36 1.36 1.36 36 1.34 1.34 1.34 60 1.19 1.10 1.22 最危险滑动面深度 20 3.00 m 3.00 m 3.00 m 36 3.00 m 3.00 m 3.00 m 60 0.99 m 1.03 m 0.99 m 3.2 工况2：两层边坡（上层渗透性大、下层渗透性小）
为了进一步说明提出的改进GA模型的有效性，本节对两层边坡（上层渗透性大、下层渗透性小）进行降雨入渗分析。第一层土饱和渗透系数为3.5 mm/h，土层厚度为0.5 m；第二层土饱和渗透系数为3 mm/h，土层厚度为2.5 m。两个土层的其他参数均相同，见表 1。图 5比较了工况2下采用上述2种模型和Richards方程数值方法计算的土体体积含水量分布。由图 5可知，提出的改进GA模型计算结果同样与Richards方程数值解吻合。表 3比较了采用不同方法计算的降雨历时20，36，60 h安全系数及最危险滑动面深度。可见，提出的改进GA模型与Richards方程数值方法计算的边坡安全系数依然吻合。如降雨历时20 h，Richards方程数值方法计算的入渗区安全系数为2.83，提出的改进GA模型结果为2.81，而修正GA模型结果为2.39，进一步体现了提出的改进GA模型的有效性及考虑客观存在的过渡层的必要性。
![]()
图 5 不同降雨历时下边坡体积含水量分布的比较（工况2）Figure 5. Comparison of distribution of water content of slope under different rainfall durations (Scenario Ⅱ)表 3 边坡安全系数和最危险滑动面深度的比较（工况2）Table 3. Comparison of factors of safety and depths of critical slip surface of slope (Scenario Ⅱ)类别 降雨历时/h Richards方程数值解 修正GA模型 提出的改进GA模型 入渗区安全系数 20 2.83 2.39 2.81 36 1.71 1.55 1.76 60 1.19 1.10 1.22 整个边坡安全系数 20 1.36 1.36 1.36 36 1.34 1.34 1.34 60 1.19 1.10 1.22 最危险滑动面深度 20 3.00 m 3.00 m 3.00 m 36 3.00 m 3.00 m 3.00 m 60 0.98 m 1.03 m 0.99 m 3.3 工况3：两层边坡（上层渗透性小、下层渗透性大）
本节进一步分析两层边坡（上层渗透性小、下层渗透性大）来验证改进GA模型的有效性。第一层土饱和渗透系数为3 mm/h，土层厚度为0.5 m；第2层土饱和渗透系数为3.5 mm/h，土层厚度为2.5 m。这两个土层的其他参数均相同，见表 1。图 6比较了采用上述2种模型和Richards方程数值方法计算的边坡体积含水量分布。表 4比较了基于不同方法计算的边坡安全系数。与工况2类似，提出的改进GA模型与Richards方程数值方法计算的安全系数一致，说明提出的改进GA模型可有效求解多层边坡降雨入渗问题。
![]()
图 6 不同降雨历时下边坡体积含水量分布的比较（工况3）Figure 6. Comparison of distribution of water content of slope under different rainfall durations (Scenario Ⅲ)表 4 边坡安全系数和最危险滑动面深度的比较（工况3）Table 4. Comparison of factors of safety and depths of critical slip surface of slope (Scenario Ⅲ)类别 降雨历时/h Richards方程数值解 修正GA模型 提出的改进GA模型 入渗区安全系数 20 2.79 2.35 2.78 36 1.72 1.53 1.74 60 1.19 1.10 1.22 整个边坡安全系数 20 1.36 1.36 1.36 36 1.34 1.34 1.34 60 1.19 1.10 1.22 最危险滑动面深度 20 3.00 m 3.00 m 3.00 m 36 3.00 m 3.00 m 3.00 m 60 1.01 m 1.04 m 1.00 m 3.4 工况4：考虑一次参数随机场典型实现的非均质边坡
工程实际中，土体饱和渗透系数ks呈空间非均质分布。本文将边坡沿垂直方向上划分为60个均质土层，厚度为d = 0.05 m[6]。基于KL级数展开方法来离散ks随机场，利用高斯型自相关函数模拟土体ks的空间自相关性，垂直自相关距离lv取0.5 m。当KL级数展开截断项数取6时，计算的期望能比率因子ε为95.67%，可以满足随机场离散精度要求（大于95%）。图 7给出了土体饱和渗透系数随机场的一次典型实现，图中颜色较深部分表示ks值较大区域，颜色较浅部分表示ks值较小区域。
![]()
图 7 饱和渗透系数随机场一次典型实现分布图Figure 7. Typical realization of random filed of hydraulic conductivity图 8比较了不同方法计算的降雨历时8，36，60 h边坡体积含水量分布。由图 8可知，在3个不同降雨历时下，提出的改进GA模型计算的体积含水量分布和Richards方程数值解同样吻合，入渗区和整个边坡的安全系数和最危险滑动面深度几乎相同（见表 5），同样说明了基于提出的改进GA模型进行多层非均质边坡降雨入渗分析的有效性。由表 5可知，降雨历时8，36 h，由这2种模型和Richards方程数值方法计算的整个边坡安全系数皆小于入渗区安全系数，此时边坡最危险滑动面均位于z = 3 m的不透水层处。然而，在降雨历时60 h时，提出的改进GA模型和Richards方程数值方法均表明边坡最危险滑动面位于z = 3 m不透水层处。而基于修正GA模型计算的整个边坡安全系数等于入渗区安全系数，即修正GA模型认为此时最危险滑动面处于入渗区（距离坡面约0.85 m），显然错误评估了边坡失稳模式。这是由于该典型实现下土层表面渗透性低，湿润锋推进缓慢，虽然湿润锋处安全系数在不断减小，但是始终大于不透水层处安全系数，故而边坡仍然沿着不透水层处发生失稳破坏。综上，采用修正GA模型会造成对边坡失稳机理的错误认识，相比之下采用提出的改进GA模型能够准确识别最危险滑动面位置，获得更合理的边坡安全系数。另外需要说明的是，提出的改进Green-Ampt模型和Richards方程数值解得到的过渡层分布存在误差，并且该误差随着降雨历时的增加而增大，见图 4（a），5，6，这是因为本文没有考虑式（15）中系数a和b随着降雨历时的增加有可能发生变化。尽管如此，两者计算的湿润锋深度随着降雨历时的增加均保持一致。
![]()
图 8 工况4不同降雨历时下边坡体积含水量分布比较Figure 8. Comparison of distribution of water content of slope under different rainfall durations (Scenario Ⅳ)表 5 边坡安全系数和最危险滑动面深度的比较（工况4）Table 5. Comparison of factors of safety and depths of critical slip surface of slope (Scenario Ⅳ)类别 降雨历时/ h Richards方程数值解 修正GA模型 提出的改进GA模型 入渗区安全系数 8 6.58 5.62 6.43 36 1.95 1.72 1.98 60 1.40 1.24 1.39 整个边坡安全系数 8 1.37 1.37 1.37 36 1.35 1.35 1.35 60 1.33 1.24 1.33 最危险滑动面深度 20 3.00 m 3.00 m 3.00 m 36 3.00 m 3.00 m 3.00 m 60 h 3.00 m 0.85 m 3.00 m 为了说明这4种工况下边坡体积含水量分布之间的差别，图 9比较了这4种工况降雨历时36，60 h的边坡含水量分布曲线。由图 9（a）可知，工况1，2，3的边坡含水量分布曲线重合，而工况4入渗区较小。这是由于工况4边坡表面渗透性小，雨水入渗率较低，进而导致入渗区小。工况1，3由于上层饱和渗透系数相同，故其在降雨历时36 h下含水量分布曲线重合。工况2与工况1，3的边坡含水量分布曲线重合的原因是，在降雨初期，降雨历时未超过积水时刻tp，雨水入渗率等于Rcosα，此时的入渗量仅与降雨历时相关，与土体饱和渗透系数大小无关。因此，工况1，2，3在降雨历时36 h入渗量相差不大，对应的边坡含水量分布曲线基本重合。相比之下，降雨历时60 h，工况4由于其自身的渗透性较低，故入渗区小；工况3下由于下层土体渗透性大，导致入渗量增大，进而入渗深度更大；工况1，2由于下层土体的饱和渗透系数一致，故这两种工况不仅入渗率和入渗量相等，而且边坡体积含水量分布也相同，见图 9（b）。
![]()
图 9 不同工况下边坡体积含水量分布的比较Figure 9. Comparison of distribution of water content of slope for different scenarios4. 结论
考虑降雨入渗作用下边坡土体中客观存在的过渡层，提出了一种适用于多层非均质边坡降雨入渗分析的改进GA模型，为合理确定不同降雨历时下多层非均质边坡的入渗区深度和土体体积含水量分布奠定了基础。同时，以无限长边坡模型为例，进行降雨入渗下边坡稳定性分析，并与采用Richards方程数值解和修正GA模型计算的边坡安全系数和最危险滑动面深度进行了对比分析，验证了提出的改进GA模型的有效性。
（1）修正GA模型局限于雨水呈矩形入渗，没有考虑雨水入渗过程中土体客观存在的过渡层，导致降雨入渗过程中单位土体重度被高估，计算的边坡安全系数偏小，造成对边坡失稳机理的错误认识。
（2）提出的改进GA模型因考虑了边坡土体入渗区内客观存在的过渡层，故计算的土体含水量分布、安全系数和滑动面分布不仅贴近客观实际，而且与Richards方程数值解吻合，能够较好地解决考虑土体饱和渗透系数空间变异性的多层非均质边坡降雨入渗问题。
（3）相比于修正GA模型，提出的改进GA模型能够准确识别非均质边坡最危险滑动面位置，为降雨型滑坡加固设计及风险防控提供了一条有效途径。如工况4中，降雨历时60 h时，修正GA模型认为最危险滑动面位于入渗区，而改进GA模型和Richards方程计算结果均证明边坡最危险滑动面位于不透水层处。
（4）工程实际中存在多种环境因素对边坡降雨入渗过程具有重要的影响，包括边坡初始状态（如含水量，土壤类型等）的不均匀分布、边坡倾角大小、覆土厚度和降雨雨型等，融入多种环境因素影响的边坡降雨入渗模型还有待进一步研究。
-
![]()
图 4 大关古滑坡工程地质剖面图1 - 1′
Figure 4. Engineering geological profile of Daguan ancient landslide
![]()
图 7 职业中学滑坡工程地质剖面图3 - 3′
Figure 7. Engineering geology profile of landslide 3 - 3′ of Daguan Vocational High School
![]()
图 8 大关职业中学滑坡地表变形特征
Figure 8. Characteristics of surface deformation of landslide Daguan Vocational High School
![]()
图 9 牌坊步行街滑坡工程地质剖面图4 - 4'
Figure 9. Engineering geology profile of Jiayuan community 4 - 4'
![]()
图 11 自然重力状态下滑坡位移和剪应变增量图
Figure 11. Displacement and shear strain increment under natural gravity
![]()
图 12 50 a一遇降雨状态下滑坡位移和剪应变增量图
Figure 12. Displacement and shear strain increment under a 50-year rainfall return period
![]()
图 13 100 a一遇降雨状态下滑坡位移和剪应变增量图
Figure 13. Displacement and shear strain increment under a 100-year rainfall return period
![]()
图 14 Ⅶ度地震烈度状态下滑坡位移和剪应变增量图
Figure 14. Displacement and shear strain increment under Ⅶ earthquake intensity
![]()
图 15 Ⅷ度地震烈度状态下滑坡位移和剪应变增量图
Figure 15. Displacement and shear strain increment under Ⅷ earthquake intensity
表 1 数值模型岩土体物理力学参数
Table 1 Physical mechanical parameters of rock mass
岩土体 体积模
量/MPa剪切模
量/MPa重度/
(kN·m-3)内摩擦
角/(°)黏聚力/
kPa含砾粉质黏土 1.47×103 1×104 23.7 27 40 含块石粉质黏土 2.0×104 1.2×104 24.6 32 60 基岩 2.5×106 2.0×106 26.3 45 120 
下载: 导出CSV
表 2 不同降雨重现期下岩土体的物理力学参数
Table 2 Physical mechanical parameters of rock mass under.different rainfall recurrence periods
降雨
工况岩土体 体积
模量/
MPa剪切
模量/
MPa重度/
(kN·m-3)内摩擦
角/(°)黏聚力/
kPa50 a
一遇含砾粉质
黏土735 5000 11.85 13.5 20 含块石粉质
黏土20000 12000 24.60 32.0 60 基岩 2.5×106 2×106 26.30 45.0 120 100 a
一遇含砾粉质
黏土735 5000 11.85 13.5 20 含块石粉质
黏土16000 9600 19.68 25.6 48 基岩 2.5×106 2×106 26.30 45.0 120
下载: 导出CSV
-
[1] 张永双, 吴瑞安, 郭长宝, 等. 古滑坡复活问题研究进展与展望[J]. 地球科学进展, 2018, 33(7): 728-740. ZHANG Yongshuang, WU Ruian, GUO Changbao, et al. Research progress and prospect on reactivation of ancient landslides[J]. Advances in Earth Science, 2018, 33(7): 728-740. (in Chinese)
[2] 何坤, 胡卸文, 马国涛, 等. 四川省盐源玻璃村特大型玄武岩古滑坡复活机制[J]. 岩土力学, 2020, 41(10): 3443-3455. HE Kun, HU Xiewen, MA Guotao, et al. The reactivated mechanism of Boli Village giant ancient basalt landslide in Yanyuan, Sichuan[J]. Rock and Soil Mechanics, 2020, 41(10): 3443-3455. (in Chinese)
[3] LOPEZ SAEZ J, CORONA C, STOFFEL M, et al. Probability maps of landslide reactivation derived from tree-ring records: Pra Bellon landslide, southern French Alps[J]. Geomorphology, 2012, 138(1): 189-202. doi: 10.1016/j.geomorph.2011.08.034
[4] 廖秋林, 李晓, 李守定, 等. 三峡库区千将坪滑坡的发生、地质地貌特征、成因及滑坡判据研究[J]. 岩石力学与工程学报, 2005, 24(17): 3146-3153. doi: 10.3321/j.issn:1000-6915.2005.17.023 LIAO Qiulin, LI Xiao, LI Shouding, et al. Occurrence, geology and geomorphy characteristics and origin of qianjiangping landslide in Three Gorges reservoir area and study on ancient landslide criterion[J]. Chinese Journal of Rock Mechanics and Engineering, 2005, 24(17): 3146-3153. (in Chinese) doi: 10.3321/j.issn:1000-6915.2005.17.023
[5] 王占巍, 赵发睿, 谢文苹, 等. 青海省高家湾滑坡的形成条件分析及稳定性评价[J]. 水土保持通报, 2020, 40(3): 81-87. WANG Zhanwei, ZHAO Farui, XIE Wenping, et al. Formation condition analysis and stability evaluation of gaojiawan landslide in Qinghai Province[J]. Bulletin of Soil and Water Conservation, 2020, 40(3): 81-87. (in Chinese)
[6] 祝建, 雷英, 赵杰. 西藏樟木口岸特大型古滑坡形成机理分析[J]. 水文地质工程地质, 2008, 35(1): 49-52. doi: 10.3969/j.issn.1000-3665.2008.01.011 ZHU Jian, LEI Ying, ZHAO Jie. Mechanism analysis of the outsized ancient landslide of Zhangmu port in Tibet[J]. Hydrogeology & Engineering Geology, 2008, 35(1): 49-52. (in Chinese) doi: 10.3969/j.issn.1000-3665.2008.01.011
[7] 黄晓虎, 易武, 龚超, 等. 开挖致使古滑坡复活变形机理研究[J]. 岩土工程学报, 2020, 42(7): 1276-1285. doi: 10.11779/CJGE202007011 HUANG Xiaohu, YI Wu, GONG Chao, et al. Reactivation and deformation mechanism of ancient landslides by excavation[J]. Chinese Journal of Geotechnical Engineering, 2020, 42(7): 1276-1285. (in Chinese) doi: 10.11779/CJGE202007011
[8] 殷跃平, 朱继良, 杨胜元. 贵州关岭大寨高速远程滑坡—碎屑流研究[J]. 工程地质学报, 2010, 18(4): 445-454. doi: 10.3969/j.issn.1004-9665.2010.04.002 YIN Yueping, ZHU Jiliang, YANG Shengyuan. Investigation of a high speed and long run-out rockslide-debris flow at Dazhai in Guanling of Guizhou Province[J]. Journal of Engineering Geology, 2010, 18(4): 445-454. (in Chinese) doi: 10.3969/j.issn.1004-9665.2010.04.002
[9] 殷跃平, 王文沛, 张楠, 等. 强震区高位滑坡远程灾害特征研究: 以四川茂县新磨滑坡为例[J]. 中国地质, 2017, 44(5): 827-841. YIN Yueping, WANG Wenpei, ZHANG Nan, et al. Long runout geological disaster initiated by the ridge-top rockslide in a strong earthquake area: a case study of the Xinmo landslide in Maoxian County, Sichuan Province[J]. Geology in China, 2017, 44(5): 827-841. (in Chinese)
[10] 殷跃平. 斜倾厚层山体滑坡视向滑动机制研究: 以重庆武隆鸡尾山滑坡为例[J]. 岩石力学与工程学报, 2010, 29(2): 217-226. YIN Yueping. Mechanism of apparent dip slide of inclined bedding rockslide: a case study of Jiweishan rockslide in Wulong, Chongqing[J]. Chinese Journal of Rock Mechanics and Engineering, 2010, 29(2): 217-226. (in Chinese)
[11] HOLTZ R D. KOVACS W D. An Introduction to Geotechnical Engineering[M]. Englewood Cliffs: Prentice Hall, 1985: 543-545.
[12] STARK T D, CHOI H, MCCONE S. Drained shear strength parameters for analysis of landslides[J]. Journal of Geotechnical and Geoenvironmental Engineering, 2005, 131(5): 575-588. doi: 10.1061/(ASCE)1090-0241(2005)131:5(575)
[13] BURDA J, HARTVICH F, VALENTA J, et al. Climate- induced landslide reactivation at the edge of the Most Basin (Czech Republic): progress towards better landslide prediction[J]. Natural Hazards and Earth System Sciences, 2013, 13(2): 361-374. doi: 10.5194/nhess-13-361-2013
[14] 杨成业, 张涛, 高贵, 等. SBAS-InSAR技术在西藏江达县金沙江流域典型巨型滑坡变形监测中的应用[J]. 中国地质灾害与防治学报, 2022, 33(3): 94-105. YANG Chengye, ZHANG Tao, GAO Gui, et al. Application of SBAS-InSAR technology in monitoring of ground deformation of representative giant landslides in Jinsha River basin, Jiangda County, Tibet[J]. The Chinese Journal of Geological Hazard and Control, 2022, 33(3): 94-105. (in Chinese)
[15] NOTTI D, GALVE J P, MATEOS R M, et al. Human-induced coastal landslide reactivation. Monitoring by PSInSAR techniques and urban damage survey (SE Spain)[J]. Landslides, 2015, 12(5): 1007-1014. doi: 10.1007/s10346-015-0612-3
[16] VAŘILOVÁ Z, KROPÁČEK J, ZVELEBIL J, et al. Reactivation of mass movements in Dessie graben, the example of an active landslide area in the Ethiopian Highlands[J]. Landslides, 2015, 12(5): 985-996. doi: 10.1007/s10346-015-0613-2
[17] GENERALI M, PIZZIOLO M. Application of a geomorphologic-heuristic model to estimate the landslides reactivation likelihood in the Emilia-Romagna Region (Italy)[M]//Engineering Geology for Society and Territory - Volume 2. Cham: Springer International Publishing, 2015: 205-209.
[18] 郭健, 许模, 赵勇, 等. 黑水河库区某古滑坡形成及复活机制[J]. 成都理工大学学报(自然科学版), 2013, 40(6): 721-728. GUO Jian, XU Mo, ZHAO Yong, et al. Formation and reactivation mechanism of an ancient landslide in Heishui reservoir of Minjiang River, China[J]. Journal of Chengdu University of Technology (Science & Technology Edition), 2013, 40(6): 721-728. (in Chinese)
[19] 工程岩体分级标准: GB/T 50218—2014[S]. 北京: 中国计划出版社, 2014. Standard for Engineering Classification of Rock Masses: GB/T 50218—2014[S]. Beijing: China Planning Press, 2014. (in Chinese)
[20] 杞明辉, 许迎杰, 杞磊, 等. 云南省极端最大日雨量时空分布特征及其重现期估算[J]. 气候变化研究快报, 2019, 8(4): 10. QI Minghui, XU Yingjie, QI Lei, et al. Temporal and spatial distributions of maximum daily precipitation and recurrence periods in Yunnan Province[J]. Climate Change Research Letters, 2019, 8(4), 561-570. (in Chinese)
[21] 朱赛楠, 魏英娟, 王平, 等. 大型单斜层状基岩滑坡变形特征与失稳机制研究: 以重庆石柱县龙井滑坡为例[J]. 岩石力学与工程学报, 2021, 40(4): 739-750. ZHU Sainan, WEI Yingjuan, WANG Ping, et al. Research on deformation characteristics and instability mechanisms of large monoclinal layered bedrock landslides: a case study of the Longjing landslide in Shizhu County, Chongqing[J]. Chinese Journal of Rock Mechanics and Engineering, 2021, 40(4): 739-750. (in Chinese)
[22] 朱赛楠, 殷跃平, 王猛, 等. 金沙江结合带高位远程滑坡失稳机理及减灾对策研究: 以金沙江色拉滑坡为例[J]. 岩土工程学报, 2021, 43(4): 688-697. doi: 10.11779/CJGE202104011 ZHU Sainan, YIN Yueping, WANG Meng, et al. Instability mechanism and disaster mitigation measures of long-distance landslide at high location in Jinsha River junction zone: case study of Sela landslide in Jinsha River, Tibet[J]. Chinese Journal of Geotechnical Engineering, 2021, 43(4): 688-697. (in Chinese) doi: 10.11779/CJGE202104011
[23] 吴瑞安, 马海善, 张俊才, 等. 金沙江上游沃达滑坡发育特征与堵江危险性分析[J]. 水文地质工程地质, 2021, 48(5): 120-128. WU Ruian, MA Haishan, ZHANG Juncai, et al. Developmental characteristics and damming river risk of the Woda landslide in the upper reaches of the Jinshajiang River[J]. Hydrogeology & Engineering Geology, 2021, 48(5): 120-128. (in Chinese)
-
期刊类型引用(5)
1. 郭旭辉，朱鸿鹄，吴冰，高宇新，胡乐乐，曹鼎峰. 基于人工神经网络的黄土含水率光纤被动感测技术研究. 岩土力学. 2025(02): 653-664 . 
百度学术
2. 高磊，袁泽，王勤，高明军. 基于OFDR技术的土体含水率模型研究. 工程地质学报. 2025(01): 96-105 . 百度学术
3. 严辉，林沛元. 深圳市岩溶地层标准贯入击数神经网络模型. 地质科技通报. 2025(02): 305-321 . 百度学术
4. 李杰，朱鸿鹄，吴冰，刘喜凤，王家琛，曹鼎峰，施斌. 下蜀土降雨入渗光纤感测及渗透系数估算研究. 工程地质学报. 2024(02): 601-611 . 百度学术
5. 朱鸿鹄. 工程地质界面：从多元表征到演化机理. 地质科技通报. 2023(01): 1-19 . 百度学术
其他类型引用(0)
-
其他相关附件
计量
- 文章访问数: 373
- HTML全文浏览量: 59
- PDF下载量: 113
- 被引次数: 5
