基于DE-LSSVM的地下洞室群支护参数智能优化

2019-09-17

长江科学院院报 2019年9期

(1.四川省交通运输厅公路规划勘察设计研究院，成都 610041；2.国电大渡河流域水电开发有限公司，成都 610016)

1 研究背景

水电能源是目前可再生和非化石能源中资源最明确、技术最成熟、最清洁和最经济的，也是世界公认的清洁能源。随着我国在降低二氧化碳排放方面面临的压力和责任越来越大，水电对我国实现低碳经济的作用和效果将更加显现。在西部高山峡谷地区修建的大型水电站中，由于受地形地质条件的限制，不可避免要在地下岩石中修建规模宏大的地下厂房洞室群。在如此庞大的地下工程中，如何用最小的开挖支护代价来保证地下空间在施工和运行期有足够的稳定性和安全性，也成为一个重要的研究课题[1-2]。

目前地下洞室群的初步支护方案主要是根据围岩分级并结合地质结构进行设计的，支护参数优化可以在此基础上采取多方案比选和智能优化方法进行[3]。多方案比选由于其选择范围有限，很难保证能在一定的全局空间中得到最优的支护参数。因此，较多的学者将神经网络和遗传算法等智能方法应用在地下洞室的支护参数优化上，建立了智能支护优化方法[4-5]。该方法通过神经网络等模型建立支护参数与评价指标之间的非线性映射关系，然后利用遗传算法、粒子群算法等优化方法以综合经济指标为目标函数，从全局空间中选取最佳的支护参数。随着统计学习理论的不断发展，支持向量机模型也逐步被应用在智能预测和模式识别等领域[6]。差异进化算法由于出众的全局优化能力也被用来搜索最优的神经网络、支持向量机等复杂模型的参数，以建立最优的预测模型[7-8]。因此，本文在已有研究成果的基础上，采用近年来迅速发展的最小二乘支持向量机模型和差异进化算法来建立智能支护优化方法，并将该方法应用在西部某大型水电站地下洞室群支护参数的优化分析中。

2 支护参数智能优化方法

2.1 差异进化算法

差异进化算法(Differential Evolution ,DE)是1995年Kenneth Price和Rainer Storn提出的一种更有效的新型优化算法，是进化算法产生以来取得的又一巨大进展。1996年该算法参加了首届进化算法大赛，在所有参赛优化算法中，DE被证明为计算速度最快的进化算法[9]，随后该算法在多个领域得到了广泛应用。由于差异进化算法也采用了“优胜劣汰、适者生存”的自然进化法则，所以也属于演化算法的一种。对于连续变量的函数优化，DE能更快、更稳定地搜索到问题的全局最优解。

差异进化算法的基本原理是：对初始种群样本中的每个个体i,从当前种群中随机选择3个点，以其中1个点为基础点，以另外2个点为参照点作一个扰动，所得点与个体i交叉后进行“自然选择”，保留较优者，实现种群的繁衍进化[10]。

假设待优化的问题为minf(x)，则DE算法步骤主要为：

(2)个体计算评价。计算每个个体Xi(t)的目标值f(Xi(t))。

(3)变异与交叉繁殖。对种群样本中的每个个体Xi(t)，随机生成3个不相同的随机整数r1,r2,r3∈(1,2, …，N)和随机整数jrand∈(1,2, …，n),如果rand[0,1]

(4)选择。如果f(xj(i)′(t))

(5)终止检验。通过迭代次数和预测误差的设定来判断是否符合终止条件，如果种群个体Xi(t+1)满足终止条件，则输出Xi(t+1)中具有最小目标值的个体作为最优解，否则转到步骤(2)。

2.2 LSSVM模型

支持向量机模型SVM是Vapnik等在统计学理论基础上提出的一种新的机器学习方法[11]。最小二乘支持向量机LSSVM是Suykens提出的一种新型支持向量机[12]，能用来解决模式识别和数据回归问题。假设训练样本数据为{xi,yi}，i=1,2,…,k，xi∈Rn为n维系统输入向量，yi∈R为向量输出值。在本文中所应用的回归问题中，可以在高维特征空间中构建拟合数据的最优函数，即

f(x)=wφ(x)+b。

(1)

式中：w为权向量；b为偏置量。在优化目标中选取误差ξi的二范数作为损失目标函数，把标准SVM的不等式约束改变为等式约束，这使得LSSVM优化问题的最小函数为[6]

(2)

约束条件为

yi-wφ(xi)=b+ξi。

(3)

式中C为惩罚因子，是一个调节常数，它能够使所求解的函数具有较好的泛化能力。

根据式(2)、式(3)可以建立Lagrange求解方程，即

(4)

式中：ai(i=1,2,…,k)是Lagrange乘子；最优的参数a和b可以通过KKT条件得到，即

(5)

