多元线性回归算法在地应力场计算中的应用

王知深 李 勇,2 朱维申 董振兴

(1.山东大学岩土与结构工程中心，山东 济南 250061； 2.山东大学土建与水利学院，山东 济南 250061)

1 概述

地应力即岩体中的初始应力状态[1]，在地下工程中，它极大的影响着岩体及围岩的力学特性[2-4]，另外，也是造成地下岩体发生变形以及破坏的重要原因[5]。地应力的形成是由多种因素造成的[6,7]，随着计算机技术的发展，数值模拟在解决地下工程问题的应用中越来越广泛，但是在地下工程数值计算中，一个稳定的初始地应力场是保证计算稳定性的前提条件[8]。因此，在工程中较为准确地通过现场测试和分析计算模拟工程考察域的初始地应力场具有重要的现实意义。

目前，常用的初始地应力场反演方法主要可以分为两种：1)位移反分析法;2)应力回归法。位移反分析法是一种间接的计算方法[9,10]，该方法主要用于没有进行初始地应力场实测或者实测资料的小范围区域的初始地应力场反演。第二种方法是应力回归法，该方法能够通过少量的地应力实测数据，高效可靠地确定地下工程较大范围内的初始地应力场的情况[11,12]。因此本文拟采用第二种方法对某地下泵站洞室围岩的三维初始地应力场进行分析。

2 工程概况

该地下泵站厂区主要建筑物包括主泵房、副厂房和安装间，安装间布置在主泵房右侧，副厂房布置在左侧，出水阀室布置在地下泵站厂房下游侧，距厂房净距20 m。

该地下泵站围岩岩体结构为层状结构，地下泵站系统以Ⅲ类围岩为主，局部分布有Ⅳ类围岩，以页岩为主。主泵室部位Ⅲ类围岩约占85%；Ⅳ类围岩约占15%。主要边墙岩性为奥陶系中统(02)灰岩、白云质灰岩夹泥灰岩；下统(01)灰岩、白云质灰岩、条带状灰岩夹泥灰岩。

3 FLAC3D的基本特点

本文拟采用大型三维数值计算分析软件快速拉格朗日分析方法FLAC3D进行计算。

FLAC3D可以模拟岩土或其他材料的力学行为，尤其在大变形问题的分析方面具有独特的优势[13]。该软件在进行数值模拟时，主要通过在计算模型上建立有限差分网格，设定材料的力学参数和本构关系，并设定初始以及边界条件来对模型进行定义。然后通过显式差分法来求解微分方程。这种计算方式在求解微分方程的过程中，中间步骤高效快速，因此可以运用该方法解决一些复杂问题。

4 三维地应力场的多元线性回归原理

(1)

假定有m个观测点，则最小二乘法的残差平方和为：

(2)

由最小二乘法可得S残为最小值的方程式为：

(3)

解得L=(L1,L2,…，Ln)T并通过叠加：

(4)

其中，j为6个应力分量的方向。

最后，通过复相关系数来对结果进行筛选[14,15]。

5 区域地应力反演计算

利用国际通用的大型三维有限元计算软件FLAC3D与回归分析程序接口来计算本次厂房区地应力场的回归，最后将各个方向上的应力分量叠加得到该区域的初始地应力场。

5.1 计算模型的建立

计算模型沿三个坐标的方位分别为：X轴S向(范围为527.6 m)，Y轴E向(范围为455 m)，Z轴与大地坐标的高程方向(从海拔0 m到地标高程)，厂房轴线方向与X轴平行。如图1，图2所示，模型不同颜色代表不同的地质情况。对计算区域共划分了518 965个单元和78 599个节点。岩体力学参数建议值见表1。

表1 岩体力学参数

类别密度/kN·m-3弹性模量/GPa泊松比抗剪强度凝聚力内摩擦角围岩Ⅲ27.80270.265.0345破碎带等Ⅳ21.250.80.230.642

5.2 实测地应力分析

在本次地应力回归计算域内，共有2个钻孔，为ZK11x-5和ZK12X-7。钻孔位置如表2所示。

在本次地应力分析中采用水压致裂法进行地应力测试。在本次计算中，两个钻孔共有13个测点，选取其中规律性较好的4个测点进行优化计算，其他测点数值因离散性较大来考虑在拟合计算中，其中4个测点的地应力位置和数值如表2所示。

表2 多元回归分析中采用的实测地应力位置及数值

5.3 多元回归分析

5.3.1计算工况

首先，将在工程现场实际监测得到的主应力，按照计算设定坐标分解为6个应力分量。然后根据自重，X和Y方向上的挤压力以及XY和YX方向上的剪切力确定工程现场的初始地应力的拟合载荷。各个计算模型工况与边界条件设置见表3。

表3 计算模型工况与边界条件设置

5.3.2系数回归

通过计算各种工况可以得知，当工况为自重+X向梯形构造应力+Y向梯形构造应力+XY剪切构造应力+YX剪切构造应力时，相关系数的取值最大。在这种工况条件下，复相关系数r=0.992 545 744，说明在这种工况下回归较好。回归系数l1=686 503.352 9，l2=0.959 324，l3=0.222 619，l4=-0.434 105,l5=0.196 296，l6=0.615 648 6。

在上述条件下，将计算结果与各个测点测得的实际应力分量进行对比，对比结果如图3，图4所示。由图3和图4可知，回归值和实测值变化趋势以及数值大小是基本吻合的。

经过计算，在包括厂房及出水阀室的中心区域，在x方向的侧压系数Kx约1.003，而在y方向Ky则为1.4左右。

6 结语

本次厂房区地应力场回归采用国际通用的大型三维有限元计算软件FLAC3D与回归分析相结合，通过叠加计算得到厂房区域内的初始地应力长分布情况,获得如下结论：

1)通过对地应力的计算值和实测值对比以及对侧压比规律的综合分析，表明得到的初始应力场能够较好地反映地形、地貌和地质构造的影响，其结果具有较好的可靠性，可作为厂房设计、施工和稳定性评价的参考依据。2)从较直观的角度来说，在包括厂房及出水阀室的中心区域，在x方向的侧压系数Kx约1.003，而在y方向Ky则为1.4左右。