基于GIS的沉陷预计分析系统及可视化研究*

2021-03-07王艳彬黄宝柱郑庆学刘义生田秀国李培现

矿山测量 2021年1期

王艳彬，黄宝柱，郭兵，郑庆学，刘义生，田秀国，李培现

(1.北京龙软科技股份有限公司，北京 100190;2.开滦(集团)有限责任公司，河北唐山 063000;3.中国矿业大学(北京) 地球科学与测绘工程学院，北京 100083)

目前，我国“三下” 压煤约为137.9亿t[1-2]，部分矿区“三下”压煤量巨大，严重制约矿区的经济社会发展，可靠的沉陷预测预报是降低矿区沉陷灾害和土地损害的关键，对指导“三下”开采实践有着重要的作用。波兰学者李特威尼申于上世纪50年代将随机介质理论引入岩层移动研究，在此基础上，我国学者刘宝琛、廖国华等[3]发展成为概率积分法。目前概率积分法已经成为我国应用最为成熟、最为广泛的预计方法。中国矿业大学、煤炭科学研究总院特采所等单位相继开发了基于概率积分法的开采沉陷预计软件。但目前大多数煤矿开采沉陷预计软件需要和CAD或其他软件平台结合，CAD不具有GIS软件的特点，对图形属性处理能力较弱，其他软件平台又和煤矿的日常生产工作联系不紧密，因此，在前人研究的基础上，依托于目前在煤矿企业应用广泛的北京龙软科技股份有限公司自主研发的LongRuan GIS平台上进行沉陷预计与分析[4-6]，旨在进行地表移动变形预计与分析的自动化处理及结果的可视化输出，降低沉陷预测工作难度，提高开采沉陷预测效率。

1 预计模型

概率积分法是一种以随机介质理论为基础的开采沉陷预计方法，概率积分法因其所用的移动和变形预计公式中含有概率积分(或真导数)而得名。其发展已较为成熟，能较精确地描述地表沉陷的分布形态[7]。目前,已成为我国应用最为广泛的倾斜煤层地表移动变形预计方法之一。

地表移动盆地内任意一点A(x,y)下沉值可以用式(1) 来计算：

(1)

式中,WA(x,y)为坐标A(x,y)的下沉值；Wcm为充分采动条件下地表最大下沉值；r为主要影响半径；D为地下开采区域。

倾斜变形i，曲率变形K，水平移动U，水平变形ε都可由下沉表达式变换求导而得。

当煤层倾角大于55°时，煤层属于急倾斜煤层。对于急倾斜煤层开采的地表下沉和变形预计，概率积分法存在明显不足，此时多采用皮尔逊Ⅲ型公式法[8]，该方法是急倾斜煤层分阶段开采地表移动变形预计的一种常用剖面函数法。

急倾斜煤层开采时走向剖面的移动和变形可按照公式(1)计算。而倾斜剖面不出现充分采动条件，在下沉盆地倾向主断面上，以底板移动边界为坐标原点、下山方向为轴的坐标系统上，倾向主断面的下沉和水平移动按式(2)和式(3)计算[3]：

W(y)=a1WmaxZa2exp (-a3Z)

(2)

u(y)=q[B-P(y)]i(y)

(3)

其中,

变形公式为(4)～(6)：

(4)

(5)

ε(y)=qBb1(1-Z)b2exp[-b3(1-Z)]

(6)

据此，可以进行地表下沉盆地倾向主断面上的全部移动变形计算。

对任意形状的工作面采用沿走向或倾向把其切割为多个近似的小矩形工作面，计算每个小矩形开采引起的同名移动变形值进行叠加,得到整个开采点引起的移动变形值[9]。

2 软件系统设计和实现

开采沉陷预计分析系统基于LongRuan GIS平台开发完成，仅使用LongRuan GIS平台即可完成煤层工作面数据预处理、沉陷预计计算和预计结果展示分析等功能。

本沉陷预计分析系统分为五大模块，即录入工作面数据模块、沉陷预计计算模块、绘制沉陷预计等值线模块、绘制沉陷变形剖面图模块、绘制损害等级云图模块。本系统中倾斜煤层开采沉陷预计采用概率积分法模型，急倾斜煤层开采倾向方向沉陷预计采用皮尔逊Ⅲ型公式法模型。

