Simulation of underground leakage process of soil on Karst slopes
-
摘要: 西南喀斯特地区石漠化现象严重,其独特的地表地下二元水文地质结构导致了表层岩溶带水土流失不仅表现为坡面水土流失,同时还存在水土向地下漏失的现象。为揭示岩溶坡面土壤地下漏失过程的规律和机理,利用物理模拟试验和CFD-DEM数值模拟相结合的方法,分析坡面土壤地下漏失的过程中土壤颗粒运移规律和特征。研究结果表明随着降雨的历程,坡面不同层位的土壤颗粒运移速度由坡顶到坡底逐渐增大的趋势转变为整个坡面土壤颗粒运移速度趋于一致,从剖面上看,土壤颗粒运移速度从降雨初期的表层 > 中层 > 底层转变为降雨后期的速度趋于一致;随着降雨历程,坡面土壤由地表流失出口和地下漏失口率先侵蚀,随后坡面形成侵蚀细沟并加剧侵蚀最终发展到坡面整体滑移的过程;随着降雨的进行,坡面颗粒与坡面的接触数先显著降低后小幅增大,裂隙口下缘的颗粒与坡面接触数逐渐增大;此外,在降雨过程中,坡面上的颗粒与坡面在Y方向的接触力整体上始终较小,表明降雨过程中颗粒与坡面的相互作用力较弱。
-
关键词:
- 石漠化 /
- 模拟试验 /
- CFD-DEM数值模拟 /
- 土壤漏失
Abstract: The karst rocky desertification is serious in the karst areas of southwest China. Its unique surface and underground dual-hydrogeological structure have led to the water and soil loss in the surface karst zones, which is manifested as slope soil loss and the phenomenon of water and soil leakage to the ground. In order to reveal the law and mechanism of underground leakage of soil on karst slopes, the physical simulation tests and CFD-DEM numerical simulations are used to analyze the law and characteristics of soil particle movement in the process of underground leakage of soil on slopes. The applicability and accuracy of the CFD-DEM numerical model are verified by comparing and analyzing the underground runoff and sediment yield process and the proportional characteristics of the runoff and sediment yield of the physical tests and the numerical model. The results show that the movement speed of soil particles in different layers on the slope surface gradually increases from the top to the bottom of the slope, and it tends to be consistent throughout the slope. From the profile, the movement speed of soil particles changes from the surface layer > the middle layer > the bottom layer at the beginning of rainfall to the same speed at the end of rainfall. With the rainfall process, the soil on the slope is first eroded from the surface runoff outlet and underground leak, and then forms erosion rills on the slope, which further aggravates the process of erosion and finally develops to the overall slope slip. During the rainfall process, the contact number between particles on the slope first significantly decreases and then slightly increases. The contact number between particles on the lower edge of the crack and the slope gradually increases. The contact force between particles and the slope in the Y direction is small, indicating that the contact force between the particles and the slope is weak during the rainfall process. -
0. 引言
高桩码头是一种典型的沿海水工建筑物。高桩码头的岸坡在地震期间可能发生整体失稳、产生较大的侧向位移进而导致结构损坏。
在进行高桩码头抗震设计时,需对岸坡地震稳定性进行验算,并可采用拟静力法[1]开展。在拟静力分析方法中,地震引起的荷载通过等效静态水平力施加,边坡的拟静力地震稳定性常用临界地震系数表示。针对边坡拟静力地震稳定性问题,Sarma[2]在1973年提出了一种使用极限平衡法确定边坡临界地震系数的方法。Chang等[3]采用极限分析法,确定了平面和对数螺旋两种破坏机制下边坡临界地震系数的严格上限解。Sloan[4-5]基于有限元法将土体离散进行极限分析,并开发了适用于边坡稳定性分析的数值极限分析程序。Loukidis等[6]使用有限元极限分析方法,评估了地震作用下边坡的稳定性。
考虑排桩稳定效果的高桩码头岸坡稳定性分析可采用有限元方法开展,但其计算较为复杂。Ito等[7-8]使用塑性变形理论计算抗滑桩受到的极限侧向力,提出了多排抗滑桩边坡的设计方法。Lee等[9]基于非耦合公式,分别考虑桩的响应和边坡稳定性,使用简化Bishop法寻找边坡的临界滑动面。Hassiotis等[10]推广摩擦圆法,将反作用力纳入到边坡稳定性分析中,从而确定了边坡临界滑动面和安全系数。Ausilio[11]采用极限分析上限法分析抗滑桩边坡稳定性,得出了给定安全系数条件下所需抗滑力的理论表达式。
针对考虑排桩加固的高桩码头岸坡稳定性,使用梯形分布抗滑力考虑排桩对岸坡的稳定作用,假设临界破坏面为穿过坡脚的对数螺旋面,基于极限分析上限法推导了岸坡临界地震系数的表达式。通过变量优化,得到临界地震系数最小上限解和对应的临界破坏面。进一步地,采用强度折减法计算给定水平地震加速度下高桩码头岸坡的安全系数,并分析不同因素对岸坡地震稳定安全系数和破坏模式的影响规律。
1. 高桩码头岸坡临界地震系数
根据极限分析上限定理,对于任何假定的破坏机制,若外力功率超过内部能量耗散速率,则土体不能承受施加的荷载。
∫STividS+∫VXividV≤∫VσijϵijdV i,j=1,2,3。 (1) 式中:Ti为作用在边界S上的牵引力;Xi为作用在体积V上的体力;vi为运动学容许的速度场;εij为与vi相容的应变速率场;σij为与Ti和Xi相对应的应力场。
对岸坡地震稳定性的分析基于以下条件开展:土体为均质、各向同性;平面应变条件;拟静力地震荷载;土体服从莫尔-库伦破坏准则及关联流动法则。
使用极限分析上限法计算临界水平地震系数,假设破坏面为穿过坡脚的对数螺旋面,破坏面方程为
r=r0e(θ−θ0)tanφ。 (2) 式中:r为对数螺旋线的半径;φ为土体内摩擦角;r0为对数螺旋线与坡顶面相交处的半径;θ0和θ为r0和r与水平线的夹角。
建立岸坡地震稳定性极限分析模型,如图 1所示。水位线与坡顶平齐,岸坡坡高为H,坡角为β。土体的总重度和有效重度分别为γ和γ′。破坏时ABC区域的土体作为刚体绕点O以角速度ω旋转,L为滑动土体顶面宽度,参数θ0和θh用于确定滑动面的具体位置,θf1至θfn分别代表码头从近海侧到陆地侧的第i排桩在对数螺旋破坏机制中所处的位置(i=1, 2, ⋯, n)。最靠近坡脚的一排桩与坡脚的距离为l,桩间距为d。
破坏模式的几何关系如下:
Hr0=[sinθhe(θh−θ0)tanφ−sinθ0], (3) Lr0=sin(β+θ0)−e(θh−θ0)tanφsin(θh+β)sinβ, (4) r0e(θh−θ0)tanφcosθh+l=r0e(θf1−θ0)tanφcosθf1,······r0e(θh−θ0)tanφcosθh+l+(n−1)d=r0e(θfn−θ0)tanφcosθfn。 (5) 速度不连续面AC上能量耗散速率为
D=cr02ω2tanφ(e2(θh−θ0)tanφ−1)。 (6) 极限状态下,作用在排桩上的侧向力可由Ito等[7]提出的公式确定。排桩对单位宽度岸坡的反作用力可由排桩受力除以桩间距得到,故可将高桩码头岸坡稳定性问题简化为平面应变问题。将排桩对岸坡稳定作用等效为梯形分布力,表达为
p=αm(a+bz), (7) 式中:a和b为使用Ito等[7]提出的公式计算出的参数;αm为群桩抗滑力系数;z为计算点距地表的深度。
近海侧到陆地侧排桩在滑动面以上产生抗滑力的有效深度z1至zn分别为:
z1=r0e(θf1−θ0)tanϕsinθf1−r0e(θh−θ0)tanϕsinθh+ltanβ ,……zn=r0e(θfn−θ0)tanϕsinθfn−r0e(θh−θ0)tanϕsinθh+[l+(n−1)d]tanβ。 (8) 近海侧到陆地侧排桩的抗滑力做功功率Wf1至Wfn分别为:
Wf1=−ωαm{(az1+bz122)[r0e(θh−θ0)tanϕsinθh−ltanβ]+(az122+bz133)} ,……Wfn=−ωαm{(azn+bz2n2)[r0e(θh−θ0)tanφsinθh−[l+(n−1)d]tanβ]+(az2n2+bz3n3)}。 (9) ABC区滑动土体重力做功功率可由OAC、OBC和OAB区土体重力的功率(W1-W2-W3)得到:
W1=γ′r30ω3(1+9tan2φ){(3tanφcosθh+sinθh)e[3(θh−θ0)tanφ]−3tanφcosθ0−sinθ0}, (10) W2=γ′r30ω6Lr0(2cosθ0−Lr0)sinθ0, (11) W3=γ′r30ω6e[(θh−θ0)tanφ][sin(θh−θ0)−Lr0sinθh]⋅{cosθ0−Lr0+cosθhe[(θh−θ0)tanφ]}。 (12) ABC区域土体单位重力加速度的水平惯性力做功功率可同理可表示为(W4-W5-W6):
W4=γr30ω3(1+9tan2φ){(3tanφsinθh−cosθh)⋅e[3(θh−θ0)tanϕ]−3tanϕsinθh+cosθ0} , (13) W5=γr30ω6Lr0sinθ0sinθ0, (14) W6=γr30ω6e[(θh−θ0)tanφ][sin(θh−θ0)−Lr0sinθh]⋅{sinθ0+sinθhe[(θh−θ0)tanφ]}。 (15) 将内部能量耗散速率与外力做功的功率相等,可得到岸坡的临界地震系数kc。
kc=D−Wf1−……−Wfn−(W1−W2−W3)(W4−W5−W6)。 (16) 前述分析表明,临界地震系数kc中有n+2个未知参数θ0, θh, θf1, …, θfn。将式(5)中的n个几何关系作为非线性约束条件,对kc进行优化,从而得到临界地震系数kc的最小上限解和破坏模式参数θ0和θh。
2. 高桩码头地震稳定性极限有限元分析
2.1 计算案例
地震作用下高桩码头岸坡断面如图 2所示,水位与坡顶地表平齐,岸坡坡高H = 10 m,坡度为1∶2。高桩码头在垂直海岸线的方向上有四排桩,近海排桩与坡脚的距离为2.5 m,排桩间距为5 m。在平行海岸线方向上,桩间距D1 = 5 m,桩间净距D2 = 3.8 m。岸坡土的黏聚力c = 5 kPa,内摩擦角φ = 30°重度γ = 20 kN/m3。群桩抗滑力系数αm取0.5。
使用本文理论方法求出临界地震系数kc= 0.2521,对应的破坏面形态参数θ0 = 1.3763、θh = 1.9883。
2.2 计算结果
为验证本文方法所得结果的合理性,使用有限元极限分析软件OPTUM G2对案例进行计算。OPTUM G2通过网格自适应迭代技术来获得更为精确的上、下限解,有限元极限分析网格如图 3所示。
通过极限分析有限元上、下限分析计算得到的临界地震系数kc分别为0.2493和0.2473。将上、下限解平均值作为精确解的估计,本文理论方法与有限元极限分析结果误差约为1.2%,这表明本文理论方法在计算岸坡临界地震系数方面的有效性。岸坡在拟静力地震荷载作用下的破坏模式可根据极限分析有限元上限分析得到的剪切耗散云图确定。理论方法得到的对数螺旋破坏面与有限元极限分析得到的剪切耗散分布如图 4所示,可看到对数螺旋破坏面与剪切耗散集中区域接近重合,证明了假设的破坏机制的合理性。
为进一步验证本文理论方法的合理性,分析了不同桩间净距与桩间距之比D2/D1、土体强度指标c和φ条件下理论分析方法与极限分析有限元方法得到的临界地震系数kc,如表 1所示。从表 1可知,理论分析方法得到的临界地震系数kc略大于极限分析上限有限元解,与本文方法上、下限解平均值误差均为1%左右,表明采用本文理论分析方法所得的岸坡临界地震系数合理。
表 1 临界地震系数对比Table 1. Comparison of critical seismic coefficientsD2/D1 c/kPa φ/(°) kc 有限元下限解 有限元上限解 本文理论方法 0.76 5 30 0.2473 0.2493 0.2521 0.80 5 30 0.2269 0.2287 0.2303 0.84 5 30 0.2107 0.2122 0.2131 0.76 2.5 30 0.1722 0.1747 0.1760 0.76 7.5 30 0.3040 0.3061 0.3101 0.76 5 27.5 0.2164 0.2181 0.2201 0.76 5 32.5 0.2810 0.2837 0.2870 3. 岸坡地震稳定性影响因素分析
强度折减法被广泛应用于边坡稳定性分析中,其通过安全系数FS折减土体强度指标c和φ,直至边坡达到临界状态,安全系数与土体强度参数的关系为
FS=cc'=tanφtanφ', (17) 式中:c和φ为土体的黏聚力和内摩擦角,c'和φ'为折减后的黏聚力和内摩擦角。
为分析土体强度参数c,φ和影响抗滑力的排桩布置参数D2/D1对稳定性和临界破坏面的影响,在岸坡拟静力稳定性分析中,使用临界地震系数kc和0.15g水平地震作用下岸坡安全系数FS评估其抗震性能。
3.1 排桩布置影响
排桩抗滑力受桩间净距与桩间距之比D2/D1影响。分别计算不同D2/D1条件下岸坡临界地震系数kc和在0.15g水平地震加速度条件下的安全系数,桩间净距与桩间距之比D2/D1分别为0.84,0.80,0.76,0.72,0.68。得到D2/D1对岸坡稳定性影响如图 5所示,临界地震系数kc和安全系数FS随D2/D1增大而减小。
不同D2/D1对应的临界破坏面如图 6所示,图 6表明随着D2/D1增大,滑动土体在上部范围缩小而在中下部范围扩大,破坏面曲率增大,但D2/D1对破坏机制的影响并不显著。
3.2 黏聚力的影响
计算不同黏聚力条件下岸坡临界地震系数kc和在0.15g水平地震加速度条件下的安全系数FS,其中,土体黏聚力c分别为1,3,5,7,9 kPa。黏聚力c的影响如图 7所示,图 7表明kc和FS随黏聚力c增大呈非线性增大,且随着黏聚力c增大,kc和FS增长速率降低。
不同黏聚力c对应的岸坡破坏面如图 8所示,从图 8可知,随着黏聚力c增大,临界破坏面以上滑动土体的范围明显增大,这表明土体黏聚力对岸坡的破坏机制有较显著的影响。
3.3 内摩擦角的影响
计算不同内摩擦角的岸坡临界地震系数kc和在0.15g水平地震加速度条件下的安全系数FS,内摩擦角φ分别为25°,27.5°,30°,32.5°,35°,内摩擦角φ对岸坡地震稳定性的影响如图 9所示,图 9表明kc和FS均随内摩擦角增大而增大。
不同内摩擦角对应的临界破坏面如图 10所示。从图 10可知随着内摩擦角增大,滑动土体在岸坡中下部范围减小,对数螺旋面曲率减小。内摩擦角对岸坡临界破坏面形态的影响在三个因素中最小。
4. 结论
针对地震条件下高桩码头岸坡稳定性提出了一种理论分析方法,分析了高桩码头岸坡的破坏模式和安全系数及随主要因素的变化规律。
(1)将排桩对码头岸坡的稳定作用简化为梯形分布力,采用对数螺旋破坏机制推导出岸坡的临界地震系数,并确定临界地震系数的最小上限解和对应的临界破坏面,建立了适用于高桩码头岸坡地震稳定性的理论分析方法。
(2)采用提出的理论分析方法分析了一个高桩码头案例的临界地震系数和破坏机制,并与极限分析有限元结果进行对比,所得结果与有限元极限分析结果相符,验证了理论分析方法的合理性。
(3)采用建立的理论分析方法,分析得出排桩布置、黏聚力、内摩擦角等对临界地震系数和破坏机制的影响规律,并通过强度折减法计算0.15g水平地震加速度条件下岸坡的安全系数。
(4)高桩码头岸坡的临界地震系数和安全系数随桩间净距与桩间距之比增大而减小,随着土体黏聚力和内摩擦角增大而增大。桩间净距与桩间距之比、土体黏聚力和内摩擦角均影响岸坡的临界破坏模式,但土体黏聚力对破坏面模式的影响更明显。
-
-
[1] LIU Q, DENG D P, YAO B J, et al. Analysis of the Karst springs' supply sources in rocky desertification area of Guanling-Huajiang, Guizhou, China[J]. Carbonates and Evaporites, 2020, 35(3): 1-11.
[2] 卢耀如, 张凤娥, 刘长礼, 等. 中国典型地区岩溶水资源及其生态水文特性[J]. 地球学报, 2006, 27(5): 393-402. doi: 10.3321/j.issn:1006-3021.2006.05.002 LU Yaoru, ZHANG Feng'e, LIU Changli, et al. Karst water resources in typical areas of China and their eco-hydrological characteristics[J]. Acta Geoscientica Sinica, 2006, 27(5): 393-402. (in Chinese) doi: 10.3321/j.issn:1006-3021.2006.05.002
[3] 魏兴萍, 谢德体, 倪九派, 等. 重庆岩溶槽谷区山坡土壤的漏失研究[J]. 应用基础与工程科学学报, 2015, 23(3): 462-473. doi: 10.16058/j.issn.1005-0930.2015.03.005 WEI Xingping, XIE Deti, NI Jiupa, et al. Soil erosion and loss on slope in Karst valley area, Chongqing with 137Cs[J]. Journal of Basic Science and Engineering, 2015, 23(3): 462-473. (in Chinese) doi: 10.16058/j.issn.1005-0930.2015.03.005
[4] 姚永慧. 中国西南喀斯特石漠化研究进展与展望[J]. 地理科学进展, 2014, 33(1): 76-84. https://www.cnki.com.cn/Article/CJFDTOTAL-DLKJ201401010.htm YAO Yonghui. Progress and prospect of Karst rocky desertification research in Southwest China[J]. Progress in Geography, 2014, 33(1): 76-84. (in Chinese) https://www.cnki.com.cn/Article/CJFDTOTAL-DLKJ201401010.htm
[5] WANG L, ZHENG F L, LIU G, et al. Seasonal changes of soil erosion and its spatial distribution on a long gentle hillslope in the Chinese Mollisol region[J]. International Soil and Water Conservation Research, 2021, 9(3): 394-404. doi: 10.1016/j.iswcr.2021.02.001
[6] 王世杰. 喀斯特石漠化概念演绎及其科学内涵的探讨[J]. 中国岩溶, 2002, 21(2): 101-105. https://www.cnki.com.cn/Article/CJFDTOTAL-ZGYR200202006.htm WANG Shijie. Concept deduction and its connotation of Karst rocky desertification[J]. Carsologica Sinica, 2002, 21(2): 101-105. (in Chinese) https://www.cnki.com.cn/Article/CJFDTOTAL-ZGYR200202006.htm
[7] 曾凌云, 汪美华, 李春梅. 基于RUSLE的贵州省红枫湖流域土壤侵蚀时空变化特征[J]. 水文地质工程地质, 2011, 38(2): 113-118. https://www.cnki.com.cn/Article/CJFDTOTAL-SWDG201102022.htm ZENG Lingyun, WANG Meihua, LI Chunmei. Study on soil erosion and its spatio-temporal change at Hongfeng Lake watershed based on RUSLE Model[J]. Hydrogeology & Engineering Geology, 2011, 38(2): 113-118. (in Chinese) https://www.cnki.com.cn/Article/CJFDTOTAL-SWDG201102022.htm
[8] 白占国, 万国江. 滇西和黔中表土中~7Be与137Cs分布特征对比研究[J]. 地理科学, 2002, 22(1): 43-48. https://www.cnki.com.cn/Article/CJFDTOTAL-DLKX200201008.htm BAI Zhanguo, WAN Guojiang. A comparative study on distribution of 7Be and 137Cs in the surface soils in the western Yunnan and the central Guizhou Provinces[J]. Scientia Geographica Sinica, 2002, 22(1): 43-48. (in Chinese) https://www.cnki.com.cn/Article/CJFDTOTAL-DLKX200201008.htm
[9] 姚邦杰, 刘琦, 任标, 等. 典型石漠化地区岩溶水系统循环演化分析[J]. 工程地质学报, 2019, 27(5): 1179-1187. https://www.cnki.com.cn/Article/CJFDTOTAL-GCDZ201905030.htm YAO Bangjie, LIU Qi, REN Biao, et al. Analysis of cyclic evolution of Karst water system in typical Karst rocky desertification area[J]. Journal of Engineering Geology, 2019, 27(5): 1179-1187. (in Chinese) https://www.cnki.com.cn/Article/CJFDTOTAL-GCDZ201905030.htm
[10] 周先齐, 徐卫亚, 钮新强, 等. 离散单元法研究进展及应用综述[J]. 岩土力学, 2007, 28(S1): 408-416. https://cpfd.cnki.com.cn/Article/CPFDTOTAL-AGLU200710007085.htm ZHOU Xianqi, XU Weiya, NIU Xinqiang, et al. A review of distinct element method researching progress and application[J]. Rock and Soil Mechanics, 2007, 28(S1): 408-416. (in Chinese) https://cpfd.cnki.com.cn/Article/CPFDTOTAL-AGLU200710007085.htm
[11] 李世海, 刘天苹, 刘晓宇. 论滑坡稳定性分析方法[J]. 岩石力学与工程学报, 2009, 28(增刊2): 3309-3324. https://www.cnki.com.cn/Article/CJFDTOTAL-YSLX2009S2008.htm LI Shihai, LIU Tianping, LIU Xiaoyu. Analysis method for landslide stability[J]. Chinese Journal of Rock Mechanics and Engineering, 2009, 28(S2): 3309-3324. (in Chinese) https://www.cnki.com.cn/Article/CJFDTOTAL-YSLX2009S2008.htm
[12] 闫钇全, 刘琦, 邓大鹏, 等. 表层岩溶裂隙带土壤地表流失/地下漏失室内模拟实验[J]. 中国岩溶, 2022, 41(2): 240-248. https://www.cnki.com.cn/Article/CJFDTOTAL-ZGYR202202008.htm YAN Yiquan, LIU Qi, DENG Dapeng, et al. Laboratory simulation study on soil surface loss and underground leakage in the epikarst fissure zone[J]. Carsologica Sinica, 2022, 41(2): 240-248. (in Chinese) https://www.cnki.com.cn/Article/CJFDTOTAL-ZGYR202202008.htm