巷道围岩初始应力场和弹性模量的区间反演.doc

巷道围岩初始应力场和弹性模量的区间反演方法摘 要 从有限元数值计算的角度对岩土工程参数的可辨识性进行了分析。在可辨识性分析的基础上，考虑观测数据的不确定性，引入不确定问题求解的区间分析思想，建立了联合采用围岩位移和应力观测数据同时识别地应力和围岩弹性模量的反演模型。并采用约束变尺度方法求解该反演模型，获得地应力和弹性模量所在区间范围。关键词 区间分析，约束变尺度方法，初始应力场，弹性模量 1 引言 输入参数（岩体力学参数、初始应力场）和岩体本构关系给不准是岩土工程数值模拟的两个瓶颈问题1。采用位移等观测数据推断岩体弹性参数或初始地应力，这类弹性问题的位移反分析，由于公式推演和计算过程的实施较为简单，且所得结果通常能较好地近似反映实际工程问题的主要特征2，70年代以来受到人们的重视2-6。实际工程中，地应力和围岩弹性参数一般均是未知的，寻找能够同时确定它们的方法更具有实用性。文3曾就同时辨识围岩弹性参数和地应力可行性进行过研究，杨林德等率先将位移和扰动应力增量量测信息同时用作反演分析的依据，提出了可同时反求初始地应力参数和地层弹性参数的位移、应力反演分析法2；文4曾在引入所谓“预估约束”的条件下同时识别弹性参数和地应力。但是，由于观测误差客观存在的不确定性，据此推演出的岩土工程参数也应为不确定的，所以，反演问题本质上是一个不确定问题，本文从该角度对此类问题进行研究。从有限元数值计算的角度，本文首先就同时辨识围岩弹性参数和地应力可行性进行了探讨。然后在可辨识性分析的基础上，采用区间分析方法的思想建立了联合利用位移和应力观测数据同时反演围岩弹性模量和地应力的区间反演计算模型。在该模型中，用有限元法建立系统的控制方程，把待求变量看作区间数，以量测信息作为设计变量的边界约束，运用约束变尺度方法求出地应力和弹性模量的最大值和最小值。2 同时辨识巷道围岩弹性模量和初始地应力的可行性分析 当巷道埋深较大时，其受力情况可以简化为图1示。巷道开挖所引起的巷道围岩位移和扰动应力可以看作是，在水平压力P=0x和垂直压力Q=0y作用下与巷道形状相同的弹性孔的位移和应力减去没有开挖巷道前岩体已经完成的位移和应力。即相当于在巷道开挖边界反向作用与地应力大小相等的边界力时，所引起的围岩的位移和应力（为作图方便，以开挖边界为矩形为例如图2示）。采用有限元计算时，假设初始地应力作用下开挖边界上的作用力为f 0，开挖引起的位移和应力增量相当于在开挖边界上反向施加等效的节点力-f 0引起的围岩围岩和应力值6。 (1) 0xy0xy0x0yP=0yQ=0xyx 图1 地应力作用下的巷道 图2 地下巷道开挖释放荷载 Fig 1 A tunnel acted by in situ stress Fig 2 Released load due to excavation of tunnel 式中： 为单元形函数矩阵； 为几何矩阵； 为岩石的自重体积力矢量，当巷道埋深较大时可近似认为 =0； 为边界力， =( )； 为单元的厚度； 是具有开挖边界的所有单元。有限元平衡方程为 (2)其中 ， ， 。节点应变和应力通过式（3）求出单元应变和应力值并采用某种光滑措施（一般相当于乘一常数矩阵）后获得。 (3) 由于弹性矩阵D与弹性模量成正比， 与初始地应力成正比，应变矩阵B仅与单元几何特性有关，由式（2）知对于任意常数 ，参数组 、 、 ， ，都可以与一个位移分布 相对应，因此，仅采用位移观测数据来确定弹性模量和初始地应力场其解是不唯一的，只能确定地应力与弹性模量的比值。再由公式（3）可见，一点的应变值与其位移值二者是相关的，同样可以有无数组（E,0x,0y,0xy）与一个应变分布 相对应，即使同时采用位移和应变观测信息也不可能改变解的唯一性。但是，由于应力分布 与位移分布 通过弹性模量 联系起来，并且由应力分布 与位移分布 可望唯一确定弹性模量 ，因此，如果有扰动应力观测信息时，则可望唯一确定所有待求变量。实际上仅采用应力观测信息，或联合采用应力、位移观测信息，或联合采用应变、应力观测信息都可望达到这种辨识的目的。当观测位移信息足够多时，哪怕只增加一个应力观测信息，便可望唯一确定围岩弹性模量和地应力。考虑到反演中采用的观测量，应该是可观测、易观测的，且观测精度愈高愈好等要求，而与应力观测相比，一般而言，位移观测更易于实施、精度高等特点，本文考虑联合使用位移和应力观测信息来同时识别位移弹性模量和初始地应力。3 基于区间分析思想的反演模型的建立及其求解方法的选择 众所周知，计算误差一直是数值分析中一个比较麻烦的问题，它来源于数据误差、截断误差和舍入误差。人们努力使计算结果能在所要求的精度内，然而，在许多问题中，往往是推测计算结果的某种精度或者使用高精度的运算以保证计算结果的精度。但是，由于计算误差的积累，可能使计算结果失去意义。巷道施工期间的监控量测，由于仪器精度和其它偶然误差，所测变形是一随机值，因而反推的围岩参数和初始地应力也是随机变量。当根据带误差观测数据采用确定性方法进行反演时，反演结果必然不稳定，而且观测数据的少量误差常可以导致反演结果有大幅度的波动2。随机理论、模糊集理论和区间分析是解决不确定性问题的三种方法6-7。用随机理论或模糊集理论求解问题时需要知道不确定参数的概率密度函数或隶属函数。而实际中这些函数往往难以确定，常常通过人为方式取舍。为了真实地反映客观实际，减少人为因素的影响，提高计算结果的可靠性，运用区间分析求解不确定问题，正是这一思想的体现，它提供了一种简便的方法。考虑到各种计算误差，作为计算结果，区间分析方法得到的是包含其精确解的一个区间和误差界，这就为人们能更准确地估计和评价所得的结果提供了一定的依据。近年来，采用区间分析求解不确定问题的方法已经有所研究。区间数与普通数不同，四则运算规则也完全不同，对于两个区间数 和 ，四则运算规定为： ， ， ， ，区间数运算容易导致区间扩张。在用区间分析求解不确定问题时，有许多新问题需要解决，如如何保证解区间不被放大或缩小。文9提出采用直接优化法进行非确定结构系统区间分析，算例表明采用直接优化方法与穷举组合法所得结果相同，而直接采用区间数运算规则所得结果区间有较大扩张，这里仍然采用这种方法处理本文问题。当把待求变量（E,0x,0y,0xy）作为一个区间数时，相应地采用有限元法要区间的区间线性结构静力学方程形式为 (4)与公式(2)不同的是式(4)中 、 、 都是区间矩阵和区间向量，其元素用区间数表示为， ， ， 。在已知位移和应力观测精度的情况下，由位移和扰动应力观测信息可得 . (5)式中： 、 分别是第 个观测点位移值和第 个应力观测点扰动应力值； 、 分别是相应观测点的位移观测值和扰动应力观测值； 和 分别相应观测点的位移观测精度和扰动应力观测精度； 、 分别是位移和扰动应力观测数量。若再考虑对于实际工程问题，一般都能给出如式(6)所示的待求变量较为宽松的上下界限，以保证问题的解与实际物理意义相符。, , , (6)从而，可以建立同时确定围岩参数和原始地应力的区间直接优化计算模型： Min / Max s.t. (7) . , , , 模型(7)是一个复杂的约束非线性优化问题，本文采用约束变尺度方法求解。约束变尺度方法的的基本思想是迭代和逼近。约束变尺度方法具有收敛快、效率高、可靠性好与整体收敛性好、适应能力强等一系列优点。求解模型(7)的约束变尺度方法程序是在Powell型约束变尺度方法的基础上，参考文献11采取的监控搜索技术、尺度矩阵保持正定等改进措施编制而成。4数值算例 某地下巷道，围岩泊桑比u=0.2，假设初始地应力场均匀分布，且剪应力0xy=0（剪应力0xy不参加反演），由观测位移值和扰动应力观测值，确定围岩弹性模量E和初始地应力分量0x ，0y。 12数值模拟时首先选取一组参数值(E,0x ,0y)=(60 000, 10, 40)MPa，用有限元进行正演计算，得到图示测线、的收敛位移值，以及观测点1 和点2的水平方向应力扰动分量（见表1）。然后把这这4个收敛位移值和扰动应力值作为观测值，并考虑观测仪器的精度确定围岩弹性模量E0和原始应力分量0x ,0y所在区间范围。用有限元法计算时，由于对称，只计算半个区域，共分为265个三角形单元。各参数上下限取 =10 000Mpa， =90 000Mpa， = =2 MPa，, = =60Mpa。计算结果如表2所示，对于某一给定初值在奔腾133微机上运行计算时间为4070s。 图 3 测点布置 Fig 3 Allocation of observations 表1 位移和扰动应力观测数据Table 1 Measured data of displacements and disturbed stresses 观 测 位 置 1 2 收 敛 位 移 /10-3m 6.0092 3.06386 -2.35280 扰动应力x /MPa 8.645100 11.02070 表2 反演识别结果Table 2 Estimated results of inverse analysis /103MPa 0x /MPa 0y /MPa 精确值 60.000 10.000 40.000 识别值(1) 59.7481262.67481 9.3340110.93187 38.4360542.40228 识别值(2) 57.0552363.10662 8.9939411.05059 36.3097643.00673 _注：识别值(1)是指位移观测精度为10-4m，扰动应力观测精度为10-2MPa时的识别结果。 识别值(2)是指位移观测精度为10-4m，扰动应力观测精度为10-1MPa时的识别结果。5结语 用区间数来描述非确定参数是一种简单有效的方法。它可以减少建模过程中的人为因素的影响，提高分析结果的可靠性10。在已知观测精度