边坡场地勘探方案预期效果评价及优化方法

Expected Effectiveness Evaluation and Optimization Methods of Slope Site Investigation Program

摘要

目前的勘探方案预期效果评价指标常未能反映物理过程且参数较难确定.此外，勘探方案优化框架中的勘探点布置策略往往依赖位置关系并需事先确定勘探范围.为解决上述问题，以不排水抗剪强度参数为例，提出并采用安全系数的均方根误差折减率期望（expected reduction rate of the root mean square error， ERRS）量化因融合参考勘探数据而导致的安全系数评估结果向参考安全系数集中效果的期望提升程度，并将其作为勘探方案预期效果评价的指标.此外，结合该指标和贪心算法构建了以优化勘探位置和数量为目的的勘探方案优化框架.ERRS指标计算过程采用乔列斯基分解中点法和改进贝叶斯更新方法离散参数完全及条件随机场实现，并基于多重二阶响应面代理模型替代确定性空间变异边坡稳定性分析，有效提高指标计算精度和效率.不排水饱和黏土边坡案例显示：提出的ERRS指标能够在无需确定复杂参数的情况下，获得与其他指标接近的评价结果；所构建的勘探方案优化框架能够在不事先确定勘探范围的情况下，得到指定勘探数量下更优的勘探点布置，进而获得更节省成本且预期效果较好的勘探方案.提出的指标和优化框架可为实际边坡工程场地勘探方案评价及优化设计提供参考.

Abstract

Current evaluation indicators for depicting the expected effectiveness of exploration program often fail to reflect physical process, and the involved parameters are difficult to determine. Besides, the exploration layout strategies in existing optimization frameworks usually rely on the positional relationship, and require determination of the investigation range in advance. To solve the above problems, taking the undrained shear strength parameter as an example, in this paper it proposes and employs the expected reduction rate of the root mean square error (ERRS) to quantify the expected degree of improvement in the concentration effect of factor of safety assessment results towards the reference factor of safety due to the incorporation of reference exploration data, which serves as an indicator to evaluate the expected effectiveness of exploration program. Additionally, a framework for optimizing exploration program is constructed using this indicator in conjunction with a greedy algorithm, aimed at optimizing the locations and number of exploration points. The ERRS indicator calculation process employs the Cholesky decomposition-based midpoint method and the improved Bayesian updating method to discretize both unconditional and conditional random fields. The multiple second-order response surface is used as a surrogate model to replace the deterministic slope stability analysis of spatially variable slopes, significantly enhancing the precision and efficiency of the indicator calculation. Example of undrained saturated clay slope shows: the proposed ERRS indicator can achieve evaluation results close to other indicators without the need to determine complex parameters; the constructed exploration program optimization framework can provide an improved arrangement of exploration points under a specified number of explorations, resulting in more cost-effective and better-performing exploration program without the need to predefine the exploration range. The proposed indicator and optimization framework can serve as a reference for evaluating and optimally designing exploration program in practical slope engineering projects.

Graphical abstract

关键词

边坡 / 场地勘探方案优化 / 空间变异性 / 贝叶斯更新 / 误差折减率期望 / 工程地质.

Key words

slope / optimization of site exploration program / spatial variability / Bayesian updating / expected reduction rate of error / engineering geology

0 引言

土体是复杂地质过程和沉积环境的产物，其参数存在固有空间变异性（韩吉伟等， 2020；孙志豪等， 2021； 仉文岗等， 2023）.为了解边坡地层性质并准确评估边坡安全性，从而给后期边坡效益成本分析提供参考，在工程现场设置勘探点进行十字板剪切或原位静力触探等岩土测试试验获得勘探数据十分重要.

勘探数据的质量与所采取的勘探方案密切相关（Cao et al.， 2016； 田密与盛小涛， 2019），勘探方案的优劣影响到边坡稳定性评估的准确性和可靠性.通过某种指标量化勘探方案的价值，进而评估勘探方案的预期效果是目前研究热点.虽然众多学者在勘探方案评价指标构建及计算方面开展了较多有益的研究工作，但是依然存在以下两方面问题.一方面是评价指标形式复杂，部分参数难以确定.例如，Jiang et al.（2020）基于信息价值理论，使用勘探数据降低工程造价和边坡失效风险的预期值作为指标定量评价边坡勘探钻孔方案的优劣，尽管该指标计算过程与边坡设计决策过程相适应，但是需要额外估计边坡失效造成的损失等参数，降低了指标的适用性.另一方面是评价指标物理意义不清晰且不完善.例如，胡金政等（2023）基于机器学习模型特征选择理论，使用边坡安全系数响应面模型的模型误差减小程度作为勘探方案优劣评价的指标，尽管该指标形式简便，但仅使用部分或不使用土体参数训练响应面模型，这样训练的响应面模型无实际物理意义.Li et al.（2016b）和Hu et al.（2023）基于安全系数样本标准差构建勘探方案评价指标，但仅考虑了勘探数据对安全系数统计特征影响的一个方面.

工程勘察中通过合理布置勘探点可以获取关键的地层数据，从而提高边坡稳定性评估的准确性，这为后期边坡设计和成本效益分析提供了重要的参考依据.然而，由于受勘探成本及现场条件的限制，在实际工程中只能进行少量的代表性勘探工作，获取有限的勘探数据.因此，开展勘探位置和勘探数量等的优化研究十分必要.目前许多研究旨在利用岩土工程剖面中实测勘探数据预测未采样位置的地层性质.例如Cai et al.（2019）提出了一种三维条件随机场方法，可仅基于有限数量的圆锥贯入试验数据估计未采样位置的各向异性土阻力；Shi and Wang（2021）提出并采用IC⁃XGBoost模型预测未采样地点的土体类型和绘制地层剖面.在上述应用中，不同位置处地层性质均具有相同的重要性.目前大多数在岩土工程设计中的勘探方案优化研究依然遵循等距或等分勘探模式（Li et al.， 2016b； Yang et al.， 2019； Zhang and Wang， 2023）.虽然常用的等距勘探点布置模式在一定程度上可行，但是这种布置模式粗略地认为场地内各处勘探价值均相等，忽视了岩土工程破坏机制对勘探位置选取的影响，导致无法最大限度发挥每个勘探点的作用.Li et al.（2016a）和Deng et al.（2017）研究表明，滑动面附近区域土体性质决定了滑动面形式，所以在这些区域布置勘探点可以更准确地估计边坡安全系数并减小安全系数标准差，但是作为最佳的勘探范围的边坡破坏机制影响区域通常无法在勘探之前有效确定.

为此，本文以不排水抗剪强度参数为例，系统研究边坡场地勘探方案预期效果评价及优化问题.首先，建立表征不排水抗剪强度空间变异性的完全随机场模型和条件随机场模型，以及边坡稳定分析的多重二阶响应面代理模型，为评价指标的构建和优化奠定基础.然后，提出采用安全系数均方根误差折减率期望（expected reduction rate of the root mean square error， ERRS）作为评价勘探方案预期效果的指标，建立基于勘探点优化布置策略的勘探方案优化框架，并给出利用该指标评价勘探方案预期效果的计算流程.最后，通过不排水饱和黏土边坡案例分析验证了本文所提的评价指标和优化框架的实用性和有效性.

1 不排水抗剪强度参数随机场模拟方法

1.1　不排水抗剪强度参数非平稳完全随机场模型

在土质边坡场地中，由于受不同沉积、后沉积、化学风化和搬运作用与荷载历史等因素的影响，不同埋深处土体特性参数不仅呈现一定的空间变异性，而且还存在明显的非平稳分布特征（蒋水华等， 2018）.例如，不排水抗剪强度的均值和标准差沿埋深方向呈现逐渐增加的趋势（Asaoka and A⁃Grivas， 1982）.本文使用非平稳随机场模型表征不排水抗剪强度的非平稳分布特征，其计算表达式如下：

s_u\left(x,z\right)=s_{u\mathrm{0}}+b\gamma z\mathrm{e}\mathrm{x}\mathrm{p}\left[w\left(x,z\right)\right]，(1)

