吴忠广等：深埋硬岩隧道围岩参数概率反演方法 81· (0°，r≥u 径向基函数.本文选择性能较好的径向基函数： 0:= (7) (0:-1,r<u /1x-x:12) 4)确定是否收敛：如果不满足要求，重复步骤 K(x,;)=exp2 (14) 1)~3),直至产生稳定序列 核函数参数σ2和惩罚因子C是影响多输出支 对于转移概率分布函数的协方差矩阵取值，参 持向量机建模精度的两个重要因素，为高效确定最 照Zhang等)研究成果，当转移概率分布函数的协 优参数值，采用粒子群搜索的方法优化（σ2，C)的取 方差矩阵为先验概率分布协方差矩阵值的0.5时， 值.粒子群算法[]主要思想是首先初始化一群随 马尔科夫链计算效率较高且能够得到合理的接受 机粒子，然后通过迭代寻找最优解.在每次迭代中， 率.对于模型残差e取值，参照Peng等)]等推荐 粒子通过跟踪两个极值来更新自己.一个是粒子本 值，假设E~N(0,0.25u),其中，u为P。、P1与EDZ 身所找到的最优解，即个体极值P,另一个是整体 三个测试值各自对应的平均值. 种群目前找到的最优解，称为全局极值g·粒子 2.2基于粒子群-多输出支持向量机的智能响应 找到这两个极值后，用下式更新自己的速度和位置： 面方法 V=wV+crand()(PBea -p)+czrand()(gBea-p) 本文采用多输出支持向量机模型建立待反演参数 (15) 与隧道监测数据之间的非线性映射关系，假设有组 P=p+V (16) 训练样本数据{x,}(i=1,2,…,n),x:∈R"为待反演 式中：V是粒子的速度，p为粒子的当前位置，c,和c2 参数，y,∈R为监测数据计算值，可建立如下关系： 是学习因子，通常在0~2间取值，rand()为(0,1) y:R"→R (8) 之间的随机数，为加权因子.针对传统粒子群算 根据多输出支持向量机理论[)，相应的映射模 法易早熟以及算法后期易在全局最优解附近产生振 型与对应的优化问题分别见式(9)与式(10)： 荡的现象，采用权重线性递减的粒子群算法[s0]解决 (x)=W(x)+b (9) 此问题.权重随算法迭代次数的变化公式为： IwI2+C()(10) 0=0m-(0s-0ea) (17) 式中：u:=le:‖=√ee,e=y-(x)rW-b, 式中：w和wn分别表示0的最大值和最小值，t W=[w2,w2,…,w]:b=[b,b2,…,b9] 为当前的迭代步数，tm为最大迭代步数，通常取 0mx=0.9,0min=0.4. 为了得到式(10)对应的最优解，采用迭代加权 在优化过程中，目标函数选用通用的均方差函 最小二乘算法(RWLS)求解2-]，进而得到下式： 数，建立的粒子群-多输出支持向量机模型作为代 L,(W,b)=g(W,b)= 替有限元方法的智能响应面用于MCMC模拟中估 a,4+T (11) 算似然函数 式中：T是与W或b无关的常数项的和； 2.3模型计算过程 概率反演B-PSO-MSVM方法具体计算流程见 0, u:≤B 下图1. a:=2C(u-E) (12) u>8 具体流程如下： 步骤一：通过室内试验与现场测量统计，获取 根据L”(W,b)相对w与驻点条件转化得到 UCS、CI/UCS和T的先验分布信息，如均值、标准 下式： 差、变异系数、分布类型等： K+D j=1,2,…,Q) 步骤二：利用均匀设计建立参数样本集，结合 Phase22数值模拟计算隧道变形及损伤区深度值，构 (13) 建多输出支持向量机学习样本； 式中：D。=diag(a1,a2,…,a1）;=[a,a2,…, 步骤三：利用粒子群方法优化模型参数，通过式 a1]'y=[yny2,…,y];p,=[B,f2,…,B]. (9)建立多输出支持向量机智能响应面模型代替数 其中，K=是核函数矩阵，Φ=[p(x1), 值模型； p(x2),…,p(x)]. 步骤四：将Bayesian方法与多输出支持向量机 常见的核函数分别有线性函数、多项式函数与 模型相结合，基于实际监测数据，构建B-粒子群-吴忠广等: 深埋硬岩隧道围岩参数概率反演方法 兹i = 兹 * , r逸u 兹{ i - 1 , r < u (7) 4)确定是否收敛;如果不满足要求,重复步骤 1) ~ 3),直至产生稳定序列. 对于转移概率分布函数的协方差矩阵取值,参 照 Zhang 等[7]研究成果,当转移概率分布函数的协 方差矩阵为先验概率分布协方差矩阵值的 0郾 5 时, 马尔科夫链计算效率较高且能够得到合理的接受 率. 对于模型残差 着 取值,参照 Peng 等[17] 等推荐 值,假设 着 ~ N(0, 0郾 25滋),其中,滋 为 P0 、P1与 EDZ 三个测试值各自对应的平均值. 2郾 2 基于粒子群鄄鄄 多输出支持向量机的智能响应 面方法 本文采用多输出支持向量机模型建立待反演参数 与隧道监测数据之间的非线性映射关系,假设有 n 组 训练样本数据{xi,y^ i}(i =1,2,…,n),xi沂R n 为待反演 参数,yi沂R 为监测数据计算值,可建立如下关系: y^:R n寅R (8) 根据多输出支持向量机理论[32] ,相应的映射模 型与对应的优化问题分别见式(9)与式(10): y^(x) = W T渍(x) + b (9) Lp (W,b) = 1 2 移 Q j = 1 椰w j椰2 + C 移 l i = 1 L(ui) (10) 式中:ui = 椰ei椰 = e T i ei,e T i = y T i - 渍( xi ) TW - b T , W = [w 1 ,w 2 ,…,w Q ];b = [b 1 ,b 2 ,…,b Q ] T . 为了得到式(10)对应的最优解,采用迭代加权 最小二乘算法(IRWLS)求解[32鄄鄄33] ,进而得到下式: Lp (W,b) = L义p (W,b) = 1 2 移 Q j = 1 椰w j椰2 + 1 2 移 l i = 1 aiu 2 i + 子 (11) 式中:子 是与 W 或 b 无关的常数项的和; ai = 0, ui臆着 2C(ui - 着) ui , ui > ì î í ïï ïï 着 (12) 根据 L义p (W,b)相对 w j 与 b j 驻点条件转化得到 下式: K + D - 1 琢 1 琢 T 1 é ë ê ê ù û ú T ú 茁 j b é ë ê ê ù û ú j ú = y j 琢 T y é ë ê ê ù û ú j ú (j = 1,2,…,Q) (13) 式中:D琢 = diag ( 琢1 , 琢2 ,…, 琢1 ); 琢 = [ 琢1 , 琢2 , …, 琢1 ] T ;yj = [yj1 ,yj2 ,…,yj1 ] T ;茁j = [茁j1 ,茁j2 ,…,茁j1 ]. 其中, K = 椎椎 T 是核函数矩阵,椎 = [渍( x1 ), 渍(x2 ),…, 渍(xl)] T . 常见的核函数分别有线性函数、多项式函数与 径向基函数. 本文选择性能较好的径向基函数: K(x,xi) = exp { - | x - xi | 2 滓 2 } (14) 核函数参数 滓 2和惩罚因子 C 是影响多输出支 持向量机建模精度的两个重要因素,为高效确定最 优参数值,采用粒子群搜索的方法优化(滓 2 ,C)的取 值. 粒子群算法[35] 主要思想是首先初始化一群随 机粒子,然后通过迭代寻找最优解. 在每次迭代中, 粒子通过跟踪两个极值来更新自己. 一个是粒子本 身所找到的最优解,即个体极值 pBest,另一个是整体 种群目前找到的最优解,称为全局极值 gBest . 粒子 找到这两个极值后,用下式更新自己的速度和位置: V = wV + c1 rand()(pBest - p) + c2 rand()(gBest - p) (15) p = p + V (16) 式中:V 是粒子的速度,p 为粒子的当前位置,c1和 c2 是学习因子,通常在 0 ~ 2 间取值,rand( )为 (0,1) 之间的随机数,w 为加权因子. 针对传统粒子群算 法易早熟以及算法后期易在全局最优解附近产生振 荡的现象,采用权重线性递减的粒子群算法[50] 解决 此问题. 权重随算法迭代次数的变化公式为: w = wmax - t(wmax - wmin ) tmax (17) 式中:wmax 和 wmin分别表示 w 的最大值和最小值,t 为当前的迭代步数, tmax 为最大迭代步数,通常取 wmax = 0郾 9,wmin = 0郾 4. 在优化过程中,目标函数选用通用的均方差函 数,建立的粒子群鄄鄄 多输出支持向量机模型作为代 替有限元方法的智能响应面用于 MCMC 模拟中估 算似然函数. 2郾 3 模型计算过程 概率反演 B鄄PSO鄄MSVM 方法具体计算流程见 下图 1. 具体流程如下: 步骤一:通过室内试验与现场测量统计,获取 UCS、CI/ UCS 和 T 的先验分布信息,如均值、标准 差、变异系数、分布类型等; 步骤二:利用均匀设计建立参数样本集,结合 Phase2 数值模拟计算隧道变形及损伤区深度值,构 建多输出支持向量机学习样本; 步骤三:利用粒子群方法优化模型参数,通过式 (9)建立多输出支持向量机智能响应面模型代替数 值模型; 步骤四:将 Bayesian 方法与多输出支持向量机 模型相结合,基于实际监测数据,构建 B鄄鄄 粒子群鄄鄄 ·81·