Pull-out bearing behavior and failure mode of scaled suction caissons
-
摘要: 鳞片桶壁吸力基础是一种从仿生角度提出的新型吸力基础，能够有效提升基础的抗拔能力，可用于锚固系泊漂浮结构及海上光伏。通过倾斜单调加载模型试验，研究了不同系泊深度处鳞片吸力基础的抗拔承载力、破坏模式。鳞片吸力基础抗拔承载力随加载角度θ增加而减小，随系泊深度增加先增大后减小，当θ=20°，40°，60°，80°时，极限抗拔承载力较传统吸力基础增长8.4%，13.6%，41.2%和28.6%；当系泊深度位于0~0.58倍基础高度时，随着加载角度增加，鳞片吸力基础由前倾逐渐过渡为后倾转动破坏，而系泊深度为0.75~1倍基础高度时，基础发生后倾转动破坏，与传统吸力基础相比，鳞片吸力基础周围土体隆起与沉陷区域范围明显扩大，由此证明：鳞片吸力基础能够调动更多土体抵抗外荷载。基于Liu建立的理论计算模型，提出了最优系泊点处鳞片吸力基础抗拔承载力计算公式，并采用模型试验结果验证方法的准确性。Abstract: The scaled suction caisson is a new type of the offshore mooring foundation, which can be applied to the floating structures and offshore photostatic power. A series of the inclined loading model tests are carried out on the scaled suction caisson, and the pull-out bearing capacity and failure mode at different mooring depths are investigated. The study shows that the pull-out bearing capacity increases with the increase of the inclined loading angle, but firstly increases, then decreases with the increase of the mooring depth. Compared with that of the traditional suction caisson, the ultimate pull-out bearing capacity of the scaled suction caisson increases by 8.4%, 13.6%, 41.2% and 28.6% when the inclined angle θ=20°, 40°, 60° and 80°. When the mooring depth locates at 0~0.58 times the foundation height, the failure mode of the scaled suction caisson changes from forward to the backward rotation with the increasing inclined loading angle. The backward rotation occurs in the scaled suction caisson when the mooring depth locates at 0.75~1 times the foundation height. Under the disturbance of the scaled suction caisson, the uplift and subsidence appear in the soils surrounding the scaled suction caisson in the loading and opposite direction, respectively. The range of the uplift and subsidence areas obviously expands compared with that of the traditional suction caisson. The results demonstrate that more soils can be mobilized by the scaled suction caisson to resist the external loads. Based on the calculation model by proposed by Liu, the formula calculating the pull-out bearing capacity of the scaled suction caisson located at optimum mooring depth is proposed. The predicted pull-out bearing capacity using the proposed method agrees well with the model test results.
-
Keywords:
- scaled suction caisson /
- model test /
- pull-out bearing capacity /
- failure mode
-
0. 引言
水利工程关乎国计民生，其在蓄洪、发电和航运等领域发挥着重大作用。但其为国民经济发展带来效益的同时，也带来了一些地质灾害隐患，如库岸边坡失稳。库岸边坡失稳不仅会造成水库功能丧失，而且下泄的水流会给下游人民生命财产造成重大损失[1]。因此，开展库岸边坡稳定性研究具有十分重要的现实意义和巨大的经济与社会效益。库岸边坡稳定性受多种因素共同影响，其中库水位降落引起的土体负孔隙水压力和抗剪强度改变对非饱和库岸边坡稳定性影响较大[2]。此外，由于天然土体受不同的沉积、后沉积条件以及不同荷载历史等缘故，土体参数呈现出空间变异性，忽略该变异性将导致不准确的边坡可靠度评估结果[3-4]。因此，有必要进行库水降落条件下考虑多参数空间变异性的库岸边坡可靠度分析研究。
在边坡可靠度分析研究方面，学者们最初采用随机变量模型研究土体参数的不确定性，但随机变量模型不能有效地描述土体参数固有的空间变异性[4]。随机场理论为土体参数空间变异性的表征提供了一条有效的途径，为此，国内外学者结合随机场理论和可靠度理论在非饱和边坡可靠度分析方面做了大量有益的研究工作[5-9]：如Cho[5]通过Karhunen-Loève（K-L）展开方法离散土体渗透系数随机场，并通过蒙特卡洛模拟（MCS）方法研究了渗透系数空间变异性对库岸边坡渗流量和渗透坡降不确定性的影响；Tan等[6]采用中点法离散了非饱和土石坝坝坡的渗透系数和水力参数随机场，并结合有限差分法和MCS方法探讨了水力参数空间变异性和相关距离对土石坝坝坡渗流量和流速的影响；Zhang等[7]研究了抗剪强度参数空间变异性、地震系数和水位波动对边坡失效概率的影响；廖文旺等[8]通过K-L展开方法离散土体渗透系数随机场，并基于一阶可靠度分析方法研究了降雨强度和降雨持时对边坡可靠度的影响；王长虹等[9]采用试验变异函数模拟土体非饱和渗透系数随机场，结合MCS方法开展了考虑渗流作用的非饱和土边坡稳定性可靠度分析。然而，以上学者的研究工作还存在以下几点不足：①在分析非饱和边坡可靠度时仅考虑单一渗透系数或水力参数的空间变异性对边坡渗流和稳定的作用，忽略了抗剪强度参数空间变异性的影响；②既有非饱和边坡可靠度分析通常只关注某时刻的静态可靠度问题，对边坡的时变可靠度问题缺乏深入研究，亟需发展有效的非饱和边坡时变可靠度分析理论和方法；③采用传统MCS方法计算考虑多参数空间变异性的非饱和边坡失效概率时计算量十分巨大，针对复杂边坡及小失效概率问题，这一缺陷更为明显，因此大多数研究仅分析了土体参数空间变异性对边坡渗流量、流速、渗透坡降及安全系数等统计特征参数的影响，而参数空间变异性对非饱和边坡失效概率影响的研究仍不充分。
为此，本文提出了库水位降落下考虑多参数空间变异性的非饱和库岸边坡时变可靠度分析方法。首先建立边坡确定性分析模型，采用随机场模型对土体参数的空间变异性进行模拟。其次，为避免由于多个关键参数经随机场离散后的随机变量维度过高导致代理模型陷入“维度灾难”问题，采用分片逆回归（SIR）方法先进行降维预处理，该方法理论严密且易于实现，在岩土工程研究中应用广泛[10]。接着，利用高斯过程回归（GPR）机器学习算法构建代理模型，GPR方法具有严格的统计理论基础、灵活的非参数推断、参数自适应获取等优点，对于小样本及非线性复杂问题具有良好的适应性[11]。进而，基于所构建的SIR-GPR模型通过MCS评估库岸边坡失效概率。在此，以水位降落作用下的非饱和库岸边坡可靠度问题为例，探讨库岸边坡在不同水位降落工况下的安全系数和失效概率变化规律，验证所提SIR-GPR-MCS方法的有效性，并基于所提方法对土体参数的敏感性进行分析。
1. 研究方法
1.1 随机场模型表征土体空间变异性
边坡岩土体参数通常表现出显著的空间变异性，随机场模型可有效表征该变异性[12-13]。在此，采用K-L展开方法离散二维随机场，该方法使用一系列随机变量和确定性空间函数对其进行表示[14]。本文选用高斯型相关函数，具体如下：
ρ（x,y）=exp(−[(|x−x'|lx)2+(|y−y'|ly)2])。 (1) 式中：(x, y)和(x', y')为随机场空间中任意两点的坐标；lx和ly分别为水平和垂直相关距离。第i个随机场Hi经过K-L展开离散为
Hi(x,y;ω)=μi+∞∑j=1σi√λifj(x,y)ξi,j(ω)，≅μi+M∑j=1σi√λifj(x,y)ξi,j(ω)。 (2) 式中：μi和σi分别为随机场的均值和标准差；ω为随机事件；ξi,j(ω)为一组独立标准正态随机变量；λi和fj(x, y)分别为协方差函数的特征值和特征函数；M为展开项数，其值取决于期望的精度、相关距离和随机场的维数等，一般通过期望能比率因子ε大于95%进行确定[15]，具体表达式为
ε=M∑i=1λi/M∑i=1λi∞∑i=1λi∞∑i=1λi。 (3) 土体抗剪强度参数之间（如有效黏聚力和有效内摩擦角之间）以及水力参数之间（如土-水特征曲线的其中两个拟合参数之间）通常呈现负相关性，且土体参数值均为正，因此采用相关对数正态随机场对参数的空间变异性进行表征[16]。式（2）K-L展开公式的相关对数正态形式为
Hi(x,y;ω)=exp[μln+M∑j=1σln√λifj(x,y)ξi,j(ω)]=exp[μln+σlnG(x,y;ω)]。 (4) 式中：μln和σln分别为相关对数正态分布随机场的均值和标准差，计算表达式为σln=[ln(1+(σi/μi)2]1/2，μln= lnµi-0.5σ2ln，G(x, y; ω)为展开项数为M的标准正态分布随机场。
1.2 渗流分析
通过K-L展开方法离散了边坡土体参数随机场后，进而对边坡非饱和渗流进行有限元分析。时变水位作用下的非饱和边坡的稳定性问题是一个典型的流固耦合过程。由达西定律和质量守恒原理推导出的饱和-非饱和渗流控制微分方程为[17]
∂∂x[kx(θ)∂H∂x]+∂∂y[ky(θ)∂H∂y]=γw∂θ∂ψ∂H∂t。 (5) 式中：θ为土体的体积含水量；kx (θ)，ky (θ)分别为土体在水平和垂直方向的渗透系数；H为总水头；γw为水的重度；t为时间，对于稳定渗流，∂H/∂t=0；ψ为土体的基质吸力；∂θ/∂ψ为土体体积含水量与基质吸力水头曲线的斜率，可由土-水特征曲线得到。土-水特征曲线的滞水效应等会影响非饱和土体抗剪强度特性及其力学行为，根据《非饱和土试验方法标准：T/CECS 1337—2023》开展相关土-水特征曲线测定试验、非饱和土三轴试验、水平土柱渗水试验及渗气试验等有利于帮助工程师准确地评估边坡土质条件[18]。为拟合或预测土-水特征曲线，目前常采用的模型有Brooks-Core模型、Van Genuchten模型和Fredlund-Xing模型。鉴于Van Genuchten模型参数化形式，对饱和、非饱和条件的全面描述及物理可解释性等优点[2]，本文采用该模型对土-水特征曲线进行拟合，其表达式为
θ=θr+θs−θr[1+(aψ)n]m。 (6) 式中：θs为土体的饱和体积含水量；θr为残余体积含水量；a，n和m均为土-水特征曲线模型拟合参数，m=1-1/n，n > 1。
渗透系数曲线表达式为[19]
ki(θ)=ks(θ−θrθs−θr)0.5[1−(1−(θ−θrθs−θr)1m)m]2。 (7) 式中：ks为土体的饱和渗透系数；i=x, y；ks，a，n统称为水力参数。由式（5）~（7）可以进行饱和-非饱和渗流分析，非饱和土体渗透系数的空间变异性可以直接通过模拟土体饱和渗透系数ks和水力模型参数a，n的空间变异性实现。
1.3 稳定性分析
渗流分析得到的孔隙水压力分布结果用于后续的边坡稳定性分析。Fredlund等[20]采用净法向应力和基质吸力双应力状态变量，结合Mohr-Coulomb强度理论，提出的非饱和土抗剪强度公式为
τ=c′+(σn−ua)tanφ′+(ua−uw)tanφb。 (8) 式中：τ为抗剪强度；c'和φ'分别为有效黏聚力和有效内摩擦角；φb为吸力内摩擦角，表示抗剪强度随基质吸力变化的速率，其计算表达式为tanφb=tanφ′[(θ− θr)/(θs−θr)]；σn为破坏面上的净法向应力；ua和uw分别为破坏面上的孔隙气压力和孔隙水压力。假设条间作用力为半正弦函数，基于Morgenstern-Price方法进行边坡稳定性分析。
1.4 代理模型
采用边坡稳定安全系数代理模型代替复杂的边坡有限元渗流及稳定性分析可以减轻计算量。非饱和边坡土体的多个关键参数随机场经过K-L展开方法离散后，通常仍然会面临维度过高的问题（特别是针对大区域边坡）。此外，边坡的极限状态函数具有高度非线性。传统的代理模型无法保证其精度，而直接进行MCS评估边坡失效概率将涉及大量复杂有限元分析，计算非常耗时[21]。为此，本节提出一种基于SIR-GPR的代理模型近似替代边坡的真实功能函数。
（1）分片逆回归降维预处理
随着回归问题维度增加，计算成本将呈指数级增长，造成统计分析难度急剧上升，该现象在统计学领域被称为“维度灾难”[10]。为避免高维度输入随机变量引起的“维度灾难”问题，本节采用分片逆回归方法进行降维预处理。其基本原理是将高维的随机变量向有效降维空间投影，得到一组相互正交的基向量[22]。有效降维空间中的信息更加简洁有效，使得回归分析过程变得简单。在SIR分析中，假设原始数据模型为
Y=f(ξ)=f(ξ1,ξ2,⋯,ξr)。 (9) 式中：ξ=[ξ1, ξ2, ..., ξr]T为r维标准正态输入随机变量向量，Y为输出向量。经SIR降维后该数据模型变为
Y=g(X1,X2,⋯,Xs)=g(ˆβT1ξ,ˆβT2ξ,⋯,ˆβTsξ,ε)。 (10) 式中：s为降维维度；ˆβi(1≤i≤s)为投影向量，这些向量构成有效降维空间；Xi(1≤i≤s)为降维后的映射变量；ε为随机扰动。式（10）成立时，原始r维空间经投影向量ˆβi投影至s(s＜r)维有效降维空间上，用预测变量的低维线性组合替换原始变量不损失任何有效信息[23]。SIR分析的主要目标就是求出该投影向量ˆβi，形成有效降维空间，完成以上降维预处理过程。SIR方法的工作步骤如下。
a）将ξi标准化：
ξi=∑−0.5ξiξi（ξi−¯ξi）， (11) 式中，∑ξiξi为样本的协方差矩阵；¯ξi为样本均值。ξi经过标准化后，均值为0，方差为1。并计算ξi对应的响应fi。
b）将fi划分为D个分片区间Md(1⩽。
c）计算各分片输入变量的均值 \overline {{\xi _d}} (1 \leqslant d \leqslant D) 。
\overline {{\xi _d}} = \frac{1}{{{n_d}}}\sum\limits_{{f_i} \in {M_d}} {\overline {{\xi _i}} } \text{，} (12) 式中，nd为第d个分组的样本数目。
d）计算条件期望Ε(ξi|fi)。
e）考虑到分片样本数不同，计算加权协方差矩阵 \hat{\boldsymbol{C}} ：
\hat{\boldsymbol{C}} = \sum\limits_d^D {\frac{{{n_d}}}{N}\overline {{\xi _d}} \overline {{\xi _d}} _{}^{\text{T}}} 。 (13) 式中，N为样本总数。
f）对分片输入变量均值 \overline {{\xi _d}} 进行主成分分析，确定 \hat{\boldsymbol{C}} 对应的特征值 {\hat \lambda _i} 及特征向量 {\hat \alpha _i}\; ：
\hat{\boldsymbol{C}}{\hat \alpha _i} = {\hat \lambda _i}{\hat \alpha _i}(1 \leqslant i \leqslant r) 。 (14) g）将特征值 {\hat \lambda _i} 排序，前s个较大特征值对应的特征向量即为有效降维空间的投影向量：
{\hat \beta _i} = {\hat \alpha _i}\Sigma _{{\xi _i}{\xi _i}}^{ - 0.5}(1 \leqslant i \leqslant s) 。 (15) h）计算分片在有效降维空间中的投影值：
{\boldsymbol{X}} = \hat \beta _i^{\text{T}} \boldsymbol{\xi} (1 \leqslant i \leqslant s)。 (16) （2）基于高斯过程回归的蒙特卡洛模拟方法
根据降维后的随机变量矩阵和安全系数矩阵建立GPR代理模型，近似替代边坡稳定性分析的复杂功能函数。N组训练样本组成的训练样本集D可表示为
D = \left\{ {({X_e},{Y_e})} \right\} = \left\{ {{\boldsymbol{X}},{\boldsymbol{Y}}} \right\}{\text{ (}}e = 1,2, \cdots ,N) 。 (17) 式中：N为训练样本数量；X为降维后的所有样本输入向量；Y为降维后所有样本输出响应向量。假设样本集存在噪声且服从标准高斯分布，据此建立回归问题的一般回归模型为
\boldsymbol{Y}=g(\boldsymbol{X})+\varepsilon, \varepsilon \sim N\left(0, \sigma_t^2 \boldsymbol{I}_n\right) 。 (18) 式中：Y=[Y1, Y2, …, YN]T为受噪声污染的N×1维观测向量；X=[X1, s, X2, s, …, XN, s]T为N×s维输入随机变量，s为随机变量的维度； \varepsilon 为相互独立且服从高斯分布的观测噪声向量； {\sigma _t}^2 为噪声方差；In为单位矩阵。高斯过程回归函数g(X)可定义为
g({\boldsymbol{X}})\sim GP[m({\boldsymbol{X}}),{\boldsymbol{K}}({\boldsymbol{X}},{\boldsymbol{X}})] \text{，} (19) m({\boldsymbol{X}}) = E\left[ {g({\boldsymbol{X}})} \right] \text{，} (20) {\boldsymbol{K}}({\boldsymbol{X}},{\boldsymbol{X}}) = E\left\{ {[g({\boldsymbol{X}}) - m({\boldsymbol{X}})]{{[g({\boldsymbol{X}}) - m({\boldsymbol{X}})]}^{\text{T}}}} \right\} 。 (21) 式中：m(X)为输入样本X的均值函数，通常为0；K(X, X)为输入样本X自身的协方差函数。由式（18），（19）可得观测目标值Y的先验分布为
\boldsymbol{Y} \sim N\left(0, \boldsymbol{K}(\boldsymbol{X}, \boldsymbol{X})+\sigma_t^2 \boldsymbol{I}_n\right) 。 (22) 进一步可得训练集输出值Y（也称观测值）和测试集输出值Y∗（也称预测值）的联合高斯先验分布为
\left[\begin{array}{l} \boldsymbol{Y} \\ \boldsymbol{Y}^* \end{array}\right] \sim N\left\{0,\left[\begin{array}{ll} \boldsymbol{K}(\boldsymbol{X}, \boldsymbol{X})+\sigma_t^2 \boldsymbol{I}_n & \boldsymbol{K}\left(\boldsymbol{X}, \boldsymbol{X}^*\right) \\ \boldsymbol{K}\left(\boldsymbol{X}, \boldsymbol{X}^*\right)^{\mathrm{T}} & \boldsymbol{K}\left(\boldsymbol{X}^*, \boldsymbol{X}^*\right) \end{array}\right]\right\} 。 (23) 式中：K(X, X*)是训练集输入值X与测试样本输入值X*之间的N×1阶协方差矩阵，K(X*, X*)是测试样本X*自身的协方差矩阵。本文选用的协方差函数型式为常用的平方指数型，其表达式为
{\boldsymbol{K}}({{\boldsymbol{X}}_p},{{\boldsymbol{X}}_q}) = {\sigma _t}^2{\delta _{pq}} + \sigma _{\text{f}}^{\text{2}}\exp \left( { - \frac{1}{{2{l^2}}}{{({{\boldsymbol{X}}_p} - {{\boldsymbol{X}}_q})}^2}} \right) 。 (24) 式中：Xp与Xq是一个样本输入值在同一维度上的两个任意元素。L， \sigma _{\text{f}}^{\text{2}} ， {\sigma _t}^2 为平方指数型协方差函数的3个超参数，一般通过极大似然法搜索出最优解[24]。对数似然函数形式为
\begin{gathered} L = \ln p({\boldsymbol{Y}}|{\boldsymbol{X}}) = - \frac{1}{2}{{\boldsymbol{Y}}^{\text{T}}}{({\boldsymbol{K}} + {\sigma _t}^2{{\boldsymbol{I}}_n})^{ - 1}}{\boldsymbol{Y}} - \\ \frac{1}{2}\ln |{\boldsymbol{K}} + {\sigma _t}^2{{\boldsymbol{I}}_n}| - \frac{N}{2}\ln 2{\mathtt{π}}。 \end{gathered} (25) 获取最优超参数后根据贝叶斯方法通过新的观测值不断更新预测值的概率分布，在新的X∗，X，Y下推断出Y∗的后验分布为
p({{\boldsymbol{Y}}^*}|{{\boldsymbol{X}}^*},{\boldsymbol{X}},{\boldsymbol{Y}})\sim N[m({{\boldsymbol{Y}}^*}),{\boldsymbol{K}}({{\boldsymbol{Y}}^*},{{\boldsymbol{Y}}^*})]。 (26) 式中，Y∗的均值m(Y∗)和方差K(Y*, Y*)的具体表达式为：
m({{\boldsymbol{Y}}^*}){\text{ = }}{\boldsymbol{K}}{({\boldsymbol{X}},{{\boldsymbol{X}}^*})^{\text{T}}}{\left[ {{\boldsymbol{K}}({\boldsymbol{X}},{\boldsymbol{X}}) + {\sigma _t}^2{{\boldsymbol{I}}_n}} \right]^{ - 1}}{\boldsymbol{Y}} \text{，} (27) \begin{aligned} & {\boldsymbol{K}}({{\boldsymbol{Y}}^*},{{\boldsymbol{Y}}^*}) = {\boldsymbol{K}}({{\boldsymbol{X}}^*},{{\boldsymbol{X}}^*}) - {\boldsymbol{K}}{({\boldsymbol{X}},{{\boldsymbol{X}}^*})^{\text{T}}} \\ & \quad\quad\quad\quad {\left[ {{\boldsymbol{K}}({\boldsymbol{X}},{\boldsymbol{X}}) + {\sigma _t}^2{{\boldsymbol{I}}_n}} \right]^{ - 1}}{\boldsymbol{K}}({\boldsymbol{X}},{{\boldsymbol{X}}^*}) 。 \end{aligned} (28) 采用以上GPR模型不断逼近边坡稳定性分析的隐式功能函数，预测边坡 {F_{\text{S}}} 。在此基础上，结合MCS方法进行时变水位作用下的边坡稳定可靠度分析，得到第t时刻边坡的失效概率为
{P_{\text{f}}}(t) = \frac{1}{{{N_{{\text{sim}}}}}}\sum\limits_{i = 1}^{{N_{{\text{sim}}}}} {I\left\{ {{F_{\text{S}}}({X_i},t) < 1} \right\}} 。 (29) 式中：Nsim为MCS抽样次数； {F_{\text{S}}} (Xi, t)为边坡在第t时刻的安全系数，Xi=(c', φ', ks, a, n)为第i组土体参数随机场实现值；Pf (t)为第t时刻边坡的失效概率； I {⋅}为指示函数， {F_{\text{S}}} (Xi, t) < 1时 I {⋅}=1，否则为0。
（3）模型精度评估
采用以上代理模型替代原有计算耗时较长的有限元模型，很大程度上节约了计算成本。但代理模型仅仅是一种近似模型，其与真实有限元仿真模拟之间仍存在一定偏差[2, 10]。通常使用精度评价指标来评估代理模型的响应预测能力：
{R^2} = 1 - \frac{{\frac{1}{N}\sum\limits_{i = 1}^N {{{({Y_i} - \bar Y)}^2}} }}{{\frac{1}{N}\sum\limits_{i = 1}^N {{{({y_i} - \bar y)}^2}} }}\text{，} (30) \sqrt {{\text{MSE}}} = \sqrt {\frac{{\sum\nolimits_{i = 1}^N {{{({Y_i} - {y_i})}^2}} }}{N}} 。 (31) 式中：R2为决定系数；\sqrt {{\text{MSE}}} 为均方根误差；{y_i}为第i个有限元模型仿真计算的 {F_{\text{S}}} ； {Y_i} 为代理模型预测的第i个边坡的 {F_{\text{S}}} ；回归模型的决定系数越接近于1或者均方根误差越接近于0，表明代理模型的精度越高。
（4）执行步骤
图 1给出了采用SIR-GPR-MCS方法进行考虑多参数空间变异性的非饱和边坡可靠度分析的流程，主要步骤如下。
![]()
图 1 考虑参数空间变异性的非饱和边坡可靠度分析流程图Figure 1. Flow chart of reliability analysis of unsaturated spatially variable reservoir bank slopesa）建立边坡确定性模型。用GEO-STUDIO软件建立非饱和边坡模型，包括定义材料属性、边界条件、土-水特征曲线和渗透系数曲线等。流固耦合分析边坡的渗流和稳定性状态，将模型文件另存为拓展名为xml的接口文件，具体见1.2和1.3节。
b）建立新的输入计算模型。首先确定参数ks，a，n，c'，φ'的统计特征，包括均值、变异系数、分布类型、相关函数、互相关系数和相关距离等。采用K-L展开方法离散参数随机场，之后通过MATLAB程序改写原xml文件中的参数，得到N个新的输入文件，具体见1.1节。
c）计算并提取结果。采用WINDOWS BAT批处理脚本调用GEO-STUDIO软件内核进行N个输入文件的计算，并输出储存着渗流分析结果和可能出现的滑动面信息的结果文件，再通过MATLAB程序提取N个最危险滑动面对应的N个 {F_{\text{S}}} 。
d）回归分析预测。采用SIR方法降维预处理随机样本矩阵，根据N个输入样本矩阵和N个输出 {F_{\text{S}}} 矩阵构建SIR-GPR代理模型。评估模型精度，满足精度要求的代理模型用于后续边坡 {F_{\text{S}}} 的预测，具体见1.4节。
e）评估边坡失效概率。采用MCS方法评估边坡Pf，具体见1.4节。
2. 算例
2.1 边坡模型与确定性分析
在此选取文献[2]中一坡高10 m，底厚10 m，坡比为1∶2的均质库岸边坡（见图 2）验证所提方法的高效性及实用性。土体饱和重度为20 kN/m3，进行该黏性土坡确定性分析所需的其他参数见表 1[2]。采用GEO-STUDIO数值计算软件进行分析，共设置稳态渗流、瞬态渗流、稳定性3个状态。初始稳态渗流分析时，库水位位于19 m处。瞬态渗流分析时边坡坡面边界条件为总水头随时间线性变化的函数，其余边界为不透水边界。总水头从初始状态19 m降落到10 m处，设置了水位降落速率分别为0.5，1，1.8，3 m/d的4种计算工况（见图 3）[2]。非饱和土体土-水特征曲线（Van Genuchten模型）所用参数{\theta _s}和{\theta _r}分别取值50%和3.5%[2]，土体对应的土-水特征曲线和渗透系数曲线如图 4（a），（b）所示。
![]()
图 2 库岸边坡尺寸和确定性模型Figure 2. Dimensions and deterministic model for reservoir bank slope表 1 土体物理力学参数Table 1. Physical and mechanical parameters of soil确定性分析 渗流 稳定 参数 ks/(10-6 m·s-1) a/kPa n c'/kPa φ'/(°) 均值 1.160 28.44 1.720 5 30 标准差 0.696 11.376 0.344 1.5 6 变异系数 0.6 0.4 0.2 0.3 0.2 互相关系数 -0.25 -0.5 ![]()
图 3 不同计算工况下库水位变化Figure 3. Variation of reservoir water level under different working conditions![]()
图 4 土-水特征曲线和渗透系数曲线Figure 4. Soil-water characteristic curve and permeability coefficient curve以水位降落速率为1.8 m/d的计算工况为例，边坡孔隙水压力分布情况如图 5（a）~（d）所示，边坡内部孔隙水压力逐渐减小，浸润线逐渐降低。随后进行稳定性分析，设置10个时间步区间（见表 2）以便观察边坡 {F_{\text{S}}} 随时间变化规律。随着水位下降， {F_{\text{S}}} 值由1.846逐渐下降至1.256；5 d后库水位达到稳定状态，此后 {F_{\text{S}}} 略微上升至1.479。水位降落阶段库岸边坡浸润线下降存在滞后效应，在边坡内部形成较大的水力梯度，且作用在坡面上的静水压力消失，导致边坡内部超孔隙水压力和渗透力指向坡外，在边坡表面形成反向渗流场，产生向外拖拽力，边坡下滑力增大， {F_{\text{S}}} 不断减小[25]。水位降落越快，库岸边坡内部渗透力越大，常在水位降落终止时 {F_{\text{S}}} 降至最低[25]。在水位稳定阶段，边坡内部浸润线由于存在滞后效应而继续下降，边坡内部应力重分布，边坡 {F_{\text{S}}} 小幅度上升。
![]()
图 5 孔隙水压力分布和最危险滑动面Figure 5. Distribution of pore water pressure and most dangerous sliding surface表 2 不同时刻库岸边坡安全系数Table 2. Safety factor of reservoir bank slope at different moments时间步序号 时间/d FS 1 0 1.846 2 1 1.544 3 2 1.358 4 3 1.258 5 5 1.256 6 8 1.316 7 11 1.355 8 15 1.393 9 22 1.441 10 30 1.479 2.2 可靠度分析
根据土体参数统计特征，利用K-L级数展开方法对其进行随机场离散。参数a与n，c'与φ'之间的互相关系数分别取为ρa, n=-0.25，ρc', φ'=-0.5[2]。在此，主要考虑土体参数ks，a，n，c'，φ'的空间变异性。水平向相关距离范围一般为2~30 m，垂直向相关距离范围一般为0.1~5 m，本文假定水平向相关距离lx=30 m，垂直向相关距离ly=3 m[26]。选用高斯型自相关函数，K-L级数展开中截断项数取10时随机场的离散误差能满足精度要求。由此，ks，a，n，c'，φ'随机场共计离散得到50个随机变量。考虑该库岸边坡的实际尺寸将模型划分为1230个三角形和四边形有限元单元。图 6（a）~（e）给出了土体参数随机场的一次典型实现，图中c'与φ'，a与n之间随机场的负相关关系也显而易见。
先采用拉丁超立方抽样（LHS）方法[26]对该边坡进行可靠度分析，作为本文所提SIR-GPR-MCS方法的对照。为保证LHS方法的计算精度，对库岸边坡在四种时变水位工况下抽取2000个样本，之后用MATLAB调用GEO-STUDIO内核对这8000个有限元模型进行FS计算，再评估边坡的Pf，得到FS均值和Pf随时间的变化规律，如图 7（a），（b）所示。水位下降初期边坡FS迅速降低，对应的Pf迅速升高，并且水位降落速率越快时FS和Pf的变化速率越快；之后FS和Pf的变化趋势变缓；当水位分别在第18天、第9天、第5天、第3天降到稳定水位时，库岸边坡FS也分别在对应的时间稳定在1.0附近，Pf稳定在0.52附近（对于本文案例，水位处于稳定阶段时，边坡FS上升和Pf下降的变化幅度较小，可忽略不计）。
![]()
图 7 不同水位降落工况下库岸边坡安全系数和失效概率随时间变化Figure 7. Safety factors and probability failure of reservoir bank slope with time for different water level plunge conditions接着，采用SIR-GPR-MCS方法进行可靠度分析。在各种时变水位工况下，本文各随机选取200组训练样本数据用于构建SIR-GPR模型。为确定满足模型精度要求的最佳降维维度，计算代理模型在不同维度r=1, 2, …, 50下的均方根误差。结果如图 8（a）所示，r=1时均方根误差为0.0237，已满足精度要求，因此本文选定SIR维度为1。接着进行1×106次MCS计算边坡Pf，结果如图 7（c）所示。图 7（b），7（c）的计算结果高度一致，表明本文所提的SIR-GPR-MCS方法用于分析考虑空间变异性的非饱和边坡时变可靠度问题较为准确。
从代理模型的精度来看，进一步评估SIR-GPR代理模型在水位降落过程中各时刻的回归分析能力。以水位降落速率为1.8 m/d的计算工况为例，SIR-GPR模型误差分析结果如图 8（b），（c）所示。水位降落前5天SIR-GPR模型的均方根误差从0.04开始持续下降，而后有略微上升并稳定在0.023附近；水位降落前5天SIR-GPR模型的决定系数在0.94~0.96波动，而后稳定在0.96附近。此现象出现的原因可能是前5天水位降落阶段边坡的FS较大，对应的Pf较小，代理模型在计算这类失效概率较小的问题时可以通过适当增加训练样本数量来进一步提升计算精度。总的来说，本文案例中两个模型误差评价指标均在合理范围内，表明SIR-GPR模型用于边坡FS的预测是较为精准的，考虑到在水位降落阶段重新确定训练样本数量构建代理模型的计算成本较大，因此本文不再展开研究。此外，随机选取200组样本分别使用有限元计算和代理模型计算FS，计算结果对比见图 9，两种方法计算的FS拟合于45°线也表明代理模型的预测性能较好。
从代理模型的计算效率来看，在配置内存为16 GB、中央处理器为11th Gen Intel (R) Core(TM) i7-11800H和主频为2.3 GHz的计算机上，采用本文所提方法计算边坡在某一特定工况下的可靠度总耗时4.5 h。其中包括前处理生成200个随机场样本耗时3.1 h，计算200个边坡样本的渗流和稳定性耗时1.2 h，后处理生成SIR-GPR模型和1×106次MCS计算耗时0.2 h。而LHS方法计算需要生成2000个随机场样本进行有限元分析再计算Pf，总耗时约32 h。因此，计算边坡在4种水位降落工况下的Pf，采用本文所提方法总耗时18 h，而LHS方法总耗时长达128 h。显然，本文所提基于SIR-GPR-MCS的非饱和库岸边坡可靠度分析方法能显著提高计算效率。
2.3 参数敏感性分析
本节将SIR-GPR-MCS方法用于具有空间变异性的边坡参数敏感性分析之中，探讨不同土体参数统计量对边坡Pf的影响规律。
（1）参数空间变异性的影响
本文土体参数变异系数（COV）取值参考文献[2]，在此，分别讨论COVks在0.3~0.9，COVa在0.2~0.8，COVn在0.1~0.4，COVc'在0.2~0.5，COVφ'在0.05~0.2之间变化对Pf的影响。当对某一参数的空间变异系数进行敏感性分析时，其余参数变异系数取值见表 1。
以水位降落速率1.8 m/d计算工况为例，边坡Pf随参数变异系数的变化规律见图 10（a）~（e）。结果为抗剪强度参数变异系数COVc'，COVφ'在较小范围内增加时，边坡Pf明显增大，表明抗剪强度参数的空间变异性显著影响边坡可靠度计算结果。而Pf随水力参数变异系数COVa，COVn， {\text{CO}}{{\text{V}}_{{k_{\text{s}}}}} 的影响相对较小。
![]()
图 10 土体参数空间变异性对库岸边坡Pf的影响Figure 10. Influences of spatial variability of soil parameters on Pf另外，图 11给出了第30天时考虑不同参数空间变异性的边坡FS概率密度函数（PDF）。如果忽略土体多参数的空间变异性（水平和垂直相关距离均设成1000 m），将得到更低更宽的PDF。PDF峰度及偏度减小，往FS小值区域偏移。可见忽略土体多参数的空间变异性可能高估边坡Pf，因此，在可靠度分析中有必要考虑多个土体参数的空间变异性。
![]()
图 11 土体参数空间变异性对FS概率密度函数的影响Figure 11. Influences of spatial variability of soil parameters on probability density function of FS（2）参数互相关性的影响
本节在参数变异系数不变的情况下，研究相关系数ρa, n和ρc', φ'对边坡Pf的影响。ρa, n和ρc', φ'的变化范围均为-0.5~0.5，当对某一互相关系数进行参数敏感性分析时，另一互相关系数的取值见表 1。
以水位降落速率为1.8 m/d的计算工况为例，边坡Pf随ρa, n、ρc', φ'的变化规律分别见图 12（a），（b）。结果表明边坡Pf与ρc', φ'正相关，c'与φ'之间的相关性显著影响边坡可靠度分析计算结果。
![]()
图 12 土体参数相关性对库岸边坡Pf的影响Figure 12. Influences of correlation among soil parameters on reservoir bank slope Pf图 10，12结果均表明，基于不同的土体参数统计量（包括变异系数及相关系数等）得到的边坡失效概率结果截然不同，因此准确地获取边坡土体参数统计量对于边坡可靠度分析及设计来说至关重要。
3. 结论
本文主要综合非饱和土理论、随机场理论、分片逆回归降维预处理、高斯过程回归分析和蒙特卡洛模拟方法，提出了同时考虑时效特性和土体多参数空间变异性的非饱和边坡SIR-GPR-MCS可靠度分析方法。探讨了库岸边坡失效概率在不同水位降落工况下的变化情况，并运用SIR-GPR-MCS方法对关键土体参数的空间变异系数和相关系数进行了敏感性分析，得到以下3点结论。
（1）所提SIR-GPR-MCS方法在求解非饱和边坡失效概率时计算精度较高，且极大地降低了计算量，为复杂边坡稳定可靠度问题提供了一种高效的评估途径。
（2）水位降落作用下库岸边坡的FS逐渐下降，相应的Pf逐渐上升，且在库水降落初期时变化迅速，而后变化趋势减缓直至Pf趋于稳定。并且对于不同水位降落工况，水位降落越快时Pf前期变化速率越快。
（3）参数敏感性分析结果表明多个土体参数的空间变异性和相关性均对库岸边坡可靠度计算结果有影响。为准确地进行库岸边坡可靠度分析，有必要同时考虑多个土体参数的空间变异性。
-
![]()
图 5 鳞片吸力基础抗拔荷载-位移曲线
Figure 5. Pull-out load-deflection curves of scaled suction caisson
![]()
图 6 鳞片吸力基础极限抗拔承载力
Figure 6. Ultimate pull-out bearing capacities of scaled suction caisson
![]()
图 8 鳞片吸力基础转动角度-位移曲线
Figure 8. Rotation angle-deflection curves of scaled suction caisson
![]()
图 9 不同倾斜荷载下鳞片吸力基础破坏模式
Figure 9. Failure modes of scaled suction caisson under different loading condition
表 1 吸力基础模型尺寸
Table 1 Dimensions of caisson models
基础形式 直径D/mm 高度H/mm 壁厚t1/
mm顶板t2/
mm鳞片高度H1/mm 鳞片宽度H1/mm 高径比H1/D1 传统吸力基础 80 160 1.5 8 — — 鳞片吸力基础 2 0.2 10 
下载: 导出CSV
-
[1] ZHANG D G, CHEN Y J, ZHANG T Y. Floating production platforms and their applications in the development of oil and gas fields in the South China Sea[J]. Journal of Marine Science and Application, 2014, 13(1): 67-75. doi: 10.1007/s11804-014-1233-2
[2] MARTINEZ A, FROST J D. The influence of surface roughness form on the strength of sand–structure interfaces[J]. Géotechnique, 2017, 7(1): 104-111. doi: 10.1680/jgele.16.00169
[3] GRAY J, LIESMANN H W. The kinetics of locomotion of the grass-snake[J]. Journal of Experimental Biology, 1950, 26(4): 354-367. doi: 10.1242/jeb.26.4.354
[4] MARTINEZ A, DEJONG J, AKIN I, et al. Bio-inspired geotechnical engineering: principles, current work, opportunities and challenges[J]. Géotechnique, 2022, 72(8): 687-705. doi: 10.1680/jgeot.20.P.170
[5] OHARA K B, MARTINEZ A. Load transfer directional of snakeskin-inspired piles during installation and pullout in sands[J]. Journal of Geotechnical and Geoenvironmental Engineering, 2022, 148(12): 04022110. doi: 10.1061/(ASCE)GT.1943-5606.0002929
[6] 李大勇, 张雨坤, 侯新宇, 等. 一种仿生鳞片式桶壁吸力基础及尺寸确定方法[P]. CN115522562A. LI Dayong, ZHANG Yukun, HOU Xinyu, et al. A suction caisson with biomimetic scale wall and size determination method[P]. CN115522562A. (in Chinese)
[7] GAO Y, QIU Y, LI B, et al. Experimental studies on the anti-uplift behavior of the suction caissons in sand[J]. Applied Ocean Research, 2013, 43: 37-45. doi: 10.1016/j.apor.2013.08.001
[8] BAI Y, LI D, ZHANG Y, et al. Bearing behavior and earth pressure distribution for modified suction caissons in sand under inclined loading[J]. Marine Georesources & Geotechnology, 2021, 39(9): 1096-1106.
[9] WANG J, LIU J, YANG Y. Model tests on failure modes and bearing capacities of suction anchors with taut mooring system[C]//The Twenty-second International Offshore and Polar Engineering Conference, OnePetro, 2012.
[10] 黎冰, 郑翔, 高玉峰, 等. 砂土中吸力式沉箱基础的最佳荷载系泊点位置模型试验研究[J]. 岩土力学, 2013, 34(9): 2521-2526. https://www.cnki.com.cn/Article/CJFDTOTAL-YTLX201309014.htm LI Bing, ZHENG Xiang, GAO Yufeng, et al. Model tests on optimal load attachment point of suction caisson foundation in sands[J]. Rock and Soil Mechanics, 2013, 34(9): 2521-2526. (in Chinese) https://www.cnki.com.cn/Article/CJFDTOTAL-YTLX201309014.htm
[11] 刘晶磊, 王建华. 软土中张紧式吸力锚破坏标准模型试验与有限元分析[J]. 岩土力学, 2013, 34(9): 2508-2514. https://www.cnki.com.cn/Article/CJFDTOTAL-YTLX201309012.htm (LIU Jinglei, WANG Jianhua. Model tests and finite element analysis of failure criterion of suction anchors with taut mooring systems in soft clay[J]. Rock and Soil Mechanics, 2013, 34(9): 2508-2514. https://www.cnki.com.cn/Article/CJFDTOTAL-YTLX201309012.htm
[12] AHMED S S, HAWLADER B C. Finite element modeling of inclined load capacity of suction caisson in sand with Abaqus/Explicit[C]//The Twenty-fourth International Ocean and Polar Engineering Conference, OnePetro, 2014.
[13] 王建华, 刘晶磊, 陈文强. 加载方向对张紧式吸力锚极限承载力的影响分析[J]. 岩土工程学报, 2012, 34(3): 385-391. http://cge.nhri.cn/article/id/14509 WANG Jianhua, LIU Jinglei, CHEN Wenqiang. Effects of loading direction on ultimate bearing capacity of suction anchors with taut mooring system[J]. Chinese Journal of Geotechnical Engineering, 2012, 34(3): 385-391. (in Chinese) http://cge.nhri.cn/article/id/14509
[14] AUBENY C, MURFF J D. Simplified limit solutions for the capacity of suction anchors under undrained conditions[J]. Ocean Engineering, 2005, 32(7): 864-877. doi: 10.1016/j.oceaneng.2004.10.006
[15] ZHAO L, GAUDIN C, O'LOUGHLIN C D, et al. Drained capacity of a suction caisson in sand under inclined loading[J]. Journal of Geotechnical and Geoenvironmental Engineering, 2019, 145(2): 04018107. doi: 10.1061/(ASCE)GT.1943-5606.0001996
[16] 黎冰, 高玉峰, 沙成明, 等. 砂土中吸力式沉箱基础的最大承载力计算方法[J]. 东南大学学报: 自然科学版, 2012, 42(6): 1201-1205. https://www.cnki.com.cn/Article/CJFDTOTAL-HJXU201602006.htm LI Bin, GAO Yufeng, SHA Chengming, et al. Calculation method for maximum bearing capacity of suction caisson foundation in sand[J]. Journal of Southeast University (Natural Science Edition), 2012, 42(6): 1201-1205. (in Chinese) https://www.cnki.com.cn/Article/CJFDTOTAL-HJXU201602006.htm
[17] LIU H X, WANG C, ZHAO Y B. Analytical study of the failure mode and pullout capacity of suction anchors in clay[J]. Ocean Systems Engineering, 2013, 3(2): 79-95. doi: 10.12989/ose.2013.3.2.079
[18] LIU H X, PENG J, ZHAO Y. Analytical study of the failure mode and pullout capacity of suction anchors in sand[J]. Ocean Systems Engineering, 2015, 5(4): 279-299. doi: 10.12989/ose.2015.5.4.279
[19] PENG J, LIU H, ZHAO Y, et al. Failure mode and pullout capacity of anchor piles in soils with cohesive and cohesionless properties[J]. Marine Georesources & Geotechnology, 2020, 38(9): 1056-1069.
[20] FRANKE E. Applicability of 1-g model tests in three cases[C]// Proceedings of the 12th International Conference on Soil Mechanics and Foundation Engineering, Rotterdam, Netherlands: A Balkema, 1989: 915-918.
[21] HONG S, VICENT S, GU K Y, et al. Effect of drainage condition on the pullout characteristics of bucket foundations in sand[J]. Ocean Engineering, 2022, 260: 111994. doi: 10.1016/j.oceaneng.2022.111994
[22] 黎冰, 郑翔, 高玉峰, 等. 砂土中吸力式沉箱基础抗拔承载特性试验研究[J]. 岩土工程学报, 2013, 35(5): 902-907. http://cge.nhri.cn/article/id/15065 (LI Bing, DENG Xiang, GAO Yufeng, et al. Model tests on pull-out capacity of suction caisson foundation in sand[J]. Chinese Journal of Geotechnical Engineering, 2013, 35(5): 902-907. http://cge.nhri.cn/article/id/15065
[23] LIU H, ZHANG W, LIU C, et al. Movement direction of drag anchors in seabed soils[J]. Applied Ocean Research, 2012, 34: 78-95. doi: 10.1016/j.apor.2011.09.007
[24] 建筑地基基础设计规范: GB 50007—2011[S]. 北京: 中国计划出版社, 2012. Code for Design of Building Foundation: GB 50007—2011[S]. Beijing: China Planning Press, 2012. (in Chinese)
-
其他相关附件
计量
- 文章访问数: 357
- HTML全文浏览量: 45
- PDF下载量: 83