消去式(5)中的w和ξ，优化问题转化为求解如下方程：

(6)

式中：y=[y1,y2，…，yk]T；a=[a1,a2,…，ak]T,Θ=[1,1，…,1]T；Ω为一个方阵，其第i行第j列元素为Ωij=φ(xi)Tφ(xj)。高维空间的内积计算可以利用原输入空间的一个核函数来进行等效变换，所以在处理非线性问题时可以不必进行非线性变换，直接采用核函数来代替内积计算，即K(xi,xj)=φ(xi)φ(xj)。在求得参数a和b后，可以得到LSSVM的回归模型，即

(7)

核函数K(xi,xj)是满足Mercer条件的任意对称函数，本文选用径向基函数，并采用差异进化算法来优化核函数的参数取值和惩罚因子C的取值。

2.3 支护优化目标函数

目标函数的选取对优化结果有重要影响，其取值直接影响优化结果的合理性。支护参数的优化目标是在一定的约束条件下使地下洞室群稳定性与经济性达到综合最佳[6]。本次优化的地下洞室群稳定性评价指标为主厂房顶拱最大位移、主变室顶拱最大位移、尾调室顶拱最大位移、塑性区体积和拉屈服区体积，经济性评价指标为锚杆支护等效总长度。由于各指标的量级和量纲有一定的差别，因此需对各个指标值进行标准化，即

pi=(xi-xi,min)/(xi,max-xi,min) 。

(8)

式中：xi为第i个指标；pi为标准化后的值。

用加权法将各指标进行综合，得到综合支护方案的综合评价指标为

(9)

式中wi为权系数，一般由专家经验确定。在本文的支护优化中，主厂房顶拱最大位移权系数为0.15，主变室顶拱最大位移权系数为0.1，尾调室最大位移权系数为0.1，塑性区体积权系数为0.15，拉屈服区体积权系数为0.1，锚杆支护等效总长度权系数为0.4，即地下洞室群稳定性指标权重为0.6，经济性指标权重为0.4。式(9)也是本文优化所采用的目标函数(即差异进化算法中的适应性函数)。约束条件为现有锚固支护方案、支护参数的可行性、施工可操作性以及计算力学原理等。

2.4 DE-LSSVM智能支护优化方法

将智能优化(差异进化算法)、计算方法(最小二乘支持向量机)和FLAC3D三维数值计算相结合，通过FLAC3D计算建立学习样本，用最小二乘支持向量机学习并建立支护参数与安全经济性评价指标之间的非线性映射关系，进而利用差异进化算法的全局寻优能力和最小二乘支持向量机的预测功能，以目标函数值作为评价指标，随机产生一组初始参数，对该组方案进行进化操作，产生下一代方案，由此进行下一步操作直至找到最小目标函数值，最小目标函数值对应的方案即为最优方案，从而实现对洞室群围岩支护参数的全局寻优。具体计算流程如图1所示。

图1 DE-LSSVM智能支护优化方法流程Fig.1 Flowchart of intelligent optimization of supporting parameters based on DE-LSSVM

3 工程实例

3.1 地下洞室群三维数值模型

某水电站地下厂房洞室群围岩主要为花岗岩，局部穿插辉绿岩脉。厂区岩体新鲜较完整，呈块状-次块状结构，岩块嵌合紧密。厂房区共安装有4台机组，每个机组的纵向长度约为34 m，所在区段的地质条件不同。如图2所示，利用FLAC3D建立数值计算模型，在分析区域内，考虑了岩脉和断层等地质结构以及全强风化层、弱风化层、微新岩体即Ⅴ类、Ⅳ类、Ⅲ类、Ⅱ类等岩体属性，岩体力学参数取值见表1。计算区域包含了主厂房、引水洞、主变室、母线洞室、尾水调压室和尾水洞。计算区域的4个侧面均采用法向约束，底部采用固定约束。围岩屈服准则采用Mohr-Coulomb模型，锚杆采用FLAC3D中的cable单元来模拟实现。

图2 三维数值计算模型Fig.2 Three-dimensional meshes of numerical calculation model of underground powerhouse

表1 地下洞室群岩体力学参数取值Table 1 Physico-mechanical parameters of rock mass

根据现场地应力测试及反演分析结果形成侧压力系数，以x向侧压力系数kx=0.6(水平上垂直地下厂房轴线方向)、y向侧压力系数ky=1.44(地下厂房轴线方向)为基本分析参数模拟初始地应力场，其中，全、强、弱风化层按自重应力场模拟，微新岩体竖直向z向应力按自重应力场施加，x向与y向应力按侧压力系数进行施加[13]。

3.2 支护参数优化方案设计