式（1）中：s_u\left(x,z\right)为不排水抗剪强度，视为非平稳随机场；s_{u\mathrm{0}}为地表处土体不排水抗剪强度；b为趋势分量，表示不排水抗剪强度均值沿埋深方向的线性变化率；\gamma为土体重度；z为土体埋深；w\left(x,z\right)为随机波动分量，表示不排水抗剪强度在均值附近随机波动的特征，视为二维平稳高斯随机场.由于以上参数取值是从文献或类似工程资料中获取，没有考虑勘探数据的影响，为此本文将式（1）称为非平稳完全随机场模型.

其中，随机波动分量w\left(x,z\right)的二维平稳高斯随机场可通过Karhunen⁃Loève（K⁃L）级数展开法、局部平均法和乔列斯基分解中点法等进行离散.其中，乔列斯基分解中点法因其计算过程简便、易编程实现的特点得到了广泛应用（蒋水华等， 2014； 胡鸿鹏等， 2024）.通过乔列斯基分解中点法（Wu et al.， 2012）获得随机波动分量的N次随机场实现{\mathit{w}}_{\mathrm{1}}~{\mathit{w}}_N的主要步骤如下：

（1）对边坡模型进行随机场单元网格剖分，提取每个单元网格的中心点坐标\left(x_i,z_i\right)，其中i=\mathrm{1},\mathrm{2},\dots ,l，l为随机场单元网格数目；

（2）基于随机波动分量的波动范围和选定的理论自相关函数计算任意两个单元网格之间的自相关系数\rho，进而构成维度为l\times l的自相关系数矩阵\mathit{\rho }，再对其进行乔列斯基分解，得到维度为l\times l的下三角矩阵{\mathit{L}}_{\mathrm{2}}；

（3）随机产生一个维度为N\times l的独立标准正态空间随机样本矩阵\mathit{\xi }；

（4）样本矩阵\mathit{\xi }右乘上三角矩阵{\mathit{L}}_{\mathrm{2}}^{\mathrm{T}}，得到维度为N\times l的独立标准正态随机场实现{\mathit{H}}^D；

（5）基于随机波动分量的均值和标准差，通过Nataf等概率变换（Li et al.，2011）将独立标准正态随机场实现{\mathit{H}}^D变换为维度为N\times l的独立正态随机场实现\mathit{H}，\mathit{H}的每一行即为随机波动分量随机场实现{\mathit{w}}_{\mathrm{1}}~{\mathit{w}}_N.

一旦获得随机波动分量的随机场实现{\mathit{w}}_{\mathrm{1}}~{\mathit{w}}_N，便可根据式（1）得到不排水抗剪强度的N次非平稳完全随机场实现{\mathit{s}}_{u_{\mathrm{1}}}~{\mathit{s}}_{u_N}，其中{\mathit{s}}_{u_N}={\left(s_{u_{N,\mathrm{1}}},s_{u_{N,\mathrm{2}}},...,s_{u_{N,l}}\right)}^T，s_{u_{N,l}}表示l个随机场单元网格处土体不排水抗剪强度的N次实现.

1.2　基于贝叶斯更新的参数条件随机场模型

通过式（1）获得的N次不排水抗剪强度完全随机场实现{\mathit{s}}_{u_{\mathrm{1}}}~{\mathit{s}}_{u_N}可视为来自不排水抗剪强度先验概率密度函数f\left({\mathit{s}}_u\right)的N个样本，因此{\mathit{s}}_{u_{\mathrm{1}}}~{\mathit{s}}_{u_N}又可称为不排水抗剪强度先验随机场实现.一旦获得勘探数据向量\mathit{d}，便可基于贝叶斯理论更新得到不排水抗剪强度后验概率密度函数f\left({\mathit{s}}_u|\mathit{d}\right)，贝叶斯更新公式如下（Jiang et al.， 2020）：

f\left({\mathit{s}}_u|\mathit{d}\right)=\frac{f\left(\mathit{d}|{\mathit{s}}_u\right)f\left({\mathit{s}}_u\right)}{f\left(\mathit{d}\right)},(2)

式（2）中：\mathit{d}=\left(d_{\mathrm{1}},d_{\mathrm{2}},...,d_{n_d}\right)为n_d个勘探数据组成的勘探数据向量；f\left({\mathit{s}}_u\right)和f\left({\mathit{s}}_u|\mathit{d}\right)分别表示不排水抗剪强度的先验和后验概率密度函数；f\left(\mathit{d}|{\mathit{s}}_u\right)为似然函数，表示在已知不排水抗剪强度为{\mathit{s}}_u的情况下获得勘探数据\mathit{d}的概率，亦可表示为似然函数L\left({\mathit{s}}_u\right)；\mathrm{1}/f(\mathit{d})为比例常数，以确保f\left({\mathit{s}}_u|\mathit{d}\right)在整个区间上的积分为1.从f\left({\mathit{s}}_u|\mathit{d}\right)获得的N次不排水抗剪强度后验随机场实现{\mathit{s}}_{u_{\mathrm{1}}}^{\mathit{d}}~{\mathit{s}}_{u_N}^{\mathit{d}}受到了勘探数据\mathit{d}的约束，因此{\mathit{s}}_{u_{\mathrm{1}}}^{\mathit{d}}~{\mathit{s}}_{u_N}^{\mathit{d}}又称为条件随机场实现.本文将式（2）称为基于贝叶斯更新的参数条件随机场模型.

对于考虑空间变异性的高维问题，式（2）没有解析解，本文采用胡鸿鹏等（2024）提出的改进BUS（bayesian updating with subset simulation）方法数值求解式（2），从而获得不排水抗剪强度后验分布及随机场实现.BUS方法基于似然函数定义一个失效区域，以建立贝叶斯更新和结构可靠度分析之间的桥梁（Straub and Papaioannou， 2015）.该方法将贝叶斯更新问题转换为等效的结构可靠度问题，再利用子集模拟方法求解该等效结构可靠度问题.然而，当勘探数据较多时，似然函数值会变得很小，导致基于似然函数构建的失效区域接受率过低，进而导致不排水抗剪强度后验样本的重复性高.为解决这一问题，提出采用改进BUS方法通过将小失效区域的构件可靠度问题转换为含多个较大失效区域的并联可靠度系统问题（胡鸿鹏等， 2024），再进行求解以提高贝叶斯更新的计算精度和效率.

改进BUS方法将似然函数分解为多个子似然函数的乘积如下：

L\left({\mathit{s}}_u\right)={\displaystyle \prod _{k=\mathrm{1}}^{n_d}}{W}_k\left({\mathit{s}}_u\right),(3)

式（3）中：L\left({\mathit{s}}_u\right)表示似然函数，其值为n_d个勘探数据采集处（对应于勘探位置所在的随机场单元）测量误差（某处的勘探数据与某处的不排水抗剪强度模拟值的差值）的联合概率密度函数.本文假设测量误差服从均值为0、标准差为2的独立联合正态分布来构建似然函数；W_k\left({\mathit{s}}_u\right)表示利用第k组勘探数据构建的子似然函数，k=\mathrm{1},\mathrm{2},\dots ,n_d，其值为第k组勘探数据采集处测量误差的概率密度函数值.

相应地，由似然函数L\left({\mathit{s}}_u\right)构建的失效区域{\Omega }_{acc}也可以分解为数个子失效区域的交集如下：

{\Omega }_{acc}={\Omega }_{acc}^{\mathrm{1}}\bigcap {\Omega }_{acc}^{\mathrm{2}}\dots \bigcap {\Omega }_{acc}^{n_d},(4)

式中：{\Omega }_{acc}^k=\left\{G_k\left({\mathit{s}}_u\right)=U_k-c_k{W}_k\left({\mathit{s}}_u\right)\le \mathrm{0}\right\}表示由W_k\left({\mathit{s}}_u\right)构建的子失效区域；U_k表示\left[\mathrm{0,1}\right]区间的均匀分布随机变量；c_k表示W_k\left({\mathit{s}}_u\right)对应的似然函数乘子，需保证对于任意不排水抗剪强度随机场实现均满足c_k{W}_k\left({\mathit{s}}_u\right)\le \mathrm{1}.

由第1.1节可知，随机场实现矩阵可以通过独立标准正态随机样本矩阵\mathit{\xi }右乘上三角矩阵{\mathit{L}}_{\mathrm{2}}^T并经过Nataf等概率变换得到.改进BUS方法可寻找对应的多次随机场实现均落在失效区域{\Omega }_{acc}内的独立标准正态随机样本矩阵.现随机生成一个独立标准正态随机样本矩阵\mathit{\xi }，其形式如下：