2.1 倾斜和急倾斜煤层工作面空间数据库的建立

为了实现压煤村庄搬迁与塌陷地治理补偿业务数据的标准化、规范化和存储共享化，建立了倾斜和急倾斜煤层工作面空间数据库。图1所示为倾斜煤层工作面数据管理界面，图2所示为急倾斜煤层工作面数据管理界面。以倾斜煤层工作面数据录入为例，本系统规定煤层工作面数据包括煤层参数、工作面参数和预计参数，煤层参数包括煤层倾角、煤层编号、煤层厚度，工作面参数包括指向下山方位角和任意形状工作面各个角点的(x,y)坐标和每个角点的采深，预计参数包括下沉系数、影响传播系数、水平移动系数、主要影响角正切和拐点偏移距。其中工作面的角点坐标可以从LongRuan GIS地图上拾取，也可以从Ascii文件中导入。倾斜和急倾斜煤层工作面空间数据库的建立为开采沉陷预测、煤层底板等高线修正提供有效的基础数据支持。

图1 倾斜煤层工作面数据管理界面

图2 急倾斜煤层工作面数据管理界面

本系统中录入工作面数据模块用于对空间数据库中倾斜、急倾斜煤层工作面数据统一管理，包括对工作面数据进行增、删、改、查操作，为了方便管理，在空间数据库中设置了第一等级项目名称和第二等级工作面名称两个等级，项目名称等级包含工作面名称等级，项目名称等级中包含具有同一属性的回采工作面数据，方便对同一项目下回采工作面数据预计地表的移动变形值，例如同一煤层的回采工作面数据或者同一年份的回采工作面数据。工作面数据预处理工作包括在矿井采掘工程平面图中拾取各个回采工作面边界坐标，通过图上钻孔煤层数据计算回采工作面的角点深度，获取回采工作面对应的煤层编号、煤层厚度和煤层倾角，以及煤矿测量人员提供的回采工作面对应的预计参数，以上信息组合成回采工作面数据，录入到空间数据库中。

2.2 不规则三角网模型存储

系统使用不规则三角网[10]模型存储沉陷预计结果数据，程序中自动预计了倾向和走向两个方向的地表移动变形值，可以预计指定方向的地表移动变形值。利用固定间距格网点上的移动变形值信息生成DTin数据并存储在指定的LongRuan GIS平台地图图层中。DTin数据是沉陷预计可视化图形绘制和预计分析的基础数据。

DTin数据中包括节点数据DTinNode和三角面数据DTinFace，数据结构分别如下：

struct DTinNode

{

int no;//>--节点编号

int marker;//>--节点类型，-1为默认【无类型】

double pos[3];//>--节点坐标X,Y,w, 【w为下沉值】

double value[12];//>--节点属性value[0-2]走向、倾向、指定方向倾斜值；value[3-5]走向、倾向、指定方向曲率值；value[6-8]走向、倾向、指定方向水平移动值；value[9-12]走向、倾向、指定方向的水平变形值。

}

struct DTinFace

{

int no;//>--三角面片编号

int node_no[3];//>--顶点编号

int nbour_no[3];//>--关联三角面片编号

int marker;//>--标识：外部布尔运算应用，默认为-1

}

DTinNode数据结构[11]中保存了三个方向的倾斜值、曲率值、水平移动值、水平变形值，共计有12个属性值，生成一次三角网通过属性的选取可以获取最多三个方向的倾斜值、曲率值、水平移动值、水平变形值，减少了数据存储量，提高了计算速度。

3 软件可视化交互

开采沉陷预计分析系统操作流程为：首先通过录入工作面数据模块在空间数据库中录入倾斜、急倾斜煤层工作面数据，通过沉陷预计计算模块选择工作面名称进行沉陷预计计算，并把预计结果保存到Longruan GIS指定图层下，再通过获取Longruan GIS指定图层下的预计结果绘制沉陷预计等值线、沉陷变形剖面图、损害等级云图、统计沉陷预计信息。系统通过可视化的交互界面，降低了开采地表变形分析基础数据获取的难度，提高了开采沉陷预测的工作效率。

