一种利用多软件联合构建复杂地质三维地应力场反演方法

1.本发明涉及地质勘查和煤矿开采技术领域，尤其涉及一种利用多软件联合构建复杂地质三维地应力场反演方法。


背景技术：

2.受资源条件的限制，我国煤炭赋存条件越发复杂。大量的矿井生产实践证明，生产过程中对地应力及围岩力学参数等基础技术数据的需要已越来越迫切但是相比较而言，煤矿地应力测量和研究工作稍落后于其他行业。目前，随采深不断增加，以静水压力理论为基础计算的地应力与实测的地应力差别较大，依此而进行的岩石力学计算已无法满足实际的需要，从而严重制约相关技术的提高。因此，近年来以地应力实测为基础的岩石力学计算在煤矿生产中大力推广应用，在此情况下研究整个矿区的地应力状况已越来越必要。
3.地应力场是地质力学与岩体力学研究的基本内容之一，不管是地面工程还是地下工程，都会直接收到地应力场的影响。在采矿领域，岩体初始地应力场是相对于采动引起的应力重新分布所形成二次应力而言的，是矿区在采掘前所具有的自然应力状态。随着采矿活动向深部转移，由于复杂地质环境和不明的应力场分布给矿山安全生产带来了隐患。通过对工程区域地应力场分布规律的分析和反演，有利于井下生产活动的开展以及对未知灾害的提前预判，为矿山生产的高效、有序、安全进行提供依据。


技术实现要素：