\mathit{\xi }=\left[\begin{array}{cccc}{\xi }_{\mathrm{1,1}}& {\xi }_{\mathrm{1,2}}& ...& {\xi }_{\mathrm{1},l}\\ {\xi }_{\mathrm{2,1}}& {\xi }_{\mathrm{2,2}}& ...& {\xi }_{\mathrm{2},l}\\ :& :& :& :\\ {\xi }_{N,\mathrm{1}}& {\xi }_{N,\mathrm{2}}& ...& {\xi }_{N,l}\end{array}\right],(5)

该矩阵维度为N\times l，其中N与随机场实现数量相对应，l与随机场单元网格数目相对应.矩阵中的第i列代表第i个随机场单元的N个独立标准正态随机样本.改进BUS方法首先利用W_{\mathrm{1}}\left({\mathit{s}}_u\right)更新矩阵\mathit{\xi }得到{\mathit{\xi }}^{\mathrm{1}}，其对应的N次随机场实现落在失效区域{\Omega }_{acc}^{\mathrm{1}}中；然后利用W_{\mathrm{2}}\left({\mathit{s}}_u\right)更新矩阵{\mathit{\xi }}^{\mathrm{1}}得到{\mathit{\xi }}^{\mathrm{2}}，其对应的N次随机场实现落在失效区域{\Omega }_{acc}^{\mathrm{1}}\bigcap {\Omega }_{acc}^{\mathrm{2}}中.以此类推，直到更新矩阵{\mathit{\xi }}^{n_d}，其对应的N次随机场实现均落在失效区域{\Omega }_{acc}^{\mathrm{1}}\bigcap {\Omega }_{acc}^{\mathrm{2}}\dots \bigcap {\Omega }_{acc}^{n_d}内.关于改进BUS方法的详细计算步骤请参考胡鸿鹏等（2024）.

2 空间变异性边坡稳定分析

一旦获得了融合勘探数据的土体参数条件随机场实现，赋值给对应的网格单元进行边坡稳定性分析.在考虑土体参数空间变异性的情况下，通常采用蒙特卡罗（Monte Carlo simulation， MCS）方法获得大量边坡样本，并对这些边坡样本进行稳定分析得到对应的安全系数，最后根据这些安全系数的统计结果评估边坡的安全性能.其中，极限平衡法因其计算原理简单，计算精度高，成为目前边坡稳定性分析的主流方法之一.GEOSTUDIO软件中SLOPE/W模块（GEO⁃SLOPE International Ltd.， 2011）调用极限平衡法中的简化Bishop法（Bishop， 1955）为空间变异性边坡稳定性分析提供了一个有效的工具.尽管如此，由于MCS方法需要调用成千上万次分析，计算效率仍较低.为此，本文进一步采用代理模型拟合边坡土体参数与安全系数之间的函数关系，代理确定性的边坡稳定性分析.

其中，多重二阶响应面（multiple second⁃order response surfaces， MSORS）是一种有效的代理模型（蒋水华， 2014），采用二阶不含交叉项多项式展开建立每条滑动面安全系数与原始空间随机变量之间的响应面函数关系，从而得到边坡多个失效模式安全系数的多重响应面.相比于其他代理模型，多重二阶响应面具有以下优势：响应面构建过程一方面无需涉及随机场实现，另一方面不依赖随机场参数取值.前者可以有效降低模型的构建难度，后者可以有效提高模型的通用性，即构建好的模型可以同时应用于完全随机场和条件随机场的边坡安全系数计算，从而极大提升边坡安全系数计算效率.假设边坡共有N_S条潜在滑动面，其中第j条潜在滑动面安全系数的响应面函数可以表示为：

\begin{array}{l}F{S}_j\left({\mathit{s}}_u\right)={\displaystyle \sum _{i=\mathrm{1}}^{N_c}}{a}_{i,j}{\psi }_{i,j}\left({\mathit{s}}_u\right)=\\ a_{\mathrm{1},j}+{\displaystyle \sum _{i=\mathrm{2}}^{l+\mathrm{1}}}{a}_{i,j}{s_u}_{{}_{i-\mathrm{1}}}+{\displaystyle \sum _{i=l+\mathrm{2}}^{\mathrm{2}l+\mathrm{1}}}{a}_{i,j}{s}_u{}_{{}_{i-l-\mathrm{1}}}^{\mathrm{2}}\end{array},（6）

式（6）中：F{S}_j\left({\mathit{s}}_u\right)为第j条潜在滑动面对应的安全系数，其中j=\mathrm{1},\mathrm{2},\dots ,N_s，N_s为潜在滑动面的数目；{\mathit{s}}_u为不排水抗剪强度随机场实现，其中{\mathit{s}}_u={\left({s_u}_{{}_{\mathrm{1}}},{s_u}_{{}_{\mathrm{2}}},\dots ,{s_u}_{{}_l}\right)}^T，l为随机变量数目，也即随机场单元网格数目；a_{i,j}为待定系数，i=\mathrm{1},\mathrm{2},\dots ,\mathrm{2}l+\mathrm{1}；{\psi }_{i,j}()为多项式展开；N_c为待定系数数目.

为建立第j条潜在滑动面安全系数的响应面，下面对式（6）中的待定系数向量{\mathit{a}}_j进行标定.首先，采用中心复合设计方法构造式（6）的\mathrm{2}l+\mathrm{1}个样本点（Bucher and Bourgund，1990），即\mathrm{2}l+\mathrm{1}个不排水抗剪强度场{{\mathit{s}}_u}_{{}_{\mathrm{1}}}~{{\mathit{s}}_u}_{{}_{\mathrm{2}l+\mathrm{1}}}与对应的\mathrm{2}l+\mathrm{1}个边坡安全系数样本F{S}_j\left({{\mathit{s}}_u}_{{}_{\mathrm{1}}}\right)~F{S}_j\left({{\mathit{s}}_u}_{{}_{\mathrm{2}l+\mathrm{1}}}\right)组合.然后，将这些样本点代入式（6）中，采用MATLAB中的lsqcurvefit函数（MathWorks Inc.， 2020）基于最小二乘法完成待定系数向量{\mathit{a}}_j的标定.其中，{{\mathit{s}}_u}_{{}_{\mathrm{1}}}~{{\mathit{s}}_u}_{{}_{\mathrm{2}l+\mathrm{1}}}中的元素按照以下规律赋值：

（1）令{{\mathit{s}}_u}_{{}_{\mathrm{1}}}中的第i个元素{s_u}_{{}_{\mathrm{1},i}}为不排水抗剪强度均值（i=\mathrm{1},\mathrm{2},\dots ,l）；

（2）令{{\mathit{s}}_u}_{{}_{\mathrm{2}}}~{{\mathit{s}}_u}_{{}_{l+\mathrm{1}}}中的第i个元素取不排水抗剪强度均值加上其标准差的k倍，其余i-1个元素取不排水抗剪强度均值（i等于{\mathit{s}}_u的下标减1）；

（3）令{{\mathit{s}}_u}_{{}_{l+\mathrm{2}}}~{{\mathit{s}}_u}_{{}_{\mathrm{2}l+\mathrm{1}}}中的第i个元素取不排水抗剪强度均值减去标准差的k倍，其余i-1个元素取不排水抗剪强度均值（i等于{\mathit{s}}_u的下标减\left[l+\mathrm{1}\right]）.

对所有N_S条潜在滑动面重复进行待定系数向量的标定，即可得到关于N_S条潜在滑动面安全系数的N_S重响应面.将随机场实现代入到N_S重响应面中即可得到边坡N_S条潜在滑动面上的N_S个安全系数，其中最小值F{S}_{\mathrm{m}\mathrm{i}\mathrm{n}}即为边坡稳定分析结果，下文简称为安全系数FS.

3 边坡场地勘探方案预期效果评价方法

3.1　安全系数均方根误差