3.1 沉陷预计计算

建立倾斜和急倾斜煤层工作面空间数据库后，利用沉陷预计计算模块，如图3所示，预计地表移动变形。沉陷预计是整个系统实现中最重要的部分，通过倾斜或急倾斜按钮读取空间数据库中的工作面数据，用户可选择工程名称以及工程名称对应的工作面名称进行组合沉陷预计，其中预计范围可通过程序自适应生成，也可手动拾取坐标的方式设置。程序中自动计算了倾向和走向两个方向的变形，还可以选择计算指定方向的变形。利用预计结果信息生成不规则三角网模型数据(以下简称DTin数据)并存储在LongRuan GIS平台指定图层下。

图3 沉陷预计计算交互界面

系统利用C++编程语言编写了基于概率积分法的预计模型算法和基于皮尔逊Ⅲ型公式法的预计模型算法，通过倾斜和急倾斜煤层工作面数据标识自动识别对应的模型算法，计算地表移动变形值。

3.2 沉陷预计等值线绘制

绘制沉陷预计等值线模块如图4所示，用于获取LongRuan GIS平台指定图层下的预计结果DTin数据，绘制下沉、倾斜、曲率、水平移动、水平变形等值线专题图及注记。通过LongRuan GIS平台中自带的“修改”菜单可以对绘制的等值线实体进行编辑、裁剪、延长、打断等操作。

图4 生成沉陷预计等值线交互界面

本系统根据间隔距离生成格网点矩阵，每个格网点(x,y)都记录了在DTinNode数据中的点号，通过选择移动变形参数、方向类型，即可获取对应的移动变形值做为格网点(x,y)的Z值，根据得到的(x,y,z)坐标值，便可通过编写的绘制等值线函数生成对应的等值线[12]专题图，如图5所示，代码为：

图5 急倾斜煤层开采下沉等值线专题图

GetPntDataByDTinNode(vSurPnt, srcDtin) ; //>-- 通过DTin数据获取格网点坐标(x,y)和对应的移动变形参数

GetIsolinesWithGridMethod_NoInter(vvIsolines, vSurPnt, iRowNum, iColNum, strIsolines, 0, 0, 0, m_dbConstPrecision, CallBack_Param)) ; //>-- 利用格网法追踪等值线-不插值提高速度

3.3 沉陷预计结果分析

LongRuan GIS平台可以查询绘制的等值线专题图中指定点的坐标、等值线的二三维长度、封闭等值线区域的面积，通过插值下沉等值线边界线面积即可以计算地表沉陷的影响范围。

系统绘制沉陷变形剖面图模块用于获取指定图层的DTin数据绘制走向和倾向的主断面图，如图6所示，以及选择剖面线绘制指定位置的断面图，方便查看建筑物、铁路、公路、河流等对应位置的下沉情况[13]。绘制损害等级云图模块可按矿区建筑物损害等级划分标准，在LongRuan GIS地图中利用倾斜、曲率、水平变形DTin数据生成临界值等值线，再通过编写的最小闭合圈算法绘制损害等级云图，由此可判断建筑物是否受到开采影响以及所受采动影响的程度[14]。

图6 走向和倾向主断面图

4 工程应用

开滦集团钱家营煤矿开采7煤，地面标高为20～23 m，采用综放开采，煤层采深为310 ～420 m，工作面位于矿井的东北部，煤层走向方位角80°～95°，煤层倾角α=7°～15°，煤层厚度4.5～5.6 m。根据邻近采区实测参数及开采工作面地质条件，选用预计参数为：下沉系数0.85，水平移动系数0.3，主要影响角正切2.5，开采影响角系数k=0.6，开采影响角90-k×α。采用本系统对开滦集团钱家营矿7煤多个倾斜煤层工作面进行沉陷预计计算，预计结果如图7所示，预计最大下沉值4.37 m，沉陷区域影响面积为6 588 870 m2，沉陷区域体积为10 590 572 m3，绘制了地表下沉等值线，确定了地面建筑物损害等级，并计算了各个等级区域的面积，与现场调查结果比对后认为，预计结果可以满足矿区生产规划要求。