4.本发明主要为解决现有问题的不足之处而提供一种利用多软件联合构建复杂地质三维地应力场反演方法。
5.为解决上述技术问题，本发明采用的一个技术方案是：一种利用多软件联合构建复杂地质三维地应力场反演方法，该方法的步骤包括：根据区域内的勘探数据、地质资料等，并通过室内实验岩样物理力学参数，建立三维有限元数值模拟模型，其三维模型基于“点——线——面——体”的建模思路，构建过程中使用了surpac、3dmine、rhino、gocad、ansys、flac3d等多种软件，实现由基础钻孔信息构建复杂地层结构，生成岩层界面等高线，并以线构面、以面构体，形成更真实的复杂的地质体模型，其中岩层界面等高线是借助surpac或3dmine等数字矿山三维软件直接生成的等高线，所生成的等高线需要进行离散化处理，生成相应的等高线拟合点数据，并经过3dmine处理后导入到gocad中，提取出等高线中拟合点的数据，其中gocad导出的等高线坐标点数据有限且较为离散，借助surfer和rhino软件对离散数据进行处理，其中数据整理包括两部分内容，一是借助surfer对点阵数据进行插值优化；二是借助rhino转换ansys直接读取格式，借助ansys软件建立三维实体模型、划分网格单元，然后导入flac3d进行数值计算。
6.其为客观反演地应力的分布规律，综合考虑地形、地质构造以及地层介质等因素，结合煤矿前期的地质勘探工作，确立合理区域范围，通过该数值模拟，更精细地表征地应力场的分布规律，其中的数值模拟是根据岩层赋存特征选择相应的力学模型及各岩层的物理
力学参数。
7.在建立的地应力反演模型，在计算时主要考虑自重应力场和边界处的水平构造应力场：
8.自重应力场在临界深度以应力状态为：
[0009][0010]
式中，σ
x
、σy和σz分别为垂直自重应力分量和两个水平自重应力分量，x为西东向，y为南北向；μ为岩体的泊松比，γ、h分别上覆岩层容重和埋深，m、s分别岩体的性质常数，σc为岩体的单轴抗压强度，mpa。
[0011]
在边界处的水平构造应力场在计算模型中，取z轴为垂直向上，x、y轴所在平面为水平面，可分为x方向上的水平均匀挤压构造运动、y方向上的水平均匀挤压构造运动、xoy水平面内的均匀剪切变形构造运动、xoz平面内的竖向均匀剪切变形构造运动、yoz平面内的竖向均匀剪切变形构造运动。
[0012]
最后通过建立三维数值计算模型，利用flac3d和应力反分析技术，反演三维地应力场，并开发切片显示软件，从而实现任意点、任意面地应力特征展示，其中反演三维地应力场是在假定构造应力服从某一分布基础上给定相应参数进行计算，通过现场实测获取某一些区域的地应力值，计算指定区域的地应力场分布并带入有限元模型中进行求解。
[0013]
区别于现有技术，本发明一种利用多软件联合构建复杂地质三维地应力场反演方法，可得到直观、可视化的整体地应力分布情况，有效提取并查看模型不同水平、不同倾斜面、不同构造处的垂直应力、最大水平主应力、最小水平主应力值及方向等信息，该方法过程简单，可以弥补现有地应力反演理论中区域地应力反演精度和效率低的不足。
附图说明
[0014]
图1为本发明的地应力场反演流程示意图。
[0015]
图2为本发明的三维地质模型构建示意图。
[0016]
图3为本发明三维地应力场计算流程示意图。
具体实施方式
[0017]
为了使本发明的目的、技术方案及优点更加清楚明白，以下结合实施例，对本发明进行进一步详细说明。应当理解，此处所描述的具体实施例仅仅用以解释本发明，并不用于限定本发明。
[0018]
针对现有技术存在的问题，本发明提供了一种利用多软件联合构建复杂地质三维地应力场反演方法，下面结合附图对本发明作详细的描述。
[0019]
如图1-图3所示，本发明实施例提供的一种利用多软件联合构建复杂地质三维地应力场反演方法包括以下步骤：
[0020]
根据工程区域测试及取样信息，建立矿区的三维地质模型，构建过程中使用了surpac、3dmine、rhino、gocad、ansys、flac3d等多种软件，实现由基础钻孔信息构建复杂地层结构，并以线构面、以面构体，形成更真实的复杂的地质体模型。借助surpac或3dmine等
数字矿山三维软件生成岩层界面等高线，生成的等高线进行离散化处理，经过3dmine处理后导入到gocad中，提取出等高线中拟合点的数据，从gocad导出的等高线坐标点数据有限且较为离散，借助surfer和rhino软件对离散数据进行处理，数据处理包括两部分内容，一是借助surfer对点阵数据进行插值优化；二是借助rhino转换ansys直接读取格式，借助ansys软件建立三维实体模型、划分网格单元，然后导入flac3d进行数值计算。
[0021]
其为能客观反演地应力的分布规律，综合考虑地形、地质构造以及地层介质等因素，结合煤矿前期的地质勘探工作，确立合理区域范围，并根据岩层赋存特征选择相应的力学模型及各岩层的物理力学参数，更精细地表征地应力场的分布规律，其中各岩层的物理力学参数主要包括弹性模量/gpa、泊松比、粘聚力/mpa、抗拉强度/mpa、内摩擦角/
°
、密度/g
·
cm-3，其相应的力学模型在划分单元时：第一，单元尺寸与应力计算的误差成正比，单元划分得越小，计算结果越精确。第二，单元划分为四角形单元，其优势在于能够真实显示岩层(特别是煤层)的变形和破坏。第三，在断层带内和煤层内，单元相对小一些，同一岩层划分在同一单元组内，尽量避免同一单元处在两种性质变化较大的岩层内。
[0022]
同时，在建立的地应力反演模型计算时主要考虑自重应力场和边界处的水平构造应力场，其自重应力场在临界深度应力状态为：
[0023][0024]
式中，σ
x
、σy和σz分别为垂直自重应力分量和两个水平自重应力分量，x为西东向，y为南北向；μ为岩体的泊松比，γ、h分别上覆岩层容重和埋深，m、s分别岩体的性质常数，σc为岩体的单轴抗压强度，mpa。
[0025]
其边界处的水平构造应力场在计算模型中，取z轴为垂直向上，x、y轴所在平面为水平面，可分为自重引起的区域应力场、x轴方向的水平挤压构造应力、y轴方向的水平挤压构造应力、水平面内均匀剪切变形构造运动、xoz平面内的竖向均匀剪切变形构造运动、yoz平面内的竖向均匀剪切变形构造运动，将影响初始地应力场形成的各种因素，以边界条件和荷载的形式施加于计算模型，并在已经建好三维有限元模型上计算值导出，然后带入地应力场反演模型之中。
[0026]
利用matlab建立回归模型，在建立过程中，可以将测区内任一点处的初始地应力值σ表述为：
[0027]
σ＝f(x，y，z，e，μ，γ，δ，p2，p3，
…
)
[0028]
式中，(x，y，z)为测点的空间坐标位置；e、μ、γ分别为岩体的弹性模量、泊松比和容重；δ为自重，p2、p3代表地质构造的影响，这些因素通过给定边界条件即可求得。
[0029]
任一测点的地应力回归计算值σ
′k可表述为：
[0030][0031]
式中，表示第k个测点处的第i种应力状态，n为总的工况数；bi对应第i种工况的多元回归系数。
[0032]
任一测点处观测值与回归计算值之间的残差平方和可表示为：
[0033][0034]
式中，sjk为第k个测点处第j个应力分量真实值与回归计算值之间的残差；为第k个测点处第j个应力分量的真实测量值；表示第i工况下第k观测点j应力分量的有限差分计算值；k为测点编号(k＝1，2，3，
…
，m)；i为应力状态或工况编号(i＝1，2，3，
…
，n)；在三维问题中，j对应六个应力分量(j＝1，2，3，
…
，6)。
[0035]
回归计算方程的系数为：
[0036]
b＝[b1，b2，b3，
…
，bn]t
[0037]
式中，(b1，b2，b3，
…
)为自变量的系数，分别对应自重力、构造应力的多元回归系数。
[0038]
任一点p的多元回归初始地应力值σ
′
p
可表述为：
[0039][0040]
根据回归系数计算初始应力场回归值，可以根据回归计算结果获得相应的计算应力分量值，进而利用模型求取该位置的地应力值。
[0041]
下面结合具体实施例对本发明的技术方案作进一步的描述。
[0042]
本发明实施要点可分为如下5步：
[0043]
第一步，确立目前分析区域，根据区域内的勘探数据等资料，建立精细的大型三维有限元数值模拟模型。
[0044]
第二步，根据地应力成因理论，研究目前区域的地应力影响因素，如自重因素、水平方向的构造因素等。通过分析这些影响因素，并进行单因素有限元计算，获得各因素在单位载荷下的有限元计算值。
[0045]
第三步，根据应力线性叠加原理，建立测点的地应力回归值计算方程，其中，以测点的地应力影响因素自重、构造应力等作为自变量，地应力的回归计算值作为因变量。
[0046]
第四步，根据最小二乘法原理，建立地应力现场实测值与有限元回归计算值两者间的残差平方和方程s，并利用多组测点数据求解回归系数b。
[0047]
第五步，将多元回归系数b带入地应力回归计算方程之中，并建立区域内指定点的回归值计算模型。
[0048]
以上所述，仅为本发明的具体实施方式，但本发明的保护范围并不局限于此，任何熟悉本技术领域的技术人员在本发明揭露的技术范围内，凡在本发明的精神和原则之内所作的任何修改、等同替换和改进等，都应涵盖在本发明的保护范围之内。