支护参数的智能优化首先要构造学习样本用于最小二乘支持向量机的学习，学习样本通过数值仿真计算来获得。大型三维数值计算是相当费时的过程，因此，正算样本方案的构造应在保证典型代表性的前提下尽可能地减少样本数，因此本文采用正交设计的方法来构造样本方案[6]。正交试验设计方法根据正交性原则来选择试验范围内的代表点，若试验包含n个因素，每个因素包含m个水平，则全面试验的试验点个数共为m个，而正交设计仅有m2个，并且具有“整齐可比性”和“均匀分散性”，非常适合于多因素和多水平的试验方案设计[14]。根据正交试验设计的优势，对最小二乘支持向量机的学习样本采用正交试验设计法进行选取，这样既可以保证模型学习、预测的准确性，又可以减少试验的次数，从而节省计算工作量。

根据设计人员建议，洞室群上部开挖的支护优化主要针对系统锚杆的设计参数进行，洞室群上部开挖顺序如图3所示。该地下厂房区围岩总体以Ⅱ类为主，依据围岩等级初步拟定的系统锚杆支护方案为：采取预应力锚杆和普通砂浆锚杆进行梅花型交叉布置，其中预应力锚杆长9 m，普通砂浆锚杆长6 m，布锚间排距为1.5 m。

图3 地下厂房洞室群上部开挖示意图Fig.3 Simplified excavation model of underground powerhouse

根据上述初步设计方案，洞室群上部开挖中优化的支护参数为预应力锚杆长度、普通砂浆锚杆长度、锚杆间排距、锚杆设计预应力。每个参数取4个水平，按正交试验设计共16种方案。依据《锚杆喷射混凝土支护技术规范》(GB 50086—2001)[15]中的规定，并根据初步支护方案来确定参数的取值范围，按正交表进行方案的设计构造，待计算的样本方案如表2所示。

3.3 支护参数优化结果分析

采用FLAC3D对样本方案进行三维数值计算，得到各个支护参数组合下地下厂房洞室群开挖支护的力学响应特征量，包括主厂房顶拱最大位移、主变室顶拱最大位移、尾调室顶拱最大位移、塑性区体积和拉屈服区体积。将每个支护方案所需的锚杆总长度结合锚杆设计预应力得到等效的锚杆支护总长度并作为支护费用的指标即经济性指标。支护样本方案对应的5个力学响应特征量和1个经济指标特征量如表3所示，以此16个计算方案的支护参数组合和响应特征量作为最小二乘支持向量机的学习特征样本的输入值和输出值，用最小二乘支持向量机建立起支护参数与力学响应特征变量和支护费用之间的非线性映射关系。

表2 计算样本方案及对应的支护参数组合
Table 2 Computation scenarios and correspondingsupporting parameters

方案编号预应力锚杆长度/m普通砂浆锚杆长度/m锚杆间排距/m锚杆设计预应力/kN1951.51502961.32003971.42504981.630051051.325061061.530071071.615081081.420091151.4300101161.6250111171.5200121181.3150131251.6200141261.4150151271.3300161281.5250

表3 计算样本方案支护效果Table 3 Supporting effects of scenarios for computation

训练过程中利用差异进化算法搜索到最优的LSSVM参数惩罚因子C=1 000，模型核函数参数σ2=0.59。为检验LSSVM模型的预测精度，以初步设计方案的支护参数作为输入，对已建立的LSSVM模型进行验证。如表4所示，LSSVM模型预测结果与FLAC3D数值计算结果相当吻合，表明LSSVM模型已正确地建立起了支护参数组合与力学响应特征量和支护费用之间的映射关系，LSSVM模型可以代替FLAC3D数值计算过程，从而为差异进化算法进行全局搜索得到最优方案的合理性提供了保证。

表4 LSSVM预测与FLAC3D计算值对比Table 4 Comparison between LSSVM prediction result and FLAC3D calculated result

以式(9)的目标函数为评价指标，对样本支护方案和设计支护方案进行评价。表3的最后一列给出了根据16个样本支护方案和设计支护方案三维数值计算结果按上述方法进行标准化后的综合评价值。其中方案2的综合评价值最小，为0.283，其次是设计方案的综合评价值为0.370。这2个方案是17个方案中相对最优的方案，仍需从全局空间搜索最优的支护参数组合。

在上述工作的基础上，结合工程实际确定各支护参数的搜索区间，通过差异进化算法进行全局寻优，得到最优的方案，支护参数(预应力锚杆长度，普通砂浆锚杆长度，锚杆间排距，锚杆设计预应力)取值分别为：9.12 m，5.91 m，1.32 m，207 kN。经验算将此支护方案参数代入FLAC3D数值计算所得目标函数值也是各方案中最小的，因此对上述方案值进行取整得到的结果即为方案2，即全局最优方案，对应的支护参数为：预应力锚杆长为9 m，普通砂浆锚杆长为6 m，布锚间排距都为1.3 m，锚杆设计预应力为200 kN。