基于二次正交试验优化的隧道围岩稳定可靠度方法.pdf

第 34 卷 第 2 期 岩 土 工 程 学 报 Vol.34 No.2 2012 年 .2 月 Chinese Journal of Geotechnical Engineering Feb. 2012 基于二次正交试验优化的隧道围岩稳定可靠度方法 苏永华 1，李 翔1，丁 云1，孙晓明2 （1. 湖南大学土木工程学院，湖南 长沙 410082；2. 深部岩土力学与地下工程国家重点实验室（北京） ，北京 100083） 摘 要通过把二次回归正交组合试验设计和显著性检验结合起来，构建出了各类因素对结构特征值影响程度的判别 方法。将该方法植入隧道围岩稳定响应面可靠性分析程序中，建立了具有识别基本因素重要性程度功能的状态函数优 化方法，消除了序列响应面状态函数建立时不能判断基本因素对结构状态重要性程度的缺陷。利用 MATAB 搜索技术 求解基于状态函数优化方法改进的响应面可靠度指标，解决了分析结果失真问题。上述两者的结合形成了具有较高效 率的复杂地层隧道围岩稳定可靠度响应面分析技术。通过某工程实例研究，演示了该技术的运行过程和操作程序，并 在效率、分析结论的准确性等方面进行了对比验证。 关键词隧道工程；可靠度；二次正交组合试验；响应面；回归 中图分类号TU43 文献标识码A 文章编号1000–4548201202–0326–07 作者简介苏永华1966– ，男，湖南涟源人，博士（后），教授，博士生导师，主要从事岩石力学和地下结构非确 定性及优化设计方面的教学与研究工作。E-mail syh5327。 Reliability degree for stability of surrounding rock of tunnels based on quadratic orthogonal experimental optimization SU Yong-hua1, LI Xiang1, DING Yun1, SUN Xiao-ming2 1. College of Civil Engineering, Hunan University, Changsha 410082, China; 2. State Key Laboratory for GeoMechanics and Deep Underground Engineering, Beijing 100083, China Abstract By integrating the quadratic regression analysis of orthogonal composite design and the statistical significance testing, a discriminant to describe the effect degree of various factors on structural eigenvalue is established. By introducing such a into the reliability analysis process of response surface for the stability of the surrounding rock of tunnels, an optimum approach for identifying the degree of importance with respect to basic ingredients for state function is developed, which can eliminate the drawback of being incapable of assessing the importance degree of basic ingredients on the structural state when establishing a series of response surfaces for the state function. A MATLAB-based search technology is utilized to solve the reliability index of the improved response surface by means of the optimum approach developed for the state function, which figures out the problem of result distortion. On the basis of the integration between both the above s, a response surface analysis technique with relatively high computational efficiency is then proposed for the reliability uation of the surrounding rock of tunnels in complicated stratum. The proposed technique is illustrated by a certain practical example, in which the running process and operating procedure are demonstrated, and then several aspects related to this technique such as the computational efficiency and accuracy of the results are also verified, showing its validity and applicability. Key words tunnel engineering; reliability; quadratic orthogonal composite design; response surface; regression 0 引 言 在地层赋存条件复杂的隧道及边坡等地质体工程 结构的可靠度分析中，由于介质的复杂性和目前数学 物理力学水平的制约，在许多情况下直接、清晰而简 明地建立起结构的解析状态方程存在无法克服的困 难。因此，随机有限元法[1-4]和基于确定性有限元的响 应面方法[5-7]就成为岩土工程稳定可靠度分析中最常 见的选择。其中基于确定性有限元的响应面方法，由 于模拟次数大大少于蒙特卡洛法，同时求解的主要工 具是确定性有限元数值试验，代价低、过程易于掌握 而吸引到更多的注意力。特别是在地下工程的围岩稳 定可靠度计算中，响应面方法正在获得越来越广泛的 ─────── 基金项目国家自然科学基金项目（51078136） ；深部岩土力学与地下 工程国家重点实验室开放基金项目（SKLGDUEK0915） ；湖南省自然 科学基金项目（09JJ3113） 收稿日期2010–11–22 第 2 期 苏永华，等. 基于二次正交试验优化的隧道围岩稳定可靠度方法 327 应用[8-10]。 在应用过程中，响应面理论也不断得到充实和更 新。如佟晓利等[11]将可靠度几何指标法引入，提高了 响应面的效率和精度；宋玉香等[9]提出了二水平试验 选点方法，拓宽了响应面在隧道结构中的应用。徐军 等[12]提出了基于有理多项式函数的响应面拟合方法； 苏永华等[13]提出了完整的基于区间变量的响应面非 概率可靠性分析方法。这些研究的重点主要在于状态 方程的拟合形式及计算精度上，无疑是非常正确的， 也有力地推动了地下工程围岩稳定可靠度方法的进 步。 本文将在前人关于响应面研究成就的基础上，针 对应用广泛的序列响应面方法中关于状态方程的形式 和可靠度指标的确定问题提出改进方法，使其在解决 复杂地下工程围岩稳定可靠度时，理论上更为严谨， 计算结果上更为有效。 1 序列响应面法隐含的前提 设隧道地层围岩与衬砌结构等物理力学表征指标 为一随机向量，表示为 XX1，X2，，Xn，结构状 态的特征值为。在目前响应面方法中，通常选用不 带交叉项的二次多项式近似表达结构特征值与基本参 数的关系[8-9]，即 fXfX1,,Xn≈hX 2 0 11 nn jjjjj jj aa Xa X    。 1 式中 n 为基本参数的个数；a0，aj， ajj为待定系数。 如果采用围岩和衬砌结构收敛变形为判断依据， 式（1）中结构特征值可取为拱顶下沉或周边收敛。 根据经验或参照规范，确定在某类岩体中的允许值 为 umax，则可建立状态方程 gXumax-hX为如下形 式 gXumax-hXumax－ 2 0 11 nn jjjjj jj aa Xa X    2 0 11 nn jjjjj jj bb Xb X    。 2 式中 j1，2，，n，n 为基本参数的数目；b0，bj， bjj为待定系数。 以式 （2） 为结构特征状态方程基本形式的可靠度 计算，目前的求解主要有最小二乘拟合法和序列迭代 法两种方案，其中序列迭代法无论在计算过程的清晰 度还是精确性方面都占有较大优势，因此在地下工程 中采用最为广泛[8, 10]。 但仔细分析序列响应面方法的思路和算法（详见 文献[14]） ，发现它必须满足两个前提首先，结构的 状态特征值为全体基本参数的函数，各基本参数对特 征值具有同样的影响程度；其次，反复迭代收敛到的 可靠度指标最小值即为所研究对象的可靠度。 很显然，这两个前提无论在理论上还是实际工程 中， 都不能够无条件地满足。 就复杂的地下结构而言， 可能相应地导致两方面的后果[15]首先是状态方程过 于庞大、复杂，计算工作量大，效率低，难于满足工 程要求；其次是得到的最小可靠度指标是全局解，而 具体实际工程可能是其中的局部解， 使分析结果失真。 因此，以解决上述问题为目标，本文首先将二次 回归正交组合试验设计和显著性检验原理结合起来， 研究基本因素对结构特征值重要程度的判别方法；然 后将该方法植入响应面状态方程构建程序中，形成对 基本因素具有过滤功能的复杂结构状态函数优化方 法；最后在可靠度指标的计算过程中，引入基于变量 范围约束的 MATLAB 搜索优化方法，提高计算效率 并使分析结果符合具体工程实际。 2 二次正交回归组合试验设计方法 在序列响应面方法中，基于计算简便，采用不带 交叉项的二次式（2）作为状态方程基本形式。本文基 于拟合精度和表达的准确性，状态方程基本形式采用 完全的二次式形式，即 2 0 11 nn jjkjkjjjj jkjj Q Xcc Xc X Xb X    。 3 式中 XX1，X2，，Xn为隧道围岩或衬砌结构基 本物理力学参数；n 为参数数目；c0，cj，ckj和 cjj为待 定系数。 2.1 二次回归正交试验水平 二次回归正交试验设计方法系正交试验设计方法 发展而来。设参数 xj的上、下边界分别为 xjup，xjlower， 算术平均值为 xj0xjup＋xjlower/2。 状态方程（3）作为二次回归形式而言，共有 n1n2/2 个待定系数，至 少需要进 行n1 n2/2 次试验。 若按照正交试验设计方法采用二水平试验抽样 点，如果全面实施则抽得试验样本点数 mc2n，导致 试验次数过多；如果 1/4 实施则 mc2n -2，试验次数远 不够；作为折中，采用 1/2 实施，即 mc2n -1，二水平 1/2 实施的试验点满足正交条件。 但 mc2n -1 次试验，样本数目还是不够，需要增 加试验样本点。零水平是一个重要的水平，所以增加 零水平的试验次数 m0。 另外，就式（3）的回归相关性而言，需要考虑基 本变量的个性影响和边界值。参照序列响应面[14]的抽 样思路，将各个基本变量边界点增加为试验点，将这 328 岩 土 工 程 学 报 2012 年 些点标记为星号“*” ，称为星号试验点，图 1 是二元 二次抽样点布置示意图。设星号试验点与中心点的距 离为 r，称为星号臂长，星号试验次数为 mr2n，总 试验次数为 mmcmrm0。 图 1 二元二次抽样点布置示意图 Fig. 1 Sketch for sampled points of dualistic variables 最后，试验表中的试验将由二水平试验、零水平 试验和星号试验三个类型的试验点组合而成。相应的 试验水平有–r，–1，0，1，r 共 5 个水平。 2.2 试验样本点及二次正交表 根据正交性条件，可以得到星号臂长度 r 的计算 式[16]为 cr0cc 2 mmm mm r   。 4 星号点作为边界点，xjr，x-jr实际上分别是 r 和–r 水平点，零水平点的相应值 xj0为 xj0[xjup＋xjlower/2]xjr＋x-jr/2 。 5 根据图 1，因素 xj的变化间距j为 jxjr – xj0/r 。 6 按照二水平正交试验理论，试验水平（规范化变 量，也称试验因素水平编码）与自然变量之间存在 一一对应关系。设自然变量 xj相应的规范化变量为 zj，xj与 zj之间的对应变换关系为 zj xj –xj0/j 。 7 由此建立起 zj与 xj之间的一一对应关系见表 1。 表 1 规范变量与自然变量对应关系 Table 1 Corresponding relation between standardized variables .and original variables 序号 规范变量 zj 自然变量 xj 1 变化间距 j 2 上星号水平 r xj r xjup 3 上水平 1 xj1 x0j ＋ j 4 零水平 0 xj0xjup＋xjlower/2xjr＋x-jr/2 5 下水平–1 x–j1 xj0 – j 6 下星号水平–r x–jr xjlower 为了配合正交试验结果的正交化操作，将以自然 变量 xj表示的方程转换为规范化变量 zj表达的方程形 式如下 ˆ y Qz 2 0 11 nn jjkjkjjjj jkjj cc zc z zc z    ， 8 式中，zz1，z2，，zn，zj为随机变量 xj的相应规 范化形式。 根据试验点组成，参照二水平正交试验表组成原 理，以基于规范化变量 zj的式（8）为拟合对象的二次 回归正交组合试验表，按照如下程序（参看实例分析 的表 3，比较直观）组构 （1）组合试验表第 1 列和第 1 栏分别为试验次数 序号和状态方程中各项名称。所以总行数、总列数相 应比试验次数和状态方程项数多 1，分别为 m1 和 n23n4/2。 （2）为计算常系数，在第 2 列安排 Z01。 （3）基本参数总个数 n，根据相应正交表 22 2 npn Lp 的前 2n-1项，可以确定二水平试验单因 素项编码，形成设计表的二水平试验编码部分（实例 分析部分表 3 中的第 1 至第 8 次试验） 。 （4） 确定星号试验次数 mr后， 先将第 1 个单因素 项依次安排两个星号试验编码“r”和“-r” ；其他依 此类推，完成星号试验编码的阶梯状分布。形成星号 试验编码部分（实例分析部分表 3 中的第 9 至第 16 次试验） 。 （5） 确定中心点试验次数 m0， 将每次试验的单因 素项编码均取为零，形成中心点试验编码部分（实例 分析部分表 3 中的第 17 次试验） 。 完成三类试验编码后，对于式（8）右边第三项交 叉项由单因素项编码相乘得到，作为影响因素排列在 表中；第四项是平方项，其编码通过对单因素项编码 取平方得到，但不满足正交性要求。为了使其满足正 交性，进行如下规范化处理。 在二次回归方程中 xj的规范化变量在式（8）中 zj的平方项为 2 j z（j1，2，，n） ，在第i次试验中 2 j z 的值为 2 ji z（j1，2，，n；i1，2，，m） 。根据 正交性定义， 2 ji z的规范化变量 ji z 的计算式为 22 1 m jijiji i zzzm     j1,2,,n 。 9 在规范化后，∑ j z 0，即具备正交性，将规范化 后的 j z 排在表中相应位置（参看后面分析实例表 3 可 以比较直观了解这一基本过程） 。 这样得到具有正交性 的二次回归组合试验设计表，然后就可以照表安排试 验。 3 状态方程确定 3.1 状态方程系数确定 第 2 期 苏永华，等. 基于二次正交试验优化的隧道围岩稳定可靠度方法 329 按照二次正交组合试验得出m次试验结果Yy1， y2，，ym后，利用试验结果和试验抽样数据，根据 回归正交试验特点，可以得出式（8）各待定系数计算 公式如下 0 1 m i i cymy    ， 10 2 11 1,2,, mm jjiiji ii cz yzjn     ， 11 2 11 mm kjkjiikji ii cz zyz z    ,1,2,,1jk kn。 12  2 11 mm jjjiiji ii czyz    。 13 3.2 各因素影响程度判别 通过正交回归得到状态方程 （8） 的具体表达式后， 则可以分析各项对状态方程的影响程度。 总平方和St为  2 2 t 111 1 mmm iii iii Syyyy m       。 14 相应地总平方和St自由度 T d1fm。 一次项偏回归平方和Sj为 22 1 m jjji i Scz    。 15 交互项偏回归平方和Skj为  2 2 1 1,2,,1 n kjkjkj i i Scz zjkkn    ；。16 平方项偏回归平方和Sjj为  2 2 1 m jjjjji i Scz    。 17 式（14）～（17）中，下标j1，2，，n。相应 各种偏回归平方和的自由度dfJ，dfKJ，dfJJ均为1。 根据上述各种偏回归平方和，得到回归平方和的 总和Sr为 rjkjjj SSSS  。 18 相应地Sr的自由度dfR为 RJKJJJ ddddffff  。 19 于是得到残差平方和Se及其自由度dfE分别为 etr SSS ， 20 ETR dddfff 。 21 对于地下岩体工程而言，显著水平可取0.1。则 回归方程中各偏回归项及Sj，Skj，Sjj及回归方程总体 和Sr的F检验统计量Fj，Fkj，Fjj，Fr分别为  eE /d jj FSSf ， 22  eE /d kjkj FSSf ， 23  eE /d jjjj FSSf ， 24  rrReE /d/dFSfSf 。 25 式（22）～（25）中，j1，2，，n；k1，2，， n–1。 对于上述各统计量，如果Fj F0.11，dfE，Fkj F0.14，12。 表 3 二次回归正交试验方案组合设计表 Table 3 Design chart of regression orthogonal combination 1 c1 d 1c1 d1c11dc1 dE12 12 c12 d2 试验 次数 Z0 Z1 Z2 Z3 Z4 Z1Z2Z1Z3Z1Z4Z2Z3Z2Z4Z3Z4 1 Z 2 Z 3 Z 4 Z 1 1 1 1 1 1 1 1 1 1 1 1 0.314 0.314 0.314 0.314 2 1 1 1 -1 -1 1 -1 -1 -1 -1 1 0.314 0.314 0.314 0.314 3 1 1 -1 1 -1 -1 1 -1 -1 1 -1 0.314 0.314 0.314 0.314 4 1 1 -1 -1 1 -1 -1 1 1 -1 -1 0.314 0.314 0.314 0.314 5 1 -1 1 1 -1 -1 -1 1 1 -1 -1 0.314 0.314 0.314 0.314 6 1 -1 1 -1 1 -1 1 -1 -1 1 -1 0.314 0.314 0.314 0.314 7 1 -1 -1 1 1 1 -1 -1 -1 -1 1 0.314 0.314 0.314 0.314 8 1 -1 -1 -1 -1 1 1 1 1 1 1 0.314 0.314 0.314 0.314 9 1 1.353 0 0 0 0 0 0 0 0 0 1.145 -0.686 -0.686 -0.686 10 1 -1.353 0 0 0 0 0 0 0 0 0 1.145 -0.686 -0.686 -0.686 11 1 0 1.353 0 0 0 0 0 0 0 0 -0.686 1.145 -0.686 -0.686 12 1 0 -1.353 0 0 0 0 0 0 0 0 -0.686 1.145 -0.686 -0.686 13 1 0 0 1.353 0 0 0 0 0 0 0 -0.686 -0.686 1.145 -0.686 14 1 0 0 -1.353 0 0 0 0 0 0 0 -0.686 -0.686 1.145 -0.686 15 1 0 0 0 1.353 0 0 0 0 0 0 -0.686 -0.686 -0.686 1.145 16 1 0 0 0 -1.353 0 0 0 0 0 0 -0.686 -0.686 -0.686 1.145 17 1 0 0 0 0 0 0 0 0 0 0 -0.686 -0.686 -0.686 -0.686 第 2 期 苏永华，等. 基于二次正交试验优化的隧道围岩稳定可靠度方法 331 表 4 相应各参数统计量 Table 4 Statistic values for various parameters 1 1 c1 d 11 统计量 Z1 Z2 Z3 Z4 Z1Z2 F 34.06 11.16 471.00 17.41 17.75 F0.11,2 8.53 8.53 8.53 8.53 8.53 1c1 1 d 1c1 1d c1d 统计量 Z1Z3 Z1Z4 Z2Z3 Z2Z4 Z3Z4 F 18.39 30.75 7.75 18.39 2.75 F0.11,2 8.53 8.53 8.53 8.53 8.53 2 1 E 2 1  2 1 c d 统计量 1 Z 2 Z 3 Z 4 Z F 23.27 2.53 16.35 5.33 F0.11,2 8.53 8.53 8.53 8.53 可知响应面函数仍然是显著的，即拟合程度高。 因此，不显著项可以剔除。在剔除不显著项后，响应 面函数为 QZ14.8870.00554Z16.005Z2133.692Z3–0.883Z4– ..1.142Z1Z2 5.821Z1Z3 – 0.883Z1Z42.778Z2Z4– 0.921Z3Z412.0285 1 Z –1.868 3 Z 。 29 将规范化变量方程还原为自然形式变量方程 gX–238.260.048E10.733 1 –32.45c11646.304d– 6.23810 -6E 1 1 0.045E1c1–0.482E1d1.513 1 d– 11.975c1d1.55710 -6 2 1 E–27.787 2 1 c 。 30 5.3 比较 按式 （26） 的求解方法， 得到极限状态方程式 （30） 的可靠度指标β为 0.824，失效概率为 20.61，最可 能失效点 MPP1 为E1， 1 ，c1，d （2830103 MPa， 26.58，0.410 -3 MPa，0.133 m） ，有限元模拟次数为 17 次。相应蒙特卡洛法（抽样在xjup和xjlower之间） 约 10000 次有限元模拟后，失效概率为 21.75，可靠 度指标为 0.78。 按照序列响应面采用不带交叉项的二次多项式， 通过 3 循环共 29 次有限元计算，可靠度指标为 0.77， 相应失效概率为 22.06， 最可能失效点 MPP2 为 （E1， 1 ，c1，d） （2760103 MPa， 24.72， 0.44810 -3 MPa， 0.128 m） ，蒙特卡洛法，有限元计算约 10000 次，失 效概率为 23.11。 采用本文方法，有限元模拟次数为蒙特卡洛模拟 法的 0.17；与蒙特卡洛法结果的绝对误差为 1.14， 相对误差为 5.24。 采用序列响应面方法，有限元次数为蒙特卡洛法 的 0.29；与蒙特卡洛法结果绝对误差为 1.05，相 对误差为 4.54。 在精度上，序列响应面法稍高；在计算效率上， 本文方法比序列响应面方法高 71。但序列响应面方 法得出的最可能失效点 MPP2 是可行域的全局解，其 中E1， 1 超出了实际的边界范围，对于该工程而言， 已经失真。因此序列响应面方法得到的可能是扩充了 范围的解答。 6 结 论 本文构建了一种响应面功能函数优化方法，解决 了序列响应面中可能导致的计算结果失真问题。在总 体上构成了具有较高效率的复杂地层隧道围岩稳定可 靠度新型分析技术。具体取得如下成果 （1） 在二次回归正交试验组合设计方法与显著性 检验结合的基础上，研究了基本参数对状态特征值影 响程度的判别准则。 （2） 将该判别准则植入隧道围岩稳定可靠度分析 状态方程的建立模式中，构建了具有自我调整功能的 响应面状态函数建立方法。 （3）研制了在变量约束条件下的基于 MATLAB 的可靠度指标搜寻技术。 （4） 通过实例分析， 展示了所研究方法的运行程 序，显示了其可行性。 参考文献 [1] 谭晓慧, 王建国, 吴礼年, 等. 边坡稳定的非线性随机有限 元加速收敛算法的研究[J]. 岩土工程学报, 2007, 297 1030–1034. TAN Xiao-hui, WANG Jian-guo, WU Li-nian, et al. Studies on accelerating convergence in nonlinear stochastic finite element analysis of slope stability[J]. Chinese Journal of Geotechnical Engineering, 2007, 297 1030–1034. in Chinese [2] 谭晓慧, 王建国. 边坡的弹塑性有限元可靠度分析[J]. 岩 土工程学报, 2007, 291 44–50. TAN Xiao-hui, WANG Jian-guo. Slope reliability analysis using elasto-plastic finite element [J]. Chinese Journal of Geotechnical Engineering, 2007, 291 44–50. in Chinese [3] 谭晓慧, 王建国, 刘新荣, 等. 边坡稳定的有限元可靠度计 算及敏感性分析[J]. 岩石力学与工程学报, 2007, 261 115–122. TAN Xiao-hui, WANG Jian-guo, LIU Xin-rong, et al. Finite element reliability computation and sensitivity analysis of slope stability[J]. Chinese Journal of Rock Mechanics and Engineering, 2007, 261 115–122. in Chinese [4] 徐 军, 刘东升, 郑颖人. 基于概率屈服准则的弹塑性随 机有限元分析[J]. 岩土工程学报, 2002, 242 225–228. XU Jun, LIU Dong-sheng, ZHENG Ying-ren. Analysis of 332 岩 土 工 程 学 报 2012 年 elastic-plastic stochastic finite element based on probabilistic yield cr