在通过乔列斯基分解中点法结合改进BUS方法获得N组不排水抗剪强度先验、后验随机场实现{{\mathit{s}}_u}_{{}_{\mathrm{1}}}~{{\mathit{s}}_u}_{{}_N}、{\mathit{s}}_{u_{\mathrm{1}}}^{\mathit{d}}~{\mathit{s}}_{u_N}^{\mathit{d}}之后，将它们赋值至边坡数值模型，即可通过多重二阶响应面法获得N个安全系数先验、后验样本F{S}_{\mathrm{1}}~F{S}_N、F{S}_{\mathrm{1}}^{\mathit{d}}~F{S}_N^{\mathit{d}}.这样{{\mathit{s}}_u}_{{}_{\mathrm{1}}}~{{\mathit{s}}_u}_{{}_N}、{\mathit{s}}_{u_{\mathrm{1}}}^{\mathit{d}}~{\mathit{s}}_{u_N}^{\mathit{d}}分别来自不排水抗剪强度先验、后验概率密度函数f\left({\mathit{s}}_u\right)、f\left({\mathit{s}}_u|\mathit{d}\right)，则F{S}_{\mathrm{1}}~F{S}_N、F{S}_{\mathrm{1}}^{\mathit{d}}~F{S}_N^{\mathit{d}}分别来自安全系数先验、后验概率密度函数f\left(FS\right)、f\left(FS|\mathit{d}\right).

现从f\left({\mathit{s}}_u\right)获得一次先验随机场实现作为不排水抗剪强度参考随机场实现{\mathit{s}}_u^{\mathrm{r}\mathrm{e}}，将其对应的安全系数视作参考安全系数F{S}^{\mathrm{r}\mathrm{e}}.从{\mathit{s}}_u^{\mathrm{r}\mathrm{e}}中提取勘探方案SIP对应的不排水抗剪强度的参考勘探数据{\mathit{d}}_{SIP}，进而可获得N组不排水抗剪强度后验随机场实现{\mathit{s}}_{u_{\mathrm{1}}}^{\mathit{d}}\left(SIP\right)~{\mathit{s}}_{u_N}^{\mathit{d}}\left(SIP\right)，以此获得的安全系数后验样本F{S}_{\mathrm{1}}^{\mathit{d}}\left(SIP\right)~F{S}_N^{\mathit{d}}\left(SIP\right)可视作融合参考勘探数据{\mathit{d}}_{SIP}之后对参考安全系数F{S}^{\mathrm{r}\mathrm{e}}的一系列估计.如果SIP具有较高的价值（勘探位置更佳、勘探数据更丰富），那么安全系数后验向量\mathit{F}{\mathit{S}}^{\mathit{d}}\left(SIP\right)将会更集中于F{S}^{\mathrm{r}\mathrm{e}}（\mathit{F}{\mathit{S}}^{\mathit{d}}\left(SIP\right)标准差减小的同时，均值更靠近F{S}^{\mathrm{r}\mathrm{e}}），即\mathit{F}{\mathit{S}}^{\mathit{d}}\left(SIP\right)对F{S}^{\mathrm{r}\mathrm{e}}的均方根误差（root mean squared error， RMSE）将更小，其计算表达式如下（孙珊珊， 2013）：

RMSE\left[\mathit{F}{\mathit{S}}^{\mathit{d}}\left(SIP\right),F{S}^{\mathrm{r}\mathrm{e}}\right]\approx \sqrt[]{\frac{\mathrm{1}}{N}{\displaystyle \sum _{i=\mathrm{1}}^N}{\left[F{S}_i^{\mathit{d}}\left(SIP\right)-F{S}^{\mathrm{r}\mathrm{e}}\right]}^{\mathrm{2}}}，（7）

式（7）中：\mathit{F}{\mathit{S}}^{\mathit{d}}\left(SIP\right)=\left[F{S}_{\mathrm{1}}^{\mathit{d}}\left(SIP\right)~F{S}_N^{\mathit{d}}\left(SIP\right)\right]表示N个安全系数后验样本构成的安全系数后验向量；F{S}^{\mathrm{r}\mathrm{e}}表示参考安全系数；RMSE\left[\mathit{F}{\mathit{S}}^{\mathit{d}}\left(SIP\right),F{S}^{\mathrm{r}\mathrm{e}}\right]表示安全系数后验向量相对于参考安全系数的均方根误差.同理，可得到安全系数先验向量\mathit{F}\mathit{S}对于参考安全系数F{S}^{\mathrm{r}\mathrm{e}}的均方根误差计算表达式如下：

RMSE\left(\mathit{F}\mathit{S},F{S}^{\mathrm{r}\mathrm{e}}\right)\approx \sqrt[]{\frac{\mathrm{1}}{N}{\displaystyle \sum _{i=\mathrm{1}}^N}{\left(F{S}_i-F{S}^{\mathrm{r}\mathrm{e}}\right)}^{\mathrm{2}}}，（8）

式（8）中：\mathit{F}\mathit{S}=\left(F{S}_{\mathrm{1}}^{}~F{S}_N^{}\right)表示N个安全系数先验样本向量； RMSE\left(\mathit{F}\mathit{S},F{S}^{\mathrm{r}\mathrm{e}}\right)表示安全系数先验向量相对于参考安全系数的均方根误差.

图1展示了安全系数先验样本、安全系数后验样本的分布示意图，红色线段代表了参考安全系数所处的位置.由于融合了勘探数据，安全系数后验样本比安全系数先验样本更集中于参考安全系数.同时由于SIP₁相比于SIP₂价值更高，因此对应的安全系数后验样本的集中效果更强.

3.2　安全系数均方根误差折减率期望

结合式（7）和式（8）可构建安全系数均方根误差折减率（reduction rate of RMSE， RRS），用来量化融合勘探方案SIP对应的勘探数据之后安全系数均方根误差的折减程度，即量化安全系数评估结果向参考安全系数集中效果的提升程度，其计算表达式如下：

RRS\left(SIP\right)=\frac{\left\{\begin{array}{l}RMSE\left(\mathit{F}\mathit{S},F{S}^{\mathrm{r}\mathrm{e}}\right)-\\ RMSE\left[\mathit{F}{\mathit{S}}^{\mathit{d}}\left(SIP\right),F{S}^{\mathrm{r}\mathrm{e}}\right]\end{array}\right\}}{RMSE\left(\mathit{F}\mathit{S},F{S}^{\mathrm{r}\mathrm{e}}\right)},（9）

式（9）中各项含义同式（7）和式（8）.考虑到不排水抗剪强度参考随机场实现{\mathit{s}}_u^{\mathrm{r}\mathrm{e}}具有不确定性，可进一步构建ERRS作为勘探方案评价指标，其量化了因融合勘探数据而导致的安全系数评估结果向参考安全系数集中效果的期望提升程度，其表达式如下：

\begin{array}{l}ERRS\left(SIP\right)\approx \frac{\mathrm{1}}{M}{\displaystyle \sum _{j=\mathrm{1}}^M}RR{S}_j\left(SIP\right)=\\ \frac{\mathrm{1}}{M}{\displaystyle \sum _{j=\mathrm{1}}^M}\frac{\left\{\begin{array}{l}RMSE\left(\mathit{F}\mathit{S},F{S}_j^{\mathrm{r}\mathrm{e}}\right)-\\ RMSE\left[\mathit{F}{\mathit{S}}_j^{\mathit{d}}\left(SIP\right),F{S}_j^{\mathrm{r}\mathrm{e}}\right]\end{array}\right\}}{RMSE\left(\mathit{F}\mathit{S},F{S}_j^{\mathrm{r}\mathrm{e}}\right)}\end{array}，（10）

式（10）中：RR{S}_j\left(SIP\right)表示第j个安全系数均方根误差折减率；F{S}_j^{\mathrm{r}\mathrm{e}}表示第j个参考安全系数；\mathit{F}{\mathit{S}}_j^{\mathit{d}}\left(SIP\right)=\left[F{S}_{j,\mathrm{1}}^{\mathit{d}}\left(SIP\right)~F{S}_{j,N}^{\mathit{d}}\left(SIP\right)\right]表示融合第j个参考随机场实现{\mathit{s}}_{u_j}^{\mathrm{r}\mathrm{e}}中的勘探数据{\mathit{d}}_{SIP}^j模拟后验随机场实现{\mathit{s}}_{u_{j,\mathrm{1}}}^{\mathit{d}}\left(SIP\right)~{\mathit{s}}_{u_{j,N}}^{\mathit{d}}\left(SIP\right)计算的N个安全系数后验向量；M表示参考安全系数的数量或参考随机场实现次数.图2展示了基于ERRS的勘探方案预期效果评价流程.在计算ERRS的过程中，通过利用乔列斯基分解中点法、改进BUS方法、多重二阶响应面法等提高勘探方案评价指标的计算精度和效率.

