0 引言

随着城市化不断推进，地下空间的建设趋势逐渐向大型化、深层化发展，随之产生了一些基坑失稳及坍塌事故。如何准确、合理地确定岩土体力学参数对于基坑动态支护设计和稳定性评价具有重要意义。由于岩土试样在取样、运输和试验过程中会受到不同程度的扰动，所以很多时候室内土工试验法得到的参数并不能完全反映现场岩土体的全部真实状态^[1]。而运用反分析方法可以较好地解决此问题，该方法也逐渐成为获取岩土力学参数的主流方法之一。由于位移容易获取，且精确可靠，所以将数值分析、智能优化算法和现场实测位移相结合是一种有效反演求取出岩土力学参数的方法^[2-3]。

近些年，许多学者在岩土力学参数的反分析方面进行了大量研究。例如文献[4-7]将改进粒子群算法、神经网络算法等引入岩土工程反分析领域，对深基坑黏弹性动态施工、深厚软土参数等问题进行了较为有效的反演分析。虽然这些智能算法取得了较好的效果，但也存在一定的局限性，比如神经网络有网络结构难以确定、依赖大样本等缺点。而极限学习机^[8](Extreme Learning Machine，ELM)作为前馈神经网络学习中一种全新的训练框架，相较于BP神经网络，ELM的学习速度和泛化能力明显提升。叶智峰^[9]利用ELM对区域地应力场进行反分析，结果表明该模型具有良好的泛化性和反演精度。张阳茁^[10]将ELM用于隧洞岩体蠕变参数反分析，建立了岩体蠕变参数与位移之间的非线性映射关系，反演参数的正算位移结果与实测位移值较为接近，验证了此方法的可行性。王述红^[11]利用模拟退火算法(SA)对遗传算法(GA)进行优化，组成GASA-ELM模型，并应用于隧道围岩力学参数反演中，结果表明其模型反演所得参数在ANSYS中正算结果与现场实测值相对误差较小。目前该方法在边坡和隧道领域的参数反演中应用比较广泛，但在深基坑方面的参数反演研究相对较少，而且上述ELM算法存在初始权值和阈值难以确定、模型泛化能力不足等问题。

为此，本研究提出了一种基于改进灰狼算法(Improved Grey Wolf Optimizer，IGWO)结合核极限学习机(Kernel-based Extreme Learning Machine，KELM)的位移反分析方法。首先利用均匀试验以基坑土层参数取值为因素、不同参数取值组合为水平数，计算不同水平下的支护结构深层水平位移值，由此组成学习样本，再使用IGWO算法对核极限学习机的正则化参数C和核参数γ进行优化，构建反演土层等效力学参数的IGWO-KELM模型，最后联合反演出符合实际情况的土层参数取值组合。将IGWO-KELM方法应用于昆明地区某基坑工程中，对基坑土层参数进行反演分析，验证了该方法的可行性。

1 基本理论 1.1 核极限学习机

ELM^[12]是一种新型单隐含层前馈神经网络算法。为了进一步增强ELM算法的泛化能力和稳定性，Huang^[13]于2010年引入核函数取代ELM中的随机映射，提出了一种基于核方法的KELM算法。具体步骤如下：

首先给定N个任意不同的训练样本，其数学模型表示如下：

$ F(\boldsymbol{x})=h(\boldsymbol{x}) \beta=\boldsymbol{H} \beta=\boldsymbol{L}, $

(1)

式中，$\boldsymbol{x}$ 为期望输入向量；$h(\boldsymbol{x})=\boldsymbol{H}$ 被称为随机特征映射矩阵；$\beta$ 为输出权重；$\boldsymbol{L}$ 为期望输出向量。

将式(1)转化为求解线性系统$\boldsymbol{H} \beta=\boldsymbol{L}$ 的最小二乘解，则输出权重$\beta$ 的公式为：

$ \beta=\boldsymbol{H}^* \boldsymbol{L}=\boldsymbol{H}^{\mathrm{T}}\left(\boldsymbol{H} \boldsymbol{H}^{\mathrm{T}}+\frac{\boldsymbol{I}}{C}\right)^{-1} \boldsymbol{L}, $

(2)

式中，$\boldsymbol{H}^*$ 表示隐含层输出矩阵$\boldsymbol{H}$ 的Moore－Penrose广义逆矩阵；$C$ 表示正则化系数；$\boldsymbol{I}$ 表示单位矩阵。