4 边坡场地勘探方案优化方法

4.1　勘探点优化布置策略

指定勘探数量下勘探点优化布置可以视作一个优化问题，即在每个勘探点的勘探位置允许取值范围之内（或称勘探位置组合所有可能取值中），找到一个能够使ERRS取最大值的勘探位置组合，构成最佳勘探方案SI{P}_{\mathrm{b}\mathrm{e}\mathrm{s}\mathrm{t}}.若预计设置5个勘探点，对应位置记为\left[L_{\mathrm{1}},L_{\mathrm{2}},L_{\mathrm{3}},L_{\mathrm{4}},L_{\mathrm{5}}\right]，则该勘探数量下勘探点优化布置的数学表达式如下（张兰， 2019）：

\begin{array}{l}\mathrm{m}\mathrm{a}\mathrm{x}ERRS\left(SIP\right),SIP=\left[L_{\mathrm{1}},L_{\mathrm{2}},L_{\mathrm{3}},L_{\mathrm{4}},L_{\mathrm{5}}\right]\\ \mathrm{s}.\mathrm{t}.\left\{\begin{array}{l}{L}_{\mathrm{1}}\in Z,l{b}_{\mathrm{1}}\le L_{\mathrm{1}}\le u{b}_{\mathrm{1}}\\ L_{\mathrm{2}}\in Z,l{b}_{\mathrm{2}}\le L_{\mathrm{2}}\le u{b}_{\mathrm{2}}\\ L_{\mathrm{3}}\in Z,l{b}_{\mathrm{3}}\le L_{\mathrm{3}}\le u{b}_{\mathrm{3}}\\ L_{\mathrm{4}}\in Z,l{b}_{\mathrm{4}}\le L_{\mathrm{4}}\le u{b}_{\mathrm{4}}\\ L_{\mathrm{5}}\in Z,l{b}_{\mathrm{5}}\le L_{\mathrm{5}}\le u{b}_{\mathrm{5}}\end{array}\right.\end{array},(11)

式（11）中：ERRS\left(SIP\right)为最优化问题的目标函数，即安全系数均方根误差折减率期望ERRS；max（）函数表示该问题为最大优化问题；SIP=\left[L_{\mathrm{1}},L_{\mathrm{2}},L_{\mathrm{3}},L_{\mathrm{4}},L_{\mathrm{5}}\right]为最优化问题的决策变量，表示由勘探位置L_{\mathrm{1}}~L_{\mathrm{5}}构成的勘探方案SIP；\mathrm{s}.\mathrm{t}.为最优化问题的约束条件，勘探位置L_{\mathrm{1}}~L_{\mathrm{5}}均为整数，分别在对应的上界u{b}_{\mathrm{1}}~u{b}_{\mathrm{5}}和下界l{b}_{\mathrm{1}}~l{b}_{\mathrm{5}}范围内进行取值.

考虑到目标函数的计算效率，采用穷举法或粒子群优化法（Kennedy and Eberhart， 1995）等智能优化算法直接从勘探位置的所有可能组合中搜寻最佳勘探位置组合以实现对式（11）的求解存在一定困难，为此本文引入贪心思想.贪心思想是一种能够得到某种度量意义下的最优解的分级处理策略，它总能做出在当前看来是最优的选择（常友渠等， 2008）.本节基于该思想将“5个勘探位置所有可能取值组合中寻得最优勘探位置组合”的最优化问题分解为“接续寻找不同勘探位置的最佳值，最后拼接成近似最优勘探位置组合”的5个子最优化问题A、B、C、D、E.

子问题A描述为：在只需布置1个勘探点的情况下，找到使得ERRS取最大值的最佳勘探位置L_{\mathrm{1}}^{\mathrm{b}\mathrm{e}\mathrm{s}\mathrm{t}}.对应的数学表达式如下：

\begin{array}{l}\mathrm{m}\mathrm{a}\mathrm{x}ERRS\left(L_{\mathrm{1}}\right)\\ \mathrm{s}.\mathrm{t}.L_{\mathrm{1}}\in Z,l{b}_{\mathrm{1}}\le L_{\mathrm{1}}\le u{b}_{\mathrm{1}}\end{array},(12)

子问题B描述为：在需要布置2个勘探点，且其中第1个勘探位置确定为L_{\mathrm{1}}^{\mathrm{b}\mathrm{e}\mathrm{s}\mathrm{t}}的情况下，找到使得ERRS取最大值的第2个最佳勘探位置L_{\mathrm{2}}^{\mathrm{b}\mathrm{e}\mathrm{s}\mathrm{t}}.对应的数学表达式如下：

\begin{array}{l}\mathrm{m}\mathrm{a}\mathrm{x}ERRS\left(\left[L_{\mathrm{1}},L_{\mathrm{2}}\right]\right)\\ \mathrm{s}.\mathrm{t}.\left\{\begin{array}{l}{L}_{\mathrm{1}}=L_{\mathrm{1}}^{\mathrm{b}\mathrm{e}\mathrm{s}\mathrm{t}}\\ L_{\mathrm{2}}\in Z,l{b}_{\mathrm{2}}\le L_{\mathrm{2}}\le u{b}_{\mathrm{2}}\end{array}\right.\end{array},(13)

以此类推，子问题E描述为：在需要布置5个勘探点，且前4个勘探位置已经确定为L_{\mathrm{1}}^{\mathrm{b}\mathrm{e}\mathrm{s}\mathrm{t}}、L_{\mathrm{2}}^{\mathrm{b}\mathrm{e}\mathrm{s}\mathrm{t}}、L_{\mathrm{3}}^{\mathrm{b}\mathrm{e}\mathrm{s}\mathrm{t}}、L_{\mathrm{4}}^{\mathrm{b}\mathrm{e}\mathrm{s}\mathrm{t}}的情况下，找到使得ERRS取最大值的第5个最佳勘探位置L_{\mathrm{5}}^{\mathrm{b}\mathrm{e}\mathrm{s}\mathrm{t}}.对应的数学表达式如下：

\begin{array}{l}\mathrm{m}\mathrm{a}\mathrm{x}ERRS\left(\left[L_{\mathrm{1}},L_{\mathrm{2}},L_{\mathrm{3}},L_{\mathrm{4}},L_{\mathrm{5}}\right]\right)\\ \mathrm{s}.\mathrm{t}.\left\{\begin{array}{l}{L}_{\mathrm{1}}=L_{\mathrm{1}}^{\mathrm{b}\mathrm{e}\mathrm{s}\mathrm{t}}\\ L_{\mathrm{2}}=L_{\mathrm{2}}^{\mathrm{b}\mathrm{e}\mathrm{s}\mathrm{t}}\\ L_{\mathrm{3}}=L_{\mathrm{3}}^{\mathrm{b}\mathrm{e}\mathrm{s}\mathrm{t}}\\ L_{\mathrm{4}}=L_{\mathrm{4}}^{\mathrm{b}\mathrm{e}\mathrm{s}\mathrm{t}}\\ L_{\mathrm{5}}\in Z,l{b}_{\mathrm{5}}\le L_{\mathrm{5}}\le u{b}_{\mathrm{5}}\end{array}\right.\end{array},(14)

子问题\mathrm{A}~\mathrm{E}对应的最佳勘探位置L_{\mathrm{1}}^{\mathrm{b}\mathrm{e}\mathrm{s}\mathrm{t}}~L_{\mathrm{5}}^{\mathrm{b}\mathrm{e}\mathrm{s}\mathrm{t}}可以分别通过绘制前几个勘探点已布置在最佳勘探位置情况下ERRS随第k个勘探位置L_k^{}的变化曲线寻找得到.最佳勘探位置L_{\mathrm{1}}^{\mathrm{b}\mathrm{e}\mathrm{s}\mathrm{t}}~L_{\mathrm{5}}^{\mathrm{b}\mathrm{e}\mathrm{s}\mathrm{t}}的组合即为勘探数量为5时的最优勘探方案SI{P}_{\mathrm{b}\mathrm{e}\mathrm{s}\mathrm{t}}=\left[L_{\mathrm{1}}^{\mathrm{b}\mathrm{e}\mathrm{s}\mathrm{t}},L_{\mathrm{2}}^{\mathrm{b}\mathrm{e}\mathrm{s}\mathrm{t}},L_{\mathrm{3}}^{\mathrm{b}\mathrm{e}\mathrm{s}\mathrm{t}},L_{\mathrm{4}}^{\mathrm{b}\mathrm{e}\mathrm{s}\mathrm{t}},L_{\mathrm{5}}^{\mathrm{b}\mathrm{e}\mathrm{s}\mathrm{t}}\right].

4.2　勘探方案优化框架

若采用贪心思想对式（11）进行求解，不仅可以得到勘探数量为5时的最优勘探位置组合，还可以得到勘探数量为1~4时的最优勘探位置组合.通过结合ERRS随L_k^{}或勘探数量变化的趋势分析，可以构建基于勘探点优化布置策略的场地勘探方案优化框架.在不确定勘探范围的情况下，该框架能够优化勘探位置和勘探数量，从而在尽量降低成本的情况下获得最大化价值的勘探方案.图3展示了利用该框架优化勘探方案的流程.

5 案例分析

5.1　边坡模型

本节采用不排水饱和黏土边坡案例说明提出的勘探方案预期效果评价和优化方法的实用性和有效性.边坡模型尺寸及随机场网格划分情况如图4所示，蒋水华等（2018）和胡金政等（2023）基于该边坡进行了勘探方案的评价和优化工作，相关计算结果可以用于对比分析.边坡模型共计划分910个单元格，每个单元格大小为2 m×0.5 m.土体重度视作常数取值为20 kN/m³.考虑土体不排水抗剪强度的空间变异性及其随埋深方向的变化趋势，通过式（1）将其模拟为非平稳随机场，随机场模型参数取值见表1（蒋水华等， 2018）.

5.2　边坡稳定性分析

图5展示了利用乔列斯基分解中点法生成的不排水抗剪强度的均值场和典型先验随机场实现的云图，并在其上标注了边坡安全系数和临界滑动面位置，其中基于均值的边坡安全系数为1.58，与蒋水华等（2018）计算的安全系数1.55基本吻合.

接着建立多重二阶响应面代理模型进行空间变异性边坡稳定性分析.模型构建过程中需要通过敏感性分析确定潜在滑动面数量N_S，根据图6所示结果将N_S确定为2 816.图7比较了使用极限平衡方法计算的安全系数F{S}_{\mathrm{L}\mathrm{E}\mathrm{M}}和多重二阶响应面法计算的安全系数F{S}_{\mathrm{M}\mathrm{S}\mathrm{O}\mathrm{R}\mathrm{S}}，数据点几乎都位于45° 线附近，决定系数R^{\mathrm{2}}为0.99，这表示所构建的代理模型具有良好的先验安全系数预测性能.

5.3　勘探方案预期效果评价

设置两个典型等距布置的勘探方案：（1）从坡顶到坡脚等距布置勘探点，勘探位置分别取21 m、27 m、33 m、39 m；（2）整个边坡区域等距布置勘探点，勘探位置分别取1 m、21 m、39 m、59 m.勘探点在边坡上的分布情况见图8.

采用ERRS按照图2所示的计算流程对以上两个勘探方案进行预期效果评价，根据图9和图10所示的敏感性分析结果，可确定需要的先验（或后验）随机场实现数量N和参考随机场实现数量M分别为5 000和30.同时根据图10的计算结果，得到勘探方案①对应的ERRS为0.66、勘探方案②对应的ERRS为0.39，表明勘探方案①具有更佳的预期效果.

5.4　勘探方案预期效果评价结果分析

为探究导致勘探方案①和勘探方案②评价结果存在差异的原因，首先从第1个参考随机场{\mathit{s}}_{u_{\mathrm{1}}}^{\mathrm{r}\mathrm{e}}中获得两个勘探方案的参考勘探数据{\mathit{d}}_{SI{P}_{①}}^{\mathrm{1}}和{\mathit{d}}_{SI{P}_{②}}^{\mathrm{1}}.然后，基于改进BUS方法获得后验随机场实现{\mathit{s}}_{u_{\mathrm{1,1}}}^{\mathit{d}}\left(SI{P}_{①}\right)~{\mathit{s}}_{u_{\mathrm{1},N}}^{\mathit{d}}\left(SIP①\right)和{\mathit{s}}_{u_{\mathrm{1,1}}}^{\mathit{d}}\left(SI{P}_{②}\right)~{\mathit{s}}_{u_{\mathrm{1},N}}^{\mathit{d}}\left(SI{P}_{②}\right)，并计算和绘制这些后验随机场实现对{\mathit{s}}_{u_{\mathrm{1}}}^{\mathrm{r}\mathrm{e}}的均方根误差云图，分别如图11a和图11b所示.图11中红色曲线表示利用参数均值边坡搜寻到的临界滑移面，其包围区域大致表征滑坡潜在影响区域.同时，计算得到两种勘探方案的第1个安全系数均方根误差折减率RR{S}_{\mathrm{1}}\left(SI{P}_{①}\right)和RR{S}_{\mathrm{1}}\left(SIP②\right)分别为0.68和0.48.观察图11a和图11b可以发现，勘探方案①对应的不排水抗剪强度均方根误差在滑坡潜在影响区域更小，而其安全系数均方根误差折减率则更高.这表明该方案更充分了解滑坡潜在影响区域土体不排水抗剪强度，将更有利于安全系数评估结果向参考安全系数集中，从而充分发挥勘探方案的性能，理想的勘探范围应为滑坡潜在影响区域.在计算RR{S}_{\mathrm{1}}\left(SI{P}_{①}\right)和RR{S}_{\mathrm{1}}\left(SI{P}_{②}\right)过程中分别利用极限平衡方法和多重二阶响应面法计算获得{\mathit{s}}_{u_{\mathrm{1,1}}}^{\mathit{d}}\left(SI{P}_{①}\right)~{\mathit{s}}_{u_{\mathrm{1,100}}}^{\mathit{d}}\left(SI{P}_{①}\right)和{\mathit{s}}_{u_{\mathrm{1,1}}}^{\mathit{d}}\left(SI{P}_{②}\right)~{\mathit{s}}_{u_{\mathrm{1,100}}}^{\mathit{d}}\left(SI{P}_{②}\right)的安全系数，并将两种方法的计算结果进行比较，如图12所示.结果表明，两种方法的计算结果基本一致，表明多重二阶响应面也具有良好的后验安全系数预测能力.

图13绘制了ERRS随单个勘探点勘探位置的变化情况，图13中蓝色实线表示边坡轮廓，红色虚线与边坡轮廓构成了大致的边坡潜在影响区域，蓝色虚线表示最优勘探位置，其位于距边坡左边界25 m处.由图13可以发现，ERRS与勘探点相对滑坡潜在滑移面的位置高度相关：在25 m之前，ERRS随着勘探位置靠近滑坡潜在影响区域从0.05逐渐上升到0.48，而在25 m之后，ERRS随着勘探位置远离滑坡潜在影响区域从0.48逐渐下降到0.这是因为勘探点与滑坡潜在影响区域中心的相对距离越近，则该勘探点越能够通过波动范围影响到该区域更大范围内的不排水抗剪强度，安全系数评估结果则越向参考安全系数集中，进而使得勘探点的ERRS结果更大.而如果勘探点与滑坡潜在影响区域中心的相对距离越远，则越难以了解该区域内的土体不排水抗剪强度，越难影响安全系数的评估，从而ERRS更小甚至为0.图13同样绘制了基于机器学习模型特征选择理论构建的勘探方案评价指标V（胡金政等， 2023）及信息价值指标E(VOI)（蒋水华等，2018）随勘探位置的变化情况.ERRS具有与V相类似的变化趋势，但是相比于后者在计算过程中仅根据边坡部分位置的土体参数训练模型预测边坡的安全系数，ERRS指标和计算过程中使用的各种方法具有更明确的物理意义.相比E(VOI)，ERRS具有更光滑的变化趋势，并且该指标在计算过程中无需对实际工程中难以确定的边坡失稳破坏后果进行估算.

5.5　勘探方案优化

采用基于勘探点优化布置策略的勘探方案优化框架，按照图3所示勘探方案优化流程进行勘探方案优化.图14a展示了优化过程中绘制的ERRS随第k个勘探位置L_k^{}的变化曲线，黑色曲线表示ERRS随第1个勘探位置L_{\mathrm{1}}^{}的变化，最终将第1个最佳勘探位置L_{\mathrm{1}}^{\mathrm{b}\mathrm{e}\mathrm{s}\mathrm{t}}取为25 m，类似地，根据其他曲线可将L_{\mathrm{2}}^{\mathrm{b}\mathrm{e}\mathrm{s}\mathrm{t}}取为33 m、L_{\mathrm{3}}^{\mathrm{b}\mathrm{e}\mathrm{s}\mathrm{t}}取为21 m、L_{\mathrm{4}}^{\mathrm{b}\mathrm{e}\mathrm{s}\mathrm{t}}取为29 m.这些变化曲线存在以下特征：（1）随着已确定的最佳勘探位置的数量越来越多，曲线的变化趋势逐渐减缓，当勘探数量为4时基本达到稳定，这是因为先确定的最佳勘探位置已经在很大程度上减小了滑坡潜在影响区域内不排水抗剪强度均方根误差，为安全系数评估结果向参考安全系数集中做出了重要贡献，而之后确定的最佳勘探位置在这一方面的贡献越来越小.（2）每条曲线都被已确定的最佳勘探位置分割成先增后减的数个部分，当前的最佳勘探位置往往不位于离上一个最佳勘探点最近的位置，例如L_{\mathrm{2}}^{\mathrm{b}\mathrm{e}\mathrm{s}\mathrm{t}}为33 m，并非距离L_{\mathrm{1}}^{\mathrm{b}\mathrm{e}\mathrm{s}\mathrm{t}}更近的29 m，这是因为两勘探点间距过小会导致勘探点折减土体不排水抗剪强度均方根误差的作用发生重叠.

图14b展示了优化过程中绘制的ERRS随勘探数量的变化曲线，随着勘探数量的增加，ERRS的增长趋势逐渐放缓，当勘探数量为4时基本稳定至0.7，据此最终将勘探数量确定为4（也可根据图14a中ERRS随L_k^{}变化的趋势将勘探数量确定为4）.图15展示了最终优化得到的勘探方案，其勘探点基本分布在滑坡潜在影响区域的中心.

5.6　勘探方案优化结果的对比分析

为进一步验证基于勘探点优化布置策略的勘探方案优化框架（本文框架）的有效性，分别基于坡顶到坡脚范围内和整个边坡区域范围内的勘探点等距布置策略搭建了勘探方案优化框架（小范围等距框架和大范围等距框架），以指导勘探方案优化设计.等距框架仅使用勘探点的相对位置关系确定指定勘探数量下的勘探位置，然后通过评价指标随勘探数量的变化情况确定勘探数量.小范围等距框架和大范围等距框架在不同勘探数量下的勘探点布置情况如图16和图17所示.

图18展示了采用这3个框架优化勘探方案过程中绘制的ERRS随勘探数量的变化情况.就同一勘探数量下勘探点布置的ERRS计算结果而言，相比于小范围等距框架和大范围等距框架，由于本文所采用的勘探点优化布置策略直接根据评价指标计算结果指导勘探点位置的布置，因此在任意勘探数量下的勘探点布置总是具有更高的ERRS计算结果.此外，尽管小范围等距框架和大范围等距框架都采用勘探点等距布置策略，但由于勘探范围的不同，二者在同一勘探数量下勘探点布置的ERRS可能存在较大差异.以勘探数量取2的情况为例，小范围等距框架在此勘探数量下的勘探点布置的ERRS为0.4，而大范围等距框架在此勘探数量下的勘探点布置的ERRS只有0.05.前者的勘探点大致落在滑坡潜在影响区域边缘，在一定程度上利用了勘探点了解该区域不排水抗剪强度的作用，有利于安全系数评估向参考安全系数集中，最终ERRS计算结果为0.4；而后者的勘探点落在滑坡潜在影响区域之外，未能了解该区域不排水抗剪强度，基本不影响安全系数的评估，使得ERRS仅为0.05.

就ERRS随勘探数量的变化趋势而言，这3个勘探方案优化框架均能使ERRS随着勘探数量的增加而总体增大.然而，由于勘探点等距布置策略的局限，后两个优化框架在优化过程中ERRS的增大趋势稍显曲折，最终优化得到的勘探方案需要包含更多的勘探数量才能达到图15所示的勘探方案的预期效果水平，这显然会增加大量的勘探成本.

6 结论

本文针对当前边坡场地勘探方案评价指标和优化设计存在的不足，提出基于ERRS的勘探方案预期效果的评价指标，构建了基于勘探点优化布置策略的勘探方案优化方法，并给出了详细评价步骤和优化流程.通过不排水饱和黏土边坡案例验证了提出方法的实用性和有效性.主要结论如下.

（1）本文提出了基于边坡ERRS的勘探方案预期效果的评价指标，物理意义明确，可有效量化因融合参考勘探数据而导致的安全系数评估结果向参考安全系数集中效果的期望提升程度，期望提升程度更大（ERRS更大），则勘探方案更优.此外，指标计算过程采用乔列斯基分解中点法和改进BUS方法离散完全及条件随机场实现，使用多重二阶响应面代理对应空间变异性边坡稳定性分析，进一步扩展了这些方法在岩土工程领域中的应用并提高了评价指标计算效率和精度.在不排水饱和黏土边坡案例中，评价结果与对滑坡潜在影响区域土体不排水抗剪强度的了解程度呈正相关，该指标可以在无需确定复杂参数的情况下获得与其他指标相一致的评价结果.

（2）提出了基于勘探点优化布置策略的边坡场地勘探方案优化方法.该方法将指定勘探数量下的勘探点布置视为以评价指标为目标函数、勘探位置组合为决策变量的最优化问题，利用贪心思想求解该问题，从而实现勘探点的优化布置，并在此基础上通过分析评价指标的变化趋势实现勘探数量的优化.针对不排水饱和黏土边坡，该框架在任意勘探数量下优化得到的勘探点布置均具有较好的预期效果，且能在无需规定勘探范围的情况下获得成本更节省且预期效果更好的勘探方案.

参考文献

[1] Asaoka,A.,A⁃Grivas,D.,1982.Spatial Variability of the Undrained Strength of Clays.Journal of the Geotechnical Engineering Division,108(5):743-756.https://doi.org/10.1061/ajgeb6.0001292

[2] Bishop,A.W.,1955.The Use of the Slip Circle in the Stability Analysis of Slopes.Géotechnique,5(1):7-17.https://doi.org/10.1680/geot.1955.5.1.7

[3] Bucher,C.G.,Bourgund,U.,1990.A Fast and Efficient Response Surface Approach for Structural Reliability Problems.Structural Safety,7(1):57-66.https://doi.org/10.1016/0167⁃4730(90)90012⁃e

[4] Cai,Y.M.,Li,J.H.,Li,X.Y.,et al.,2019.Estimating Soil Resistance at Unsampled Locations Based on Limited CPT Data.Bulletin of Engineering Geology and the Environment,78(5):3637-3648.https://doi.org/10.1007/s10064⁃018⁃1318⁃2

[5] Cao,Z.J.,Wang,Y.,Li,D.Q.,2016.Quantification of Prior Knowledge in Geotechnical Site Characterization.Engineering Geology,203:107-116.https://doi.org/10.1016/j.enggeo.2015.08.018

[6] Chang,Y.Q.,Xiao,G.Y.,Zeng,M.,2008.Exploration of Greedy Algorithm.Journal of Chongqing Electric Power College,13(3):40-42,47(in Chinese with English abstract).

[7] Deng,Z.P.,Li,D.Q.,Qi,X.H.,et al.,2017.Reliability Evaluation of Slope Considering Geological Uncertainty and Inherent Variability of Soil Parameters.Computers and Geotechnics,92:121-131.https://doi.org/10.1016/j.compgeo.2017.07.020

[8] GEO⁃SLOPE International Ltd.,2011.Stability Modeling with SLOPE/W 2007 Version.Fourth Edition.Translated by CnTech.Metallurgical Industry Press,Beijing (in Chinese).

[9] Han,J.W.,Liu,X.M.,Yang,Y.H.,et al.,2020.Study on the Horizontal Bearing Capacity of Single Pile Foundation Considering Spatial Variability of Soil Strength.Hydro⁃Science and Engineering,(6):108-114(in Chinese with English abstract).

[10] Hu,H.P.,Jiang,S.H.,Chen,D.,et al.,2024.Probabilistic Back Analysis of Slope Parameters and Reliability Evaluation Using Improved Bayesian Updating Method.Rock and Soil Mechanics,45(3):835-845 (in Chinese with English abstract).

[11] Hu,J.Z.,Zhang,J.,Huang,H.W.,et al.,2023.Value of Information Assessment and Optimization of Slope Boreholes.Earth Science,48(5):1977-1988 (in Chinese with English abstract).

[12] Hu,J.Z.,Zheng,J.G.,Zhang,J.,et al.,2023.Bayesian Framework for Assessing Effectiveness of Geotechnical Site Investigation Programs.ASCE⁃ASME Journal of Risk and Uncertainty in Engineering Systems,Part A:Civil Engineering,9(1):04022054.https://doi.org/10.1061/ajrua6.0001278

[13] Jiang,S.H.,2014.A Non⁃Intrusive Stochastic Method for Slope Reliability in Hydroelectricity Engineering (Dissertation).Wuhan University,Wuhan (in Chinese with English abstract).

[14] Jiang,S.H.,Li,D.Q.,Zhou,C.B.,et al.,2014.Slope Reliability Analysis Considering Effect of Autocorrelation Functions.Chinese Journal of Geotechnical Engineering,36(3):508-518 (in Chinese with English abstract).

[15] Jiang,S.H.,Liu,X.,Yao,R.Z.,et al.,2018.Optimization Design Approach for Layout Scheme of Slope Boreholes Based on Bayesian Updating and Value of Information Analysis.Chinese Journal of Geotechnical Engineering,40(10):1871-1879 (in Chinese with English abstract).

[16] Jiang,S.H.,Papaioannou,I.,Straub,D.,2020.Optimization of Site⁃Exploration Programs for Slope⁃Reliability Assessment.ASCE⁃ASME Journal of Risk and Uncertainty in Engineering Systems,Part A:Civil Engineering,6(1):04020004.https://doi.org/10.1061/ajrua6.0001042

[17] Kennedy,J.,Eberhart,R.,1995.Particle Swarm Optimization.Proceedings of ICNN’95 - International Conference on Neural Networks.Perth,WA,Australia.1942-1948.https://doi.org/10.1109/icnn.1995.488968

[18] Li,D.Q.,Chen,Y.F.,Lu,W.B.,et al.,2011.Stochastic Response Surface Method for Reliability Analysis of Rock Slopes Involving Correlated Non⁃Normal Variables.Computers and Geotechnics,38(1):58-68.https://doi.org/10.1016/j.compgeo.2010.10.006

[19] Li,D.Q.,Qi,X.H.,Cao,Z.J.,et al.,2016a.Evaluating Slope Stability Uncertainty Using Coupled Markov Chain.Computers and Geotechnics,73:72-82.https://doi.org/10.1016/j.compgeo.2015.11.021

[20] Li,Y.J.,Hicks,M.A.,Vardon,P.J.,2016b.Uncertainty Reduction and Sampling Efficiency in Slope Designs Using 3D Conditional Random Fields.Computers and Geotechnics,79:159-172.https://doi.org/10.1016/j.compgeo.2016.05.027

[21] MathWorks Inc.,2020.MATLAB Version:9.9.0.1467703(R2020b).The MathWorks Inc,Natick,Massachusetts.

[22] Shi,C.,Wang,Y.,2021.Development of Subsurface Geological Cross⁃Section from Limited Site⁃Specific Boreholes and Prior Geological Knowledge Using Iterative Convolution XGBoost.Journal of Geotechnical and Geoenvironmental Engineering,147(9):04021082.https://doi.org/10.1061/(asce)gt.1943⁃5606.0002583

[23] Straub,D.,Papaioannou,I.,2015.Bayesian Updating with Structural Reliability Methods.Journal of Engineering Mechanics,141(3):04014134.https://doi.org/10.1061/(asce)em.1943⁃7889.0000839

[24] Sun,S.S.,2013.Statistical Inference Methods Based on Item Response Theory and Cognitive Diagnosis(Dissertation).Northeast Normal University,Changchun(in Chinese with English abstract).

[25] Sun,Z.H.,Tan,X.H.,Sun,Z.B.,et al.,2021.Reliability of Spatially Variable Earth Slopes Based on the Upper Bound Analysis.Rock and Soil Mechanics,42(12):3397-3406 (in Chinese with English abstract).

[26] Tian,M.,Sheng,X.T.,2019.Method for Determining Minimum Test Data Quantity for Geotechnical Engineering Investigation.Rock and Soil Mechanics,40(Supp.1):400-408(in Chinese with English abstract).

[27] Wu,S.H.,Ou,C.Y.,Ching,J.,et al.,2012.Reliability⁃Based Design for Basal Heave Stability of Deep Excavations in Spatially Varying Soils.Journal of Geotechnical and Geoenvironmental Engineering,138(5):594-603.https://doi.org/10.1061/(asce)gt.1943⁃5606.0000626

[28] Yang,R.,Huang,J.S.,Griffiths,D.V.,et al.,2019.Optimal Geotechnical Site Investigations for Slope Design.Computers and Geotechnics,114:103111.https://doi.org/10.1016/j.compgeo.2019.103111.

[29] Zhang,L.,2019.Research on Quantum Behaved Particle Swarm Optimization Algorithm and Its Applications for Single⁃Objective and Multi⁃Objective Optimization Problems (Dissertation).Northwestern Polytechnical University,Xi’an(in Chinese with English abstract).

[30] Zhang,L.,Wang,L.,2023.Optimization of Site Investigation Program for Reliability Assessment of Undrained Slope Using Spearman Rank Correlation Coefficient.Computers and Geotechnics,155:105208.https://doi.org/10.1016/j.compgeo.2022.105208

[31] Zhang,W.G.,Meng,F.S.,He,C.J.,et al.,2023.Stability Analysis of Slope with Spatial Variability in Karst Area under Building Load.Journal of Disaster Prevention and Mitigation Engineering,43(2):316-323 (in Chinese with English abstract).

[32] 常友渠,肖贵元,曾敏,2008.贪心算法的探讨与研究.重庆电力高等专科学校学报,13(3):40-42,47.

[33] GEO-SLOPE International Ltd.,2011.边坡稳定性分析软件SLOPE/W用户指南.中仿科技公司,译.北京:冶金工业出版社.

[34] 韩吉伟,刘晓明,杨月红,等,2020.考虑土体强度空间变异性的单桩水平承载力研究.水利水运工程学报 (6):108-114.

[35] 胡鸿鹏,蒋水华,陈东,等,2024.基于改进贝叶斯更新方法的边坡参数概率反分析及可靠度评估.岩土力学,45(3):835-845.

[36] 胡金政,张洁,黄宏伟,等,2023.边坡勘察钻孔信息价值评价及优化布置方法.地球科学,48(5):1977-1988.

[37] 蒋水华,2014.水电工程边坡可靠度非侵入式随机分析方法(博士学位论文).武汉:武汉大学.

[38] 蒋水华,李典庆,周创兵,等,2014.考虑自相关函数影响的边坡可靠度分析.岩土工程学报,36(3):508-518.

[39] 蒋水华,刘贤,尧睿智,等,2018.基于贝叶斯更新和信息量分析的边坡钻孔布置方案优化设计方法.岩土工程学报,40(10):1871-1879.

[40] 孙珊珊,2013.项目反应理论与认知诊断的统计推断方法(博士学位论文).长春:东北师范大学.

[41] 孙志豪,谭晓慧,孙志彬,等,2021.基于上限分析的空间变异土质边坡可靠度.岩土力学,42(12):3397-3406.

[42] 田密,盛小涛,2019.岩土工程最小勘探数据量确定方法.岩土力学,40(增刊1):400-408.

[43] 张兰,2019.单目标和多目标量子粒子群优化算法的研究及应用(博士学位论文).西安:西北工业大学.

[44] 仉文岗,孟凡胜,何昌杰,等,2023.岩溶区空间变异性边坡建筑荷载作用下稳定分析.防灾减灾工程学报,43(2):316-323.

基金资助

国家自然科学基金项目(42272326)

国家自然科学基金项目(52222905)

国家自然科学基金项目(521791035)

国家自然科学基金项目(52408371)

江西省自然科学基金项目(20232ACB204031)

江西省自然科学基金项目(20242BAB24001)

江西省自然科学基金项目(S20255146)