利用Mercer’s条件来定义核矩阵, 并用核矩阵$\boldsymbol{\varOmega}_{\text {ELM }}$ 替代ELM中的随机矩阵$\boldsymbol{H} \boldsymbol{H}^{\mathrm{T}}$，KELM核矩阵定义如下：

$ \left\{\begin{array}{l} \boldsymbol{\varOmega}_{\mathrm{ELM}(i, j)}=h\left(\boldsymbol{x}_i\right) h\left(\boldsymbol{x}_j\right)=K\left(\boldsymbol{x}_i, \boldsymbol{x}_j\right) \\ \boldsymbol{\varOmega}_{\mathrm{ELM}(i, j)}=\boldsymbol{H H}^{\mathrm{T}} \end{array}, \right. $

(3)

式中，$\boldsymbol{x}_i, \boldsymbol{x}_j$为试验输入向量；$K\left(\boldsymbol{x}_i, \boldsymbol{x}_j\right)$ 是高斯核函数。

将式(3)代入式(2)得到KELM模型的预测函数，如式(4)所示：

$ \boldsymbol{F}(x)=\left[\begin{array}{c} K\left(\boldsymbol{x}, \boldsymbol{x}_1\right) \\ \vdots \\ K\left(\boldsymbol{x}, \boldsymbol{x}_N\right) \end{array}\right]^{\mathrm{T}}\left(\boldsymbol{\varOmega}_{\mathrm{ELM}}+\frac{\boldsymbol{I}}{C}\right)^{-1} \boldsymbol{L}。$

(4)

从上述推导过程可以看出，KELM用稳定的核映射代替了ELM的随机初始化的问题，由此来增强模型的稳定性和泛化能力。文献[14]总结发现核参数C和γ对KELM的算法性能影响较大，通过优化C和γ参数能够有效地改善预测结果。因此，将引入IGWO算法来优化C和γ值，以提高模型的泛化性能。

1.2 灰狼算法及改进

灰狼算法是由Mirjalili^[15]受到了自然界中灰狼的等级制度与狩猎活动的启发而提出的一种新型智能优化算法。在该算法中，灰狼种群按照社会等级高低分为：α，β，δ和ω狼。定义适应度最优个体为α狼，次优个体为β狼，第3优个体为δ狼，其余为底层ω狼，即α，β，δ狼均为1只，且狼群总个体数为N。

在基础的灰狼算法中，搜索过程是由每次迭代中3只适应度最佳α，β，δ狼引领至目标解的搜索空间，这容易使狼群过早聚集，造成算法全局和局部搜索不平衡，跳出局部解的能力不足等缺点。为了克服这种缺陷，通过采用一种新的运动搜索策略，即基于维度学习的狩猎(Dimension Learning-Based Hunting，DLH)搜索策略方法^[16]，来改进灰狼算法。

IGWO算法包括3个阶段：信息初始化、信息传播、选择和更新，如下所示。

(1) 信息初始化阶段：将N只灰狼随机分布在区间(l_i，u_j)中，如式(5)所示：

$ \begin{gathered} X_{i j}=l_j+\operatorname{rand}_j[0, 1] \times\left(u_j-l_j\right), \\ i \in[1, N], j \in[1, D], \end{gathered} $

(5)

式中，第$t$ 次迭代中第$i$ 只狼的位置向量$\boldsymbol{X}_i(t)= \left\{x_{i 1}, x_{i 2}, \cdots, x_{i D}\right\}; ~ D$ 表示问题维度；$N$ 为灰狼矩阵的行数。

(2) 信息传播阶段：在IGWO算法结果中由2位候选解组成，第一个候选者X_i－GWO(t+1) 通过基础GWO搜索计算出，另一个候选者X_i－DLH(t+1)采用DLH搜索策略方法计算出。

在DLH搜索策略中，首先要用当前位置X_i(t)和候选位置X_i－DLH(t+1)计算出相对距离R_i(t)，如式(6)所示：

$ R_i(t)=\left\|X_i(t)-X_{i-\text { GWO }}(t+1)\right\| 。$

(6)

再以R_i(t)为半径，构造出个体灰狼X_i(t)向量的邻域N_i(t)，如式(7)所示：

$ \begin{gathered} N_i(t)=\left\{X_j(t) \mid D_i\left(X_i(t), X_j(t)\right) \leqslant R_i(t), \right. \\ \left.X_j(t) \in { pop }\right\}, \end{gathered} $

(7)

式中，$D_i$ 表示$X_i(t)$ 与$X_j(t)$ 的欧氏距离，$p o p$ 代表整个狼群。构建了个体的信息邻域以后，再通过多重邻域学习法，从$N_i(t)$ 内随机选取灰狼位置$X_{\{n, d\}}(t)$ 和狼群矩阵pop内随机灰狼位置$X_{\{r, d\}}(t)$来计算候补位置$X_{\{i-\mathrm{DLH}, d\}}(t+1)$，基于维度信息共享的狼群候选解公式如下：

$ X_{i-\mathrm{DLH}, d}(t+1)=X_{i, d}(t)+\operatorname{rand}\left[X_{n, d}(t)-X_{r, d}(t)\right] 。$

(8)

(3) 选择和更新阶段：通过对狼群进行迭代搜索，比较2个不同候选者的适应度值，选出较优个体，如式(9)所示：

$ X_i(t+1)=\left\{\begin{array}{cc} X_{i-\mathrm{GWO}}(t+1), & f\left|X_{i-\mathrm{GWO}}\right|<f\left|X_{i-\mathrm{DLH}}\right| \\ X_{i-\mathrm{DLH}}(t+1), & \text { 其他 } \end{array}\right. $

(9)

1.3 均匀试验

均匀试验^[17]是由方开泰与王元教授在1978年提出的一种试验设计方案，其设计的试验方法能够以较少的试验点获得最多的试验信息。对于有k个因素X₁，X₂, …，X_k的试验，其均匀试验的设计步骤如下：

(1) 选择每个试验因素的取值范围：[X_imin, X_imax](i=1，2, …，k)，其中X_imin，X_imax分别为第i个因素的最小值和最大值。

(2) 确定因素的试验水平，平分为n个水平，即：

$ X_{i j}=X_{i \min }+\frac{j-1}{n-1}\left(X_{i \max }-X_{i \min }\right), $

(10)

式中，i=1，2，…，k，为因素序号；j=1，2，…，n，为水平序号。

(3) 根据上述试验选择合适的均匀设计表Un (r^k)，其中U表示均匀表代号；n表示试验次数；r表示因素水平数，一般与n相等；k表示均匀表纵列数。

2 基于IGWO-KELM模型的反演流程

由于岩土体是一种非连续性、多相性、变异性的复合材料，具有不确定性特征，传统的显示函数表达式无法有效地表示支护结构位移与岩土力学参数之间的非线性关系。所以利用IGWO算法优化KELM模型建立土层抗剪强度参数与支护结构深层水平位移间的响应关系。根据文献[18]的方法，计算出工况基坑坑顶以下某深度处桩身水平位移实测值的平均值$D_{i=h}$ 和标准差$\sigma$，确定不同深度处实测值的反演区间$\left[D_{i=h}-\sigma, D_{i=h}+\sigma\right]$，在区间内独立随机抽样得到待反演的实测位移值。将智能算法的预测位移值与实测位移值的误差平方和作为目标函数，当目标函数被优化到最小值时，即可得到反演参数最优值。反演目标函数可表示为：

$ F\left(x_{i}\right)=\sum\limits_{i=1}^{n}\left[y\left(x_{i}\right)-y_{i}\right]^{2}, $

(11)

式中，x_i为待反演参数，y (x_i)和y_i分别为某一工况下第i个点的预测位移值和实测位移值；n为该工况的测点数量。参数反分析主要步骤如下所示。

步骤1：根据工程勘测资料，估计待反演的抗剪强度参数的初始取值范围。

步骤2：构造学习样本。将均匀试验设计的参数组合水平作为输入样本，根据《建筑基坑支护技术规程》(JGJ120—2012)计算的深层水平位移作为输出样本。

步骤3：初始化种群数量pop、最大迭代次数iter_max、灰狼搜索区间的范围[X_min，X_max]。

步骤4：评价适应度值，计算当前灰狼候选解X_i－GWO(t+1)和X_i－DLH(t+1)，将候选解代入适应度函数计算适应值，从小到大依次排序X_α，X_β，X_δ。

步骤5：判断当前最优适应度值f_best是否小于期望值f_wish(优化目标所需的收敛精度)，若是，则输出最优C和γ；否则，返回步骤3。

步骤6：构建IGWO-KELM模型。

步骤7：用模型的预测位移值与现场实测位移值的误差平方和作为反演目标函数。

步骤8：将随机抽样获得的实测位移值输入公式(11)中，收敛达到最小值时，即反演出最优土层参数值。

3 工程案例

本研究以昆明草海隧道深基坑工程为例，选取开挖全线中AM21剖面进行分析。该剖面采用SMW工法和“钻孔灌注桩+旋喷桩止水帷幕”2种支护形式，基坑深13.42 m，基坑安全等级为一级。基坑围护结构采用Φ850@600的SMW工法桩，其中支护桩桩长为29.9 m，型钢类型为H700×300×13×24，型钢长度为27 m，采用隔一插二的形式进行布置，顶部采用1 200×800 mm的冠梁连接。依据《建筑基坑支护技术规程》(JGJ120—2012)的计算规则和工程地质勘察报告的相关数据设置合理的施工工况，基坑施工工况如表 1所示。图 2给出了具体基坑支护布置和地层分布情况。基坑土层物理力学指标和抗剪强度取值范围分别见表 2和表 3。

本研究所选的深层水平位移监测点为AM21-CX5，基坑监测从2013年6月9日至2013年7月10日，共计32天。其中施工期间工况3、工况5、工况7相对稳定的监测结果，如图 3~5所示。在MATLAB软件中运用方幂法构建均匀设计表，对上述因素各取40个水平，选用U40(4018)均匀设计表安排试验。根据《建筑基坑支护技术规程》(JGJ120—2012)计算出不同试验水平下深层水平位移值，选择试验水平参数组合作为输入样本，以工况3，5，7地面以下的3 m(测点1)、5 m(测点2)、7 m(测点3)、9 m(测点4)、11 m(测点5)、13 m(测点6)、15 m(测点7)17 m(测点8)处深层水平位移值作为输出样本。选取前30个样本作为模型训练数据，余下的10个样本作为模型测试数据。由于篇幅有限，只列出了工况3，5，7部分样本输入、输出值，见表 4~6。本算法在MATLAB软件中进行编程运行，通过联合多个位移输出点进行反演计算，即先分工况建立每个输出点的KELM预测函数，再以多个位移预测函数和多个实测位移值的误差平方和编程反演目标函数。

4 预测模型评价 4.1 评价指标

试验中，选择平均百分比误差和均方误差根作为模型误差评价指标。

$ e_{\mathrm{MAP}}=\frac{1}{n} \sum\limits_{i=1}^{n} \frac{\left|y_{i}-\hat{y}_{i}\right|}{\left|\hat{y}_{i}\right|} \times 100 \%, $

(12)

$ e_{\mathrm{RMS}}=\sqrt{\frac{1}{n} \sum\limits_{i=1}^{n}\left(y_{i}-\hat{y}_{i}\right)^{2}}, $

(13)

式中，e_MAP为平均百分比误差；e_RMS为均方误差根；n为训练数据数目；y_i为预测数据; $\hat{y}_{i}$ 为真实数据。

为对比IGWO-KELM和GWO-KELM的预测性能，选取工况5测点2的样本数据进行预测分析，先划分训练集和测试集，以训练样本的均方误差根作为优化核参数的目标函数，代入最优核参数建立响应面模型，最后计算真实数据和预测数据的误差指标。平均适应度与最优适应度越接近说明该模型计算过程越稳定。

4.2 核参数寻优

将GWO和IGWO算法的种群数量pop和迭代次数T_max分别设置为30和100，KELM的核参数范围设置为：C∈ [1, 400]，γ∈ [1, 100]。分别利用IGWO、GWO算法耦合KELM核极限学习机，以训练样本的均方根误差作为适应度函数对模型内核参数进行寻优，其GWO和IGWO算法分别对KELM模型核参数C和γ寻优的适应度曲线和结果见图 6和图 7。通过对比图 6和图 7发现，在相同种群数量和迭代次数下，模型GWO-KELM和IGWO-KELM都能在迭代过程中以较快的速度进行收敛，分别在迭代60次和5次后趋于稳定，但IGWO算法在内核参数搜寻过程中的收敛速度明显优于GWO算法。通过图 7可以看出，模型经过IGWO优化后适应度曲线上下起伏较小，说明该模型的计算过程稳定性较好，其最佳适应度值2.01×10^―3也远小于另一个模型的最佳适应度值(GWO-KELM: 4.50×10^―2)。得到优化后的KELM模型后，代入工况5测点2的输入样本计算出预测位移值，再与测试样本的真实位移值进行比较，预测评价结果见表 7。从表 7可以看出，相比于GWO-KELM模型，IGWO-KELM模型的预测值与实际值比较接近。上述表明，采用经过改进的灰狼算法可以更好更快地找到最佳优化值，对KELM模型的预测精确性和寻优稳定性具有较大提升。

5 土层等效力学参数反演

研究表明IGWO-KELM模型具有较好的预测精度和适应性。因此基于该方法对工况3，5，7分别建立IGWO-KELM模型。以预测函数联合多个测点的实测位移值编写反演目标函数，在对应工况测点的实测位移值区间内$\left[D_{i=h}-\sigma, D_{i=h}+\sigma\right]$重复20次独立随机抽样组合，再次调用IGWO算法反演出解集，最后依据解集极值区间[X_min, X_max]得到等效力学参数区间。工况3开挖面以上测点为3，5，7 m；工况5开挖面以上测点为3，5，7，9，11 m；工况7开挖面以上测点为：3，5，7，9，11，13 m，3个工况全部测点为3，5，7，9，11，13，15，17 m，其部分反演迭代曲线见图 8~11。由于篇幅有限，以下只给出了工况3反演的解集，见表 8和表 9。按照2种取点方式分别整理出3个工况的反演结果，见表 10和表 11。

从图 8~11可以看出，在IGWO优化后，以开挖面以上测点位移反演的适应度值分别为0.68和0.89，结果整体比全部测点的好。从表 8和表 9可以看出，工况3开挖范围内(开挖至7.52 m)的第1~3层土的参数反演极值范围比较稳定，但由第4层③₂粉土开始，开挖范围外的土层参数极值区间出现波动较大的情况。通过分析表 10和表 11，发现在工况5和工况7开挖面以上取点反演时，都表现出开挖面以上的土层(①₁填土、③₁泥炭质土、③_1-1黏土、③₂粉土)结果范围比较收敛，而从开挖面以下的③₃黏土开始，土层参数反演区间浮动较大，与上述工况3分析出的规律基本一致。由此可知该方法对开挖面以上土层的抗剪强度参数反演时，反演结果比较可靠，而对开挖面以下土层的抗剪强度参数反演时，结果不太理想。

分析其主要原因是基坑土方开挖对基坑开挖深度范围内的土体影响比较显著，而对基坑开挖面以下嵌固端部分土体影响较小，导致开挖范围内支护桩位移变形比嵌固端的明显且准确。对比表 10和表 11可知，以开挖面以上测点位移反演的结果范围从整体上要比以全部测点位移反演的结果范围更小，得到的交集区间也相对更加可靠。因此选择以开挖面以上测点位移反演的交集结果作为最终反演结果，见表 12。

通过表 12的对比，发现原始设计参数和最终反演参数结果存在一定差异，原始设计参数取值存在不合理之处。将表 12中的最终土层参数反演分析结果取区间中值，代入同一基坑下不同剖面模型中计算工况3，5，7处基坑坑顶以下测点4，8，12 m处的桩体水平位移，并与对应位置处统计得到位移实测值的均值进行对比，如表 13所示。对比结果表明，在3个不同工况下，深层桩体水平位移测点计算值与对应实测位移均值的绝对误差较小，验证了该算法的可行性。

6 结论

(1) 本研究基于均匀试验设计得到学习样本，再利用一种信息共享搜索策略改进的GWO算法对KELM模型进行核参数优化，得到最优核参数C和γ，再用优化后的模型进行预测计算。对比GWO-KELM模型，IGWO-KELM模型的2种误差指标RMSE和MAPE分别降低了37.89%和33.33%，可知IGWO算法能够有效提高KELM模型的预测精度。

(2) 基于IGWO-KELM的深基坑土层等效力学参数反演分析方法可以较好地反演分析动态施工的土层抗剪强度参数。但该方法对工况开挖面以上土层的抗剪强度参数反演时结果比较准确，对工况开挖面以下土层的抗剪强度参数反演时结果范围浮动较大，对于基坑工况开挖面以下土层的抗剪强度参数反演方法还需进一步研究。

(3) 由于观测仪器、外界条件等诸多因素的影响，监测数据中出现系统误差和粗差的可能性较大。本研究在监测数据点的容差区间里进行随机抽样反演计算，并验证其正算结果与实测位移值比较接近，可知采用多点联合的方法对含有误差的监测数据反演时具有一定可行性。