CN115115180A - 一种基于多参量分析的矿区地表塌陷风险识别与预测方法 - Google Patents

一种基于多参量分析的矿区地表塌陷风险识别与预测方法 Download PDF

Info

Publication number
CN115115180A
CN115115180A CN202210574837.1A CN202210574837A CN115115180A CN 115115180 A CN115115180 A CN 115115180A CN 202210574837 A CN202210574837 A CN 202210574837A CN 115115180 A CN115115180 A CN 115115180A
Authority
CN
China
Prior art keywords
area
risk
subsidence
mining
monitoring
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Pending
Application number
CN202210574837.1A
Other languages
English (en)
Inventor
朱权洁
谷雷
刘晓云
梁娟
于正兴
尹永明
郑贵强
滑帅
张震
隋龙琨
刘晓辉
王大仓
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
North China Institute of Science and Technology
Original Assignee
North China Institute of Science and Technology
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by North China Institute of Science and Technology filed Critical North China Institute of Science and Technology
Priority to CN202210574837.1A priority Critical patent/CN115115180A/zh
Publication of CN115115180A publication Critical patent/CN115115180A/zh
Pending legal-status Critical Current

Links

Images

Classifications

    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06Q—INFORMATION AND COMMUNICATION TECHNOLOGY [ICT] SPECIALLY ADAPTED FOR ADMINISTRATIVE, COMMERCIAL, FINANCIAL, MANAGERIAL OR SUPERVISORY PURPOSES; SYSTEMS OR METHODS SPECIALLY ADAPTED FOR ADMINISTRATIVE, COMMERCIAL, FINANCIAL, MANAGERIAL OR SUPERVISORY PURPOSES, NOT OTHERWISE PROVIDED FOR
    • G06Q10/00—Administration; Management
    • G06Q10/06—Resources, workflows, human or project management; Enterprise or organisation planning; Enterprise or organisation modelling
    • G06Q10/063—Operations research, analysis or management
    • G06Q10/0635—Risk analysis of enterprise or organisation activities
    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01V—GEOPHYSICS; GRAVITATIONAL MEASUREMENTS; DETECTING MASSES OR OBJECTS; TAGS
    • G01V1/00—Seismology; Seismic or acoustic prospecting or detecting
    • G01V1/28—Processing seismic data, e.g. analysis, for interpretation, for correction
    • G01V1/288—Event detection in seismic signals, e.g. microseismics
    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01V—GEOPHYSICS; GRAVITATIONAL MEASUREMENTS; DETECTING MASSES OR OBJECTS; TAGS
    • G01V1/00—Seismology; Seismic or acoustic prospecting or detecting
    • G01V1/28—Processing seismic data, e.g. analysis, for interpretation, for correction
    • G01V1/30—Analysis
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06Q—INFORMATION AND COMMUNICATION TECHNOLOGY [ICT] SPECIALLY ADAPTED FOR ADMINISTRATIVE, COMMERCIAL, FINANCIAL, MANAGERIAL OR SUPERVISORY PURPOSES; SYSTEMS OR METHODS SPECIALLY ADAPTED FOR ADMINISTRATIVE, COMMERCIAL, FINANCIAL, MANAGERIAL OR SUPERVISORY PURPOSES, NOT OTHERWISE PROVIDED FOR
    • G06Q50/00—Systems or methods specially adapted for specific business sectors, e.g. utilities or tourism
    • G06Q50/02—Agriculture; Fishing; Mining
    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01V—GEOPHYSICS; GRAVITATIONAL MEASUREMENTS; DETECTING MASSES OR OBJECTS; TAGS
    • G01V2210/00—Details of seismic processing or analysis
    • G01V2210/60—Analysis

Abstract

本发明涉及基于多参量分析的矿区地表塌陷风险识别与预测方法及应用，通过微震监测反演监测采空区上覆岩层移动变形位置，进而确定地表沉陷潜在沉陷区域，并进行风险等级划分，对采空区地表划分风险区域，分别对高、中、低风险区域选取合适的监测手段展开监测，可以实现节约人力、物力，可实现多参量监测手段的最大利用率，减小数据冗余，且区域小范围内的监测精度高。本发明既能实时反映现场情况，又可从机理上认识灾害发生过程；既能实现时间上的超前预测，又可达到空间上的总体把握，提高预警预测的时间、空间准确性，实现预警工作的实时性、超前性，达到预警方法的普适性，能够大大节省现场人力、物力资源的浪费。

Description

一种基于多参量分析的矿区地表塌陷风险识别与预测方法
技术领域
本发明属于矿山灾害监测预警领域，尤其涉及一种基于多参量分析的矿区地表塌陷风险识别与预测方法。
背景技术
煤矿地下开采引起岩层移动和地表沉陷的现象和过程，又称“矿山岩层及地表移动”。矿产被采出以后，采区周围岩体内部原有的力学平衡状态受到了破坏，使岩层发生了移动、变形和破坏。当开采面积达到一定范围之后，移动和破坏将波及到地表。沉陷区域对矿区及矿区周边的农业、工业、生态环境及民众生活均会造成巨大的影响与破坏。随着经济社会的发展，煤炭开采量会越开越大，沉陷区给环境、土地、水资源、建筑物等造成的灾害会越来越大。
现有的监测技术手段，在井上通常是利用测绘工程专业领域内的测量仪器手段对采空区地表沉陷区域展开监测，井下通常利用测斜仪、伸长仪等对岩层展开监测，井上井下的监测关系并不紧密，通常存在数据来源广，较为分散，缺乏统一的数据交互标准规范，不利于对多存储来源数据进行综合性、关联性分析。传统的监测手段中，通常是对地表已沉陷区域通过沉陷监测，建立模型，沉陷预测的手段，对地表沉陷展开监测预警，无法预知沉陷潜在风险区域，并展开监测，因此造成了监测资源的浪费。
发明内容
本发明要解决的技术问题是设计一种基于多参量分析的矿区地表塌陷风险识别与预测方法及应用，解决目前大面积采空区地表沉陷灾害监测手段单一、监测费用高昂、数据冗余度大等问题。
为解决上述技术问题，本发明的基于多参量分析的矿区地表塌陷风险识别与预测方法包括如下步骤：
步骤S1：通过微震监测系统获取微震事件并进行分析，对作用在采空区上覆岩层的数据进行整理，主要步骤如下：
步骤S11：通过微震处理分析软件对微震事件进行反演，获得震源点坐标；
步骤S12：通过基于LOF的K-means聚类方法对微震事件按照空间坐标x、y、z进行聚类，获得微震数据集Q0；
步骤S13：对聚类后的微震数据集按照走向、倾向绘制微震事件密度云图、微震事件能量云图，通过对比微震事件发生的空间范围、密集程度对采空区上覆岩层移动进行识别，判断事件发生的密集区域；
步骤S2：地表潜在沉陷风险区域辨识和风险等级划分，包括如下步骤：
步骤S21：地表沉陷潜在风险区域定位：将上覆岩层震源点坐标(Xn，Yn，Zn)投影至采空区地表，得到地表潜在沉陷点(Xn，Yn，Z)，n为上覆岩层震源点数；
步骤S22：地表沉陷潜在风险区域范围：对已确定的地表潜在沉陷风险区域点进行划分，通过“点连接成面”的方式辨识地表沉陷风险区域；
步骤S23：根据工作面的走向，岩层移动角确定采空区地表沉陷风险区域，同时结合地下煤层的开采计划、开采深度、开采厚度、地下采空区范围、区域地下微震事件聚集程度因素对划分好的地表潜在沉陷风险区域进行风险等级划分，将地表沉陷灾害风险等级划分为高风险区、中风险区和低风险区；采空区地表沉陷的区域半径计算公式为：
Figure BDA0003661715850000021
r为地面影响区半径，H为开采深度，β为岩层移动角；
岩层移动角的获取方式包括《导水裂缝带高度的钻孔冲洗液漏失量观测方法》中所述的井下钻探法；《建筑物、水体、铁路及主要井巷煤柱留设与压煤开采规程》中所给出的特定场景下的边界角、裂缝角；邓志毅等在《厚松散层条件下边界角和移动角求取方法》中表述的计算方法。
一般来说，在开采深度、开采厚度等因素一定的情况下，地表沉陷随着开采计划的进行，地下采空区范围不断扩大，地表沉陷区域也在不断扩大，通过辨识风险区域地下微震事件聚集程度确定地表沉陷风险等级的大小。高风险等级：地表风险区域所对应的下方基岩中微震事件密集程度高，且大能量微震事件在所有微震事件中占比高；中风险等级：地表风险区域所对应的下方基岩中微震事件密集程度较高，且中小能量微震事件在所有微震事件中占比高，大能量微震事件数量少；低风险等级：地表风险区域所对应的下方基岩中微震事件密集程度一般，且小能量微震事件在所有微震事件中占比高，中能量事件较少，大能量事件数量可以忽略不计。
地表下沉的主要原因在于地表下方基岩发生运动，且厚表土层由于强度较低，难以形成稳定的结构，原本稳定的地表会随着岩层的运动而快速下沉，造成地表沉陷灾害，所以岩层破裂是地表移动的前兆。微震监测技术可通过探测井下微小震动确定上覆岩层破坏位置、并根据能量大小预知破坏情况，因此，采用高精度检波器拾取岩层破裂震动信号，采用优化定位方法确定震源位置，能够追踪工作面开采过程中岩层破裂与运动规律，实现地表潜在沉陷风险区域辨识和风险等级划分。
步骤S3：确定矿山地表沉陷区域监测台网布设，综合监测数据并进行预测；高风险区域地表沉陷具有沉陷范围广、沉陷深度大、地下岩层破坏严重、危险性大的特点，选用D-InSAR技术、测斜仪对高风险区域进行井上-井下联合监测；中风险区域地表沉陷表现为大面积的地裂缝、塌陷沟等破坏形式，因此在监测时要特别注意中风险区向高风险区的转变，现场监测方式选用GPS监测手段与全站仪监测手段相结合，实现全天候监测，并时时关注地下开采所造成的微震事件，避免地表沉陷向高风险等级转变所造成的损失；低风险区域地表沉陷具有受灾面积小、风险等级提高可能性小的特点，在监测时选用高精度水准测量；
由于本发明通过微震监测对地表潜在沉陷区域进行提前预测，因此可以充分利用矿区已有的监测手段，对地表沉陷潜在区域提前布设监测台站，达到节约人力、物力、时间的要求。一般矿区沉陷监测手段主要有：
(1)常规大地测量方法，包括：水准测量、GPS测量、三角高程测量等。
(2)特殊测量方法，包括：倾斜仪、伸长仪、应力、应变测量、静力水准测量等。
(3)摄影及遥感测量方法，包括：航空、地面摄影测量、D-InSAR、PS-InSAR等。
本发明中，采空区地表沉陷监测网站布设不局限于某一手段，可根据矿区现有手段灵活布设，选用矿区常见测量手段进行实例讲解，包括水准测量、GPS测量、D-InSAR技术测斜仪。现对常用手段进行分析，具体内容如下：
(1)水准测量点位的布设由于采空区地下岩层破坏严重，点位布设的稳定性需要特别注意稳定性问题，存在选点困难、经济代价高、沉陷中心点误差大等困难。
(2)GPS测量具有测站间无需通视、全天候、三维信息展示、测量范围大等优点，但是其精度受环境影响较大，通常以厘米级、亚毫米级精度为目标，因此在实际应用时需要与其他手段相结合，成本较高。
(3)D-InSAR技术具有全天候、无接触、低成本、监测范围大(100km×100km)等特点，随着遥感监测的不断进步，D-InSAR技术对的精度可达到毫米级。
(4)测斜仪是通过逐段测量其测斜管与铅垂线的夹角，利用几何关系计算得到测斜管周线与铅垂线的水平偏移量。
上述监测方案中采用多种手段对采空区地表沉陷展开监测，且监测方案及手段并不局限于以上方式。
步骤S4：地表潜在沉陷风险区合理性验证并进行智能决策：当确定好采动造成的地表潜在沉陷风险区后，通过理论研究、数值模拟手段对地表沉陷区域进行验证，主要包括：
(1)通过理论计算、数值模拟手段分别获取正常开采情况过程中上覆岩层破裂高度，并与微震聚类分析结果进行对比；
(2)通过数值计算方法对未来开采计划进行模拟，获取超前数值模拟结果，并对数值模拟地表下沉值、理论分析预测值、监测手段预测值进行对比分析；
若二者满足要求，则需对地表潜在沉陷风险区做出合理规划；若二者相差很大，则需重新对地表潜在沉陷风险区重新划分。
进一步的，步骤S1中微震监测系统基于煤体破裂定位方法的分析，将2N个检波器平均分为2组，每组N个，检波器通过人工方式随工作面的推进而向前移动，每组检波器通过通讯电缆连接自采集站，采集站通过光纤连接至传输站，传输站通过网络将信息传输至中央控制记录系统。
进一步的，步骤S11中，由于P波在岩体中传播速度最快，而且初至时间易于识别，所以在一般情况下宜采用P波定位。采用此法定位时，假设岩层是均匀速度模型P波传播速度为已知，同时要在至少4个以上不同地点布设监测台站，P1，P2，P3，P4为4个监测点q点为震源的位置，S1，S2，S3，S4为震源点和监测点之间的距离，通过计算列出求解震源的方程组为：
Si＝V(ti-t0)；
Figure BDA0003661715850000041
进一步的，步骤S12中，基于LOF的K-means聚类方法的基本数学描述为：
步骤S121：计算坐标集合中的LOF值，并进行归一化处理，剔除原始震源坐标集合中的离群点，防止离群点所带来的误差，构建震源点坐标集合Q1；
LOF值计算公式为：
Figure BDA0003661715850000042
式中：lrdk(xi)是指对象xi的局部可达密度；
Figure BDA0003661715850000043
式中：
Figure BDA0003661715850000044
为对象xi的归一化值；
Figure BDA0003661715850000045
为对象xi的正则化值；S(xi)为进行正则化计算的对象，即LOF值；bases为基准值，当使用LOF算法模型时，bases＝baseLOF＝1为高斯误差函数；μs为所有对象xi正则化值的均值；σs为所有对象xi正则化值的标准差；
步骤S122：在剩余坐标集合中选取LOF值最小的点X1作为首个初始聚类中心，计算剩余各点到聚类中心的距离，将符合要求的归为一类；
步骤S123：剔除步骤S122中X1点及其聚类后构建震源点坐标集合Q2；
步骤S124：重复步骤S122选出X2作为第二个初始聚类中心点，进入步骤S123，以此往复，直到选出K个初始聚类中心；
步骤S125：利用K-means聚类方法，使用步骤S124中的K个聚类中心点，对震源坐标集合Q1进行划分，并迭代选出最优的聚类中心。
进一步的，步骤S3中，结合现场矿区提供的资料，利用开源软件GRASS GIS对所监测点位的坐标信息进行可视化展示，实现数据的融合处理，提供的资料包括矿山开采计划、CAD图纸和所在区域高精度地图。进行多源数据融合处理主要是为了避免数据在获取及使用过程中面临如下问题：①数据来源多、规范不统一、管理条例低；②数据较为分散，缺乏统一的数据交互标准规范，不利于对多存储来源数据进行综合性、关联性分析。
进一步的，步骤S3中，开采沉陷预测采用差分进化混合野狗优化算法(diffentialevolution Dingo Optimization Algorithm，DEDOA)和SVR预测模型，DEDOA-SVM算法包括如下步骤：
(1)输入训练数据和测试数据；
(2)数据归一化处理；
(3)确定差分野狗优化算法初始参数；
(4)确定野狗种群初始位置；
(5)利用下述公式分别求解种群围攻和单独捕捉到猎物时野狗的位置；
Figure BDA0003661715850000051
其中，t代表当前的迭代次数，
Figure BDA0003661715850000052
是野狗新位置；na是在[2,SizePop/2]的逆序中生成的随机整数，其中SizePop是野狗种群规模；
Figure BDA0003661715850000053
是将围攻野狗的子集，其中
Figure BDA0003661715850000054
是随机生成的野狗种群；
Figure BDA0003661715850000055
是当前野狗的位置；
Figure BDA0003661715850000056
是上一次迭代中发现的最佳野狗；βi是在[-2，2]内均匀生成的随机数，它是一个比例因子，可改变野狗轨迹的大小；
Figure BDA0003661715850000057
其中，
Figure BDA0003661715850000058
是野狗的新位置；
Figure BDA0003661715850000059
是上一次迭代中发现的最佳野狗；β2是在[-1，1]区间内均匀生成的随机数；r1是在从1到最大搜索代理(野狗)大小的区间内生成的随机数；
Figure BDA0003661715850000061
是随机选择的第r1个野狗，其中i≠r1；
(6)根据下述公式求解野狗随机搜寻腐肉时的位置；
Figure BDA0003661715850000062
其中，
Figure BDA0003661715850000063
是野狗的新位置；r1是在从1到最大搜索代理(野狗)大小的区间内生成的随机数；
Figure BDA0003661715850000064
是随机选择的第r1个野狗，其中i≠r1；σ是是算法随机生成的二进制数，σ∈{0,1}；
(7)根据下述公式(14)求解野狗种群中的生存率，目的是踢出弱小野狗，提高种群的适应能力；
Figure BDA0003661715850000065
其中，fitnessmax是当前一代中最差的适应度值；fitnessmin是当前一代中最佳的适应度值；公式中的生存向量包含[0，1]区间内的归一化适应度；
(8)更新种群；
(9)判断迭代次数，若迭代次数达到最大迭代次数择输出最优值c和g，否则返回步骤(4)继续进行；
(10)将寻优结果c和g代入到SVR模型中进行训练并预测；
(11)输出沉降预测值，用绝对误差、相对误差作为评价指标对其预测精度进行分析。
进一步的，步骤S4中，理论计算方法为：
垮落带：
Figure BDA0003661715850000066
其中，h1表示垮落带高度；m表示所采煤层的厚度；k欧式岩石碎胀系数；α表示所采煤层的倾角；
裂隙带(全部垮落)：
Figure BDA0003661715850000067
其中，h2表示裂隙带高度；n表示煤分层层数；M表示累计采厚；
Wmax＝kqm cosα，
其中，Wmax为地表最大下沉值；k为影响调整系数；α为岩层倾角；q为下沉系数。
进一步的，还包括地表沉陷灾害可视化模拟与智能决策，包括如下步骤：
(1)首先，利用开源软件GRASS GIS对数据进行可视化展示，主要操作为将获取的D-InSAR差分干涉影像与开采平面、水准测量数据以及井下开采CAD图纸数据导入GRASSGIS信息库，并对数据进行坐标转换，将所有数据处于同一坐标系下进行叠加分析；
(2)其次，通过GRASS GIS对D-InSAR差分干涉影像区域进行“重分类”，通过设置不同阈值，通过颜色变化来区分潜在沉陷区域风险等级大小；
(3)最后，结合三维重建算法，建立矿山地表三维模型；
通过上述方法获取矿区地表沉降量、预测值，实现矿区地表潜在沉陷风险区域风险等级划分。
进一步的，地表沉陷灾害可视化模拟与智能决策方法基于步骤S3中监测获取的多源数据、采空区环境背景专题库、建筑设施和道路信息构建数据库；基于步骤S4中的数值模拟模型以及预测模型构建模型库；基于专家咨询、查阅文献、行业标准和设计规程构建知识库；基于模型库识别问题、确认问题和统筹数据库提出解决办法构建方法库。
本发明还提供一种基于多参量分析的矿区地表塌陷风险识别与预测系统，包括：
一个或多个处理器；以及
一个或多个存储器，其中存储有计算机可执行程序，当由所述处理器执行所述计算机可执行程序时，执行前述基于多参量分析的矿区地表塌陷风险识别与预测方法。
本发明的基于多参量分析的矿区地表塌陷风险识别与预测方法通过微震监测反演监测采空区上覆岩层移动变形位置，进而确定地表沉陷潜在沉陷区域，并进行风险等级划分，对采空区地表划分风险区域，分别对高、中、低风险区域选取合适的监测手段展开监测，可以实现节约人力、物力，可实现多参量监测手段的最大利用率，减小数据冗余，且区域小范围内的监测精度高。传统地表沉陷风险等级评估，一般需要等煤层充分采动后，地表沉陷趋于稳定后在对风险等级进行评估，而本发明先对地下尚未开采区域或采动影响范围较小的区域展开沉陷风险评估，再针对不同的区域进行监测预测，具有重要的创新意义。本发明既能实时反映现场情况，又可从机理上认识灾害发生过程；既能实现时间上的超前预测，又可达到空间上的总体把握，提高预警预测的时间、空间准确性，实现预警工作的实时性、超前性，达到预警方法的普适性，能够大大节省现场人力、物力资源的浪费。
附图说明
下面结合附图对本发明的具体实施方式做进一步阐明。
图1为本发明的基于多参量分析的矿区地表塌陷风险识别与预测方法的总体流程示意图；
图2为实施例中微震监测系统布置示意图；
图3为微震事件追踪岩石破裂移动轨迹；
图4为震动场、位移场演化示意图；
图5为震源点定位示意图；
图6为基于LOF的K-means聚类方法流程图；
图7为微震事件走向密度云图；
图8为采空区地表沉陷监测技术路线图；
图9为本发明中DEDOA-SVR预测模型算法流程图；
图10为本发明中采空区地表潜在沉陷风险区域辨识与风险等级划分示意图；
图11(a)为实施例中微震走向微震事件云图；
图11(b)实施例中微震倾向微震事件云图；
图12为实施例中地表潜在沉陷风险等级示意图。
具体实施方式
结合图1，本发明的基于多参量分析的矿区地表塌陷风险识别与预测方法包括如下步骤：
步骤S1：通过微震监测系统获取微震事件并进行分析，对作用在采空区上覆岩层的数据进行整理，主要步骤如下：
步骤S11：通过微震处理分析软件对微震事件进行反演，获得震源点坐标；
步骤S12：通过基于LOF的K-means聚类方法对微震事件按照空间坐标x、y、z进行聚类，获得微震数据集Q0；
步骤S13：对聚类后的微震数据集按照走向、倾向绘制微震事件密度云图、微震事件能量云图，通过对比微震事件发生的空间范围、密集程度对采空区上覆岩层移动进行识别，判断事件发生的密集区域；微震事件走向密度云图如图7所示；
步骤S2：地表潜在沉陷风险区域辨识和风险等级划分，包括如下步骤：
步骤S21：地表沉陷潜在风险区域定位：将上覆岩层震源点坐标(Xn，Yn，Zn)投影至采空区地表，得到地表潜在沉陷点(Xn，Yn，Z)，n为上覆岩层震源点数；
步骤S22：地表沉陷潜在风险区域范围：对已确定的地表潜在沉陷风险区域点进行划分，通过“点连接成面”的方式辨识地表沉陷风险区域；
步骤S23：根据工作面的走向，岩层移动角确定采空区地表沉陷风险区域，同时结合地下煤层的开采计划、开采深度、开采厚度、地下采空区范围、区域地下微震事件聚集程度因素对划分好的地表潜在沉陷风险区域进行风险等级划分，将地表沉陷灾害风险等级划分为高风险区、中风险区和低风险区；采空区地表沉陷的区域半径计算公式为：
Figure BDA0003661715850000091
r为地面影响区半径，H为开采深度，β为岩层移动角；
岩层移动角的获取方式包括《导水裂缝带高度的钻孔冲洗液漏失量观测方法》中所述的井下钻探法；《建筑物、水体、铁路及主要井巷煤柱留设与压煤开采规程》中所给出的特定场景下的边界角、裂缝角；邓志毅等在《厚松散层条件下边界角和移动角求取方法》中表述的计算方法。
一般来说，在开采深度、开采厚度等因素一定的情况下，地表沉陷随着开采计划的进行，地下采空区范围不断扩大，地表沉陷区域也在不断扩大，通过辨识风险区域地下微震事件聚集程度确定地表沉陷风险等级的大小。高风险等级：地表风险区域所对应的下方基岩中微震事件密集程度高，且大能量微震事件在所有微震事件中占比高；中风险等级：地表风险区域所对应的下方基岩中微震事件密集程度较高，且中小能量微震事件在所有微震事件中占比高，大能量微震事件数量少；低风险等级：地表风险区域所对应的下方基岩中微震事件密集程度一般，且小能量微震事件在所有微震事件中占比高，中能量事件较少，大能量事件数量可以忽略不计。
地表下沉的主要原因在于地表下方基岩发生运动，且厚表土层由于强度较低，难以形成稳定的结构，原本稳定的地表会随着岩层的运动而快速下沉，造成地表沉陷灾害，所以岩层破裂是地表移动的前兆。微震监测技术可通过探测井下微小震动确定上覆岩层破坏位置、并根据能量大小预知破坏情况，因此，采用高精度检波器拾取岩层破裂震动信号，采用优化定位方法确定震源位置，能够追踪工作面开采过程中岩层破裂与运动规律，实现地表潜在沉陷风险区域辨识和风险等级划分。
步骤S3：确定矿山地表沉陷区域监测台网布设，综合监测数据并进行预测；高风险区域地表沉陷具有沉陷范围广、沉陷深度大、地下岩层破坏严重、危险性大的特点，选用D-InSAR技术、测斜仪对高风险区域进行井上-井下联合监测；中风险区域地表沉陷表现为大面积的地裂缝、塌陷沟等破坏形式，因此在监测时要特别注意中风险区向高风险区的转变，现场监测方式选用GPS监测手段与全站仪监测手段相结合，实现全天候监测，并时时关注地下开采所造成的微震事件，避免地表沉陷向高风险等级转变所造成的损失；低风险区域地表沉陷具有受灾面积小、风险等级提高可能性小的特点，在监测时选用高精度水准测量；采空区地表沉陷监测技术路线图如图8所示；
由于本发明通过微震监测对地表潜在沉陷区域进行提前预测，因此可以充分利用矿区已有的监测手段，对地表沉陷潜在区域提前布设监测台站，达到节约人力、物力、时间的要求。一般矿区沉陷监测手段主要有：
(1)常规大地测量方法，包括：水准测量、GPS测量、三角高程测量等。
(2)特殊测量方法，包括：倾斜仪、伸长仪、应力、应变测量、静力水准测量等。
(3)摄影及遥感测量方法，包括：航空、地面摄影测量、D-InSAR、PS-InSAR等。
本发明中，采空区地表沉陷监测网站布设不局限于某一手段，可根据矿区现有手段灵活布设，选用矿区常见测量手段进行实例讲解，包括水准测量、GPS测量、D-InSAR技术测斜仪。现对常用手段进行分析，具体内容如下：
(1)水准测量点位的布设由于采空区地下岩层破坏严重，点位布设的稳定性需要特别注意稳定性问题，存在选点困难、经济代价高、沉陷中心点误差大等困难。
(2)GPS测量具有测站间无需通视、全天候、三维信息展示、测量范围大等优点，但是其精度受环境影响较大，通常以厘米级、亚毫米级精度为目标，因此在实际应用时需要与其他手段相结合，成本较高。
(3)D-InSAR技术具有全天候、无接触、低成本、监测范围大(100km×100km)等特点，随着遥感监测的不断进步，D-InSAR技术对的精度可达到毫米级。
(4)测斜仪是通过逐段测量其测斜管与铅垂线的夹角，利用几何关系计算得到测斜管周线与铅垂线的水平偏移量。
上述监测方案中采用多种手段对采空区地表沉陷展开监测，且监测方案及手段并不局限于以上方式。
步骤S4：地表潜在沉陷风险区合理性验证并进行智能决策：当确定好采动造成的地表潜在沉陷风险区后，通过理论研究、数值模拟手段对地表沉陷区域进行验证，主要包括：
(1)通过理论计算、数值模拟手段分别获取正常开采情况过程中上覆岩层破裂高度，并与微震聚类分析结果进行对比；
(2)通过数值计算方法对未来开采计划进行模拟，获取超前数值模拟结果，并对数值模拟地表下沉值、理论分析预测值、监测手段预测值进行对比分析；
若二者满足要求，则需对地表潜在沉陷风险区做出合理规划；若二者相差很大，则需重新对地表潜在沉陷风险区重新划分。
本实施例优选地，步骤S1中微震监测系统基于煤体破裂定位方法的分析，将2N个检波器平均分为2组，每组N个，检波器通过人工方式随工作面的推进而向前移动，每组检波器通过通讯电缆连接自采集站，采集站通过光纤连接至传输站，传输站通过网络将信息传输至中央控制记录系统。本实施例中，微震监测系统布置如图2所示，本实施例中，取N＝10。
微震监测技术追踪岩层移动轨迹并对地表沉陷展开研究的原理为：
(1)当地下煤层开采后，采空区顶板岩层在自重应力及上覆岩层的作用下，产生向下的移动和弯曲。当其内部拉应力超过岩层的抗拉强度时，直接顶板首先断裂破碎相继冒落而老顶岩层则以梁板的形式沿层面法向方向移动、弯曲进而产生断裂离层。随着工作面不断向前推进，受到采动影响的岩层范围不断扩大。岩体内的应力状态主要与上覆岩层的状态密切相关，局部矿体的开采，引起岩体内部应力重分布，直至达到平衡，这个过程伴随着岩层移动和破坏。
(2)据经验可知，在矿山开采过程中，煤岩层在周围岩体应力作用下发生破坏时会伴随微震事件的发生，如图3所示。上覆岩层的弯曲、垮落、挤出、滑移、崩落等现象在发生到一定程度后，岩层变形会扩展到地表，表现为地裂缝、塌陷坑、沉陷沟槽、塌陷盆地等不同程度的地表沉陷现象。
(3)对微震事件进行监测，能够准确的反映地表沉陷区域，并可利用监测结果的大小、集合、位置来反映地表沉陷的风险区域。
具体可以概括为，煤矿开采损害的发展是一个缓慢的过程，随着发展高度的变化，其物理场也在不断变化。主要包括：
(1)震动场演化。随着开采范围的扩大和时间的推移，采场覆岩破裂高度不断向上发展，岩层破裂产生震动，即震源出现的高度不断变大，如图4所示。岩层移动边界与覆岩中震动场的演化特征有关。
(2)位移场演化。低位岩层破裂是高位岩层移动的前提，岩层破裂是地表移动的前提。随着开采范围的扩大和时间的推移，地表移动范围在不断发展，岩层破裂高度也在不断发展，并最终交汇，如图4所示。因此，通过追踪岩层破裂轨迹、地表下沉边界移动轨迹，能够准确得到顶板运动与地表沉陷一体化的开采损害模型。
本实施例优选地，步骤S11中，由于P波在岩体中传播速度最快，而且初至时间易于识别，所以在一般情况下宜采用P波定位。采用此法定位时，假设岩层是均匀速度模型P波传播速度为已知，同时要在至少4个以上不同地点布设监测台站，P1，P2，P3，P4为4个监测点q点为震源的位置，S1，S2，S3，S4为震源点和监测点之间的距离，震源点定位如图5所示，通过计算列出求解震源的方程组为：
Si＝V(ti-t0)；
Figure BDA0003661715850000121
本实施例优选地，如图6所示，步骤S12中，基于LOF的K-means聚类方法的基本数学描述为：
步骤S121：计算坐标集合中的LOF值，并进行归一化处理，剔除原始震源坐标集合中的离群点，防止离群点所带来的误差，构建震源点坐标集合Q1；
LOF值计算公式为：
Figure BDA0003661715850000122
式中：lrdk(xi)是指对象xi的局部可达密度；
Figure BDA0003661715850000123
式中：
Figure BDA0003661715850000124
为对象xi的归一化值；
Figure BDA0003661715850000125
为对象xi的正则化值；S(xi)为进行正则化计算的对象，即LOF值；bases为基准值，当使用LOF算法模型时，bases＝baseLOF＝1为高斯误差函数；μs为所有对象xi正则化值的均值；σs为所有对象xi正则化值的标准差；
步骤S122：在剩余坐标集合中选取LOF值最小的点X1作为首个初始聚类中心，计算剩余各点到聚类中心的距离，将符合要求的归为一类；
步骤S123：剔除步骤S122中X1点及其聚类后构建震源点坐标集合Q2；
步骤S124：重复步骤S122选出X2作为第二个初始聚类中心点，进入步骤S123，以此往复，直到选出K个初始聚类中心；
步骤S125：利用K-means聚类方法，使用步骤S124中的K个聚类中心点，对震源坐标集合Q1进行划分，并迭代选出最优的聚类中心。
本实施例优选地，步骤S3中，结合现场矿区提供的资料，利用开源软件GRASS GIS对所监测点位的坐标信息进行可视化展示，实现数据的融合处理，提供的资料包括矿山开采计划、CAD图纸和所在区域高精度地图。进行多源数据融合处理主要是为了避免数据在获取及使用过程中面临如下问题：①数据来源多、规范不统一、管理条例低；②数据较为分散，缺乏统一的数据交互标准规范，不利于对多存储来源数据进行综合性、关联性分析。
本实施例优选地，步骤S3中，开采沉陷预测采用差分进化混合野狗优化算法(diffential evolution Dingo Optimization Algorithm，DEDOA)和SVR预测模型，DEDOA-SVR预测模型算法流程图如图9所示：
1、DEDOA算法介绍
差分进化算法是根据生物进化理论提出来的一种全局搜索算法，该算法具有结构简单、参数少、鲁棒性强的优点，主要通过变异、交叉和选择3个步骤来进行。野狗优化算法是根据野狗群搜索捕猎活动所得到的一种种群智能优化算法。因为野外条件的恶劣，野狗群必须紧密合作，充分交流沟通，形成良好的合作机制才能保证狩猎成功。由于野狗群之间存在严格的等级分配关系，和良好的沟通能力，在优化过程中利用这些特点确定最优的参数值。本发明通过利用DE算法对DOA算法的搜索猎物过程中的全局搜索能力进行全局搜索，从而避免了野狗算法的局部最优问题。
(1)种群初始化
野狗种群在搜索边界内随机初始化：
Figure BDA0003661715850000131
其中，qbi表示个体上边界；mbi表示个体的下边界，randi是[0，1]之间的随机数。
(2)围攻
围攻行为顾名思义，会有数只野狗向着目标进行攻击。野狗算法中则是随机选择一定数量的野狗在全局最优位置附近进行搜索并得到一个新位置。
Figure BDA0003661715850000132
其中，t代表当前的迭代次数，
Figure BDA0003661715850000133
是野狗新位置；na是在[2,SizePop/2]的逆序中生成的随机整数，其中SizePop是野狗种群规模；
Figure BDA0003661715850000134
是将围攻野狗的子集，其中
Figure BDA0003661715850000135
是随机生成的野狗种群；
Figure BDA0003661715850000136
是当前野狗的位置；
Figure BDA0003661715850000137
是上一次迭代中发现的最佳野狗；βi是在[-2，2]内均匀生成的随机数，它是一个比例因子，可改变野狗轨迹的大小。
(3)攻击
野狗通常捕猎小猎物，直到单独捕获为止。式(11)模拟了这种行为：
Figure BDA0003661715850000141
其中，
Figure BDA0003661715850000142
是野狗的新位置；
Figure BDA0003661715850000143
是上一次迭代中发现的最佳野狗；β2是在[-1，1]区间内均匀生成的随机数；r1是在从1到最大搜索代理(野狗)大小的区间内生成的随机数；
Figure BDA0003661715850000144
是随机选择的第r1个野狗，其中i≠r1。
(4)清扫
清扫行为被定义为当野狗在它们的栖息地随意行走时找到腐肉吃的行为，用于模拟这种行为的算式为：
Figure BDA0003661715850000145
其中，
Figure BDA0003661715850000146
是野狗的新位置；r1是在从1到最大搜索代理(野狗)大小的区间内生成的随机数；
Figure BDA0003661715850000147
是随机选择的第r1个野狗，其中i≠r1；σ是是算法随机生成的二进制数，σ∈{0,1}。
(5)野狗存活率
在DOA中，野狗的存活率值计算式为：
Figure BDA0003661715850000148
其中，fitnessmax是当前一代中最差的适应度值；fitnessmin是当前一代中最佳的适应度值，该算式的生存向量包含[0，1]区间内的归一化适应度。
下式通过步骤(4)应用于低生存率，例如，生存率值等于或小于0.3。
Figure BDA0003661715850000149
其中，
Figure BDA00036617158500001410
是将被更新的生存率较低的野狗，r1和r2是在从1到搜索代理(野狗种群)的最大大小的区间内生成的随机数，r1≠r2；
Figure BDA00036617158500001411
和
Figure BDA00036617158500001412
是随机选择的第r1，r2个野狗；
Figure BDA00036617158500001413
是上一次迭代中发现的最佳野狗；σ是是算法随机生成的二进制数，σ∈{0,1}。
2、SVR模型建立如下：
假设有一组样本D＝{(x1,y1),(x2,y2),(x3,y3),…(xn,yn)},(yi∈R,i＝1,i＝2,3,…,n)，希望找到一个函数
Figure BDA00036617158500001414
式中，ω是一个特征空间向量；x→Rn；b为位移项。
引入损失函数参数ε，并使其符合式以下两个算式的条件：
Figure BDA0003661715850000151
Figure BDA0003661715850000152
如果函数无法使用损失函数参数ε估计所有数据，则误差ζ必须添加。上述两式应表示为如下形式：
Figure BDA0003661715850000153
Figure BDA0003661715850000154
式中，c是惩罚因子；ζ是误差。
通过引入拉格朗日函数和对偶变量，找到函数f(x)。但是它只能解决线性回归问题，开采沉陷预测是非线性回归问题，因此需要引入SVM核函数k(xi，xj)解决非线性回归问题。SVM核函数可以处理维数递增的复杂性问题，利用扩展定理避免高维核函数表达式，本发明是为了获得使用少量样本数据来预测开采沉陷的最佳回归模型。核函数的形式主要有3种，分别是径向基函数(radial basis function，RBF)、多项式和S形核函数，与其他核函数相比，RBF核函数具有参数少、计算速度快、可以降低数据中的噪声干扰的优势。因此，RBF核函数更适用于此回归模型：
k(xi,xj)＝e(-gamma||xi-xj||2)；
Figure BDA0003661715850000155
的惩罚因子c和核参数gamma(此后表示为g)是SVM回归模型的重要参数。提高SVM回归的准确性利用DEDOA算法，该算法用于模型参数优化。本文确定g和c的初始值，以g和c作为野狗的位置向量，使野狗在全局范围内搜索g和c的最佳值，最后建立预测方程。
DEDOA-SVM算法包括如下步骤：
(1)输入训练数据和测试数据；
(2)数据归一化处理；
(3)确定差分野狗优化算法初始参数；
(4)确定野狗种群初始位置；
(5)利用下述公式分别求解种群围攻和单独捕捉到猎物时野狗的位置；
Figure BDA0003661715850000161
其中，t代表当前的迭代次数，
Figure BDA0003661715850000162
是野狗新位置；na是在[2,SizePop/2]的逆序中生成的随机整数，其中SizePop是野狗种群规模；
Figure BDA0003661715850000163
是将围攻野狗的子集，其中
Figure BDA0003661715850000164
是随机生成的野狗种群；
Figure BDA0003661715850000165
是当前野狗的位置；
Figure BDA0003661715850000166
是上一次迭代中发现的最佳野狗；βi是在[-2，2]内均匀生成的随机数，它是一个比例因子，可改变野狗轨迹的大小；
Figure BDA0003661715850000167
其中，
Figure BDA0003661715850000168
是野狗的新位置；
Figure BDA0003661715850000169
是上一次迭代中发现的最佳野狗；β2是在[-1，1]区间内均匀生成的随机数；r1是在从1到最大搜索代理(野狗)大小的区间内生成的随机数；
Figure BDA00036617158500001610
是随机选择的第r1个野狗，其中i≠r1；
(6)根据下述公式求解野狗随机搜寻腐肉时的位置；
Figure BDA00036617158500001611
其中，
Figure BDA00036617158500001612
是野狗的新位置；r1是在从1到最大搜索代理(野狗)大小的区间内生成的随机数；
Figure BDA00036617158500001613
是随机选择的第r1个野狗，其中i≠r1；σ是是算法随机生成的二进制数，σ∈{0,1}；
(7)根据下述公式(14)求解野狗种群中的生存率，目的是踢出弱小野狗，提高种群的适应能力；
Figure BDA00036617158500001614
其中，fitnessmax是当前一代中最差的适应度值；fitnessmin是当前一代中最佳的适应度值；公式中的生存向量包含[0，1]区间内的归一化适应度；
(8)更新种群；
(9)判断迭代次数，若迭代次数达到最大迭代次数择输出最优值c和g，否则返回步骤(4)继续进行；
(10)将寻优结果c和g代入到SVR模型中进行训练并预测；
(11)输出沉降预测值，用绝对误差、相对误差作为评价指标对其预测精度进行分析。
本实施例优选地，步骤S4中，理论计算方法为：
垮落带：
Figure BDA0003661715850000171
其中，h1表示垮落带高度；m表示所采煤层的厚度；k欧式岩石碎胀系数；α表示所采煤层的倾角；
裂隙带(全部垮落)：
Figure BDA0003661715850000172
其中，h2表示裂隙带高度；n表示煤分层层数；M表示累计采厚；
Wmax＝kqm cosα，
其中，Wmax为地表最大下沉值；k为影响调整系数；α为岩层倾角；q为下沉系数。
本实施例优选地，还包括地表沉陷灾害可视化模拟与智能决策，包括如下步骤：
(1)首先，利用开源软件GRASS GIS对数据进行可视化展示，主要操作为将获取的D-InSAR差分干涉影像与开采平面、水准测量数据以及井下开采CAD图纸数据导入GRASSGIS信息库，并对数据进行坐标转换，将所有数据处于同一坐标系下进行叠加分析；
(2)其次，通过GRASS GIS对D-InSAR差分干涉影像区域进行“重分类”，通过设置不同阈值，通过颜色变化来区分潜在沉陷区域风险等级大小；
(3)最后，结合三维重建算法，建立矿山地表三维模型；
通过上述方法获取矿区地表沉降量、预测值，实现矿区地表潜在沉陷风险区域风险等级划分。
本实施例优选地，地表沉陷灾害可视化模拟与智能决策方法基于步骤S3中监测获取的多源数据、采空区环境背景专题库、建筑设施和道路信息构建数据库；基于步骤S4中的数值模拟模型以及预测模型构建模型库；基于专家咨询、查阅文献、行业标准和设计规程构建知识库；基于模型库识别问题、确认问题和统筹数据库提出解决办法构建方法库。
本实施例还提供一种基于多参量分析的矿区地表塌陷风险识别与预测系统，包括：
一个或多个处理器；以及
一个或多个存储器，其中存储有计算机可执行程序，当由所述处理器执行所述计算机可执行程序时，执行前述基于多参量分析的矿区地表塌陷风险识别与预测方法。
结合图10，本实施例的采空区地表潜在沉陷风险区域辨识与风险等级划分示意图，通过微震监测获取采空区上覆岩层微震事件，经过聚类、反演、投影，获取地表沉陷区域，通过地表下微震事件密度、能量大小对潜在沉陷区域风险等级进行划分。
本实施例中，微震走向和倾向微震事件云图如图11(a)和图11(b)所示，微震事件在x轴主要集中在3887000-3888000附近，在y轴主要集中在36491600附近，可辨识地表对应区域为潜在沉陷区域，进一步研究地表潜在沉陷风险等级如图12所示，Ⅰ标识区域为高风险区，此处潜在沉陷面积较大，可采用遥感卫星监测，Ⅱ标识区为中风险区，沉陷面积中等水平，可采用GPS测量手段，Ⅲ区域为低风险区，风险区域面积较小，由于此处有地表河流，可采用精密水准测量监测风险等级上升趋势。同时可在上述沉陷风险区对应地下开采部位布设伸长仪等监测手段，形成联合监测机制。
对上述监测所获取的数据采用DEDOA-SVR预测模型进行预测。
将联合监测结果导入GRASS GIS中进行分析，并进行智能决策。
在以上的描述中阐述了很多具体细节以便于充分理解本发明。但是以上描述仅是本发明的较佳实施例而已，本发明能够以很多不同于在此描述的其它方式来实施，因此本发明不受上面公开的具体实施的限制。同时任何熟悉本领域技术人员在不脱离本发明技术方案范围情况下，都可利用上述揭示的方法和技术内容对本发明技术方案做出许多可能的变动和修饰，或修改为等同变化的等效实施例。凡是未脱离本发明技术方案的内容，依据本发明的技术实质对以上实施例所做的任何简单修改、等同变化及修饰，均仍属于本发明技术方案保护的范围内。

Claims (10)

1.一种基于多参量分析的矿区地表塌陷风险识别与预测方法，其特征在于：包括如下步骤：
步骤S1：通过微震监测系统获取微震事件并进行分析，对作用在采空区上覆岩层的数据进行整理，主要步骤如下：
步骤S11：通过微震处理分析软件对微震事件进行反演，获得震源点坐标；
步骤S12：通过基于LOF的K-means聚类方法对微震事件按照空间坐标x、y、z进行聚类，获得微震数据集Q0；
步骤S13：对聚类后的微震数据集按照走向、倾向绘制微震事件密度云图、微震事件能量云图，通过对比微震事件发生的空间范围、密集程度对采空区上覆岩层移动进行识别，判断事件发生的密集区域；
步骤S2：地表潜在沉陷风险区域辨识和风险等级划分，包括如下步骤：
步骤S21：地表沉陷潜在风险区域定位：将上覆岩层震源点坐标(Xn，Yn，Zn)投影至采空区地表，得到地表潜在沉陷点(Xn，Yn，Z)，n为上覆岩层震源点数；
步骤S22：地表沉陷潜在风险区域范围：对已确定的地表潜在沉陷风险区域点进行划分，通过“点连接成面”的方式辨识地表沉陷风险区域；
步骤S23：根据工作面的走向，岩层移动角确定采空区地表沉陷风险区域，同时结合地下煤层的开采计划、开采深度、开采厚度、地下采空区范围、区域地下微震事件聚集程度因素对划分好的地表潜在沉陷风险区域进行风险等级划分，将地表沉陷灾害风险等级划分为高风险区、中风险区和低风险区；采空区地表沉陷的区域半径计算公式为：
Figure FDA0003661715840000011
r为地面影响区半径，H为开采深度，β为岩层移动角；
步骤S3：确定矿山地表沉陷区域监测台网布设，综合监测数据并进行预测；高风险区域地表沉陷具有沉陷范围广、沉陷深度大、地下岩层破坏严重、危险性大的特点，选用D-InSAR技术、测斜仪对高风险区域进行井上-井下联合监测；中风险区域地表沉陷表现为大面积的地裂缝、塌陷沟等破坏形式，因此在监测时要特别注意中风险区向高风险区的转变，现场监测方式选用GPS监测手段与全站仪监测手段相结合，实现全天候监测，并时时关注地下开采所造成的微震事件，避免地表沉陷向高风险等级转变所造成的损失；低风险区域地表沉陷具有受灾面积小、风险等级提高可能性小的特点，在监测时选用高精度水准测量；
步骤S4：地表潜在沉陷风险区合理性验证并进行智能决策：当确定好采动造成的地表潜在沉陷风险区后，通过理论研究、数值模拟手段对地表沉陷区域进行验证，主要包括：
(1)通过理论计算、数值模拟手段分别获取正常开采情况过程中上覆岩层破裂高度，并与微震聚类分析结果进行对比；
(2)通过数值计算方法对未来开采计划进行模拟，获取超前数值模拟结果，并对数值模拟地表下沉值、理论分析预测值、监测手段预测值进行对比分析；
若二者满足要求，则需对地表潜在沉陷风险区做出合理规划；若二者相差很大，则需重新对地表潜在沉陷风险区重新划分。
2.根据权利要求1所述的基于多参量分析的矿区地表塌陷风险识别与预测方法，其特征在于：步骤S1中微震监测系统基于煤体破裂定位方法的分析，将2N个检波器平均分为2组，每组N个，检波器通过人工方式随工作面的推进而向前移动，每组检波器通过通讯电缆连接自采集站，采集站通过光纤连接至传输站，传输站通过网络将信息传输至中央控制记录系统。
3.根据权利要求1所述的基于多参量分析的矿区地表塌陷风险识别与预测方法，其特征在于：步骤S11中采用P波定位，假设岩层是均匀速度模型P波传播速度为已知，同时要在至少4个以上不同地点布设监测台站，P1，P2，P3，P4为4个监测点q点为震源的位置，S1，S2，S3，S4为震源点和监测点之间的距离，通过计算列出求解震源的方程组为：
Si＝V(ti-t0)；
Figure FDA0003661715840000021
4.根据权利要求1所述的基于多参量分析的矿区地表塌陷风险识别与预测方法，其特征在于：步骤S12中，基于LOF的K-means聚类方法的基本数学描述为：
步骤S121：计算坐标集合中的LOF值，并进行归一化处理，剔除原始震源坐标集合中的离群点，防止离群点所带来的误差，构建震源点坐标集合Q1；
LOF值计算公式为：
Figure FDA0003661715840000022
式中：lrdk(xi)是指对象xi的局部可达密度；
Figure FDA0003661715840000023
式中：
Figure FDA0003661715840000024
为对象xi的归一化值；
Figure FDA0003661715840000025
为对象xi的正则化值；S(xi)为进行正则化计算的对象，即LOF值；bases为基准值，当使用LOF算法模型时，bases＝baseLOF＝1为高斯误差函数；μs为所有对象xi正则化值的均值；σs为所有对象xi正则化值的标准差；
步骤S122：在剩余坐标集合中选取LOF值最小的点X1作为首个初始聚类中心，计算剩余各点到聚类中心的距离，将符合要求的归为一类；
步骤S123：剔除步骤S122中X1点及其聚类后构建震源点坐标集合Q2；
步骤S124：重复步骤S122选出X2作为第二个初始聚类中心点，进入步骤S123，以此往复，直到选出K个初始聚类中心；
步骤S125：利用K-means聚类方法，使用步骤S124中的K个聚类中心点，对震源坐标集合Q1进行划分，并迭代选出最优的聚类中心。
5.根据权利要求1所述的基于多参量分析的矿区地表塌陷风险识别与预测方法，其特征在于：步骤S3中，结合现场矿区提供的资料，利用开源软件GRASS GIS对所监测点位的坐标信息进行可视化展示，实现数据的融合处理，提供的资料包括矿山开采计划、CAD图纸和所在区域高精度地图。
6.根据权利要求1所述的基于多参量分析的矿区地表塌陷风险识别与预测方法，其特征在于：步骤S3中，开采沉陷预测采用差分进化混合野狗优化算法和SVR预测模型，DEDOA-SVM算法包括如下步骤：
(1)输入训练数据和测试数据；
(2)数据归一化处理；
(3)确定差分野狗优化算法初始参数；
(4)确定野狗种群初始位置；
(5)利用下述公式分别求解种群围攻和单独捕捉到猎物时野狗的位置；
Figure FDA0003661715840000031
其中，t代表当前的迭代次数，
Figure FDA0003661715840000032
是野狗新位置；na是在[2,SizePop/2]的逆序中生成的随机整数，其中SizePop是野狗种群规模；
Figure FDA0003661715840000033
是将围攻野狗的子集，其中
Figure FDA0003661715840000034
是随机生成的野狗种群；
Figure FDA0003661715840000035
是当前野狗的位置；
Figure FDA0003661715840000036
是上一次迭代中发现的最佳野狗；βi是在[-2，2]内均匀生成的随机数，它是一个比例因子，可改变野狗轨迹的大小；
Figure FDA0003661715840000037
其中，
Figure FDA0003661715840000041
是野狗的新位置；
Figure FDA0003661715840000042
是上一次迭代中发现的最佳野狗；β2是在[-1，1]区间内均匀生成的随机数；r1是在从1到最大搜索代理大小的区间内生成的随机数；
Figure FDA0003661715840000043
是随机选择的第r1个野狗，其中i≠r1；
(6)根据下述公式求解野狗随机搜寻腐肉时的位置；
Figure FDA0003661715840000044
其中，
Figure FDA0003661715840000045
是野狗的新位置；r1是在从1到最大搜索代理大小的区间内生成的随机数；
Figure FDA0003661715840000046
是随机选择的第r1个野狗，其中i≠r1；σ是是算法随机生成的二进制数，σ∈{0,1}；
(7)根据下述公式(14)求解野狗种群中的生存率，目的是踢出弱小野狗，提高种群的适应能力；
Figure FDA0003661715840000047
其中，fitnessmax是当前一代中最差的适应度值；fitnessmin是当前一代中最佳的适应度值；公式中的生存向量包含[0，1]区间内的归一化适应度；
(8)更新种群；
(9)判断迭代次数，若迭代次数达到最大迭代次数择输出最优值c和g，否则返回步骤(4)继续进行；
(10)将寻优结果c和g代入到SVR模型中进行训练并预测；
(11)输出沉降预测值，用绝对误差、相对误差作为评价指标对其预测精度进行分析。
7.根据权利要求1所述的基于多参量分析的矿区地表塌陷风险识别与预测方法，其特征在于：步骤S4中，理论计算方法为：
垮落带：
Figure FDA0003661715840000048
其中，h1表示垮落带高度；m表示所采煤层的厚度；k欧式岩石碎胀系数；α表示所采煤层的倾角；
裂隙带：
Figure FDA0003661715840000049
其中，h2表示裂隙带高度；n表示煤分层层数；M表示累计采厚；
Wmax＝kqmcosα，
其中，Wmax为地表最大下沉值；k为影响调整系数；α为岩层倾角；q为下沉系数。
8.根据权利要求1所述的基于多参量分析的矿区地表塌陷风险识别与预测方法，其特征在于：还包括地表沉陷灾害可视化模拟与智能决策，包括如下步骤：
(1)首先，利用开源软件GRASS GIS对数据进行可视化展示，主要操作为将获取的D-InSAR差分干涉影像与开采平面、水准测量数据以及井下开采CAD图纸数据导入GRASS GIS信息库，并对数据进行坐标转换，将所有数据处于同一坐标系下进行叠加分析；
(2)其次，通过GRASS GIS对D-InSAR差分干涉影像区域进行“重分类”，通过设置不同阈值，通过颜色变化来区分潜在沉陷区域风险等级大小；
(3)最后，结合三维重建算法，建立矿山地表三维模型；
通过上述方法获取矿区地表沉降量、预测值，实现矿区地表潜在沉陷风险区域风险等级划分。
9.根据权利要求8所述的基于多参量分析的矿区地表塌陷风险识别与预测方法，其特征在于：地表沉陷灾害可视化模拟与智能决策方法基于步骤S3中监测获取的多源数据、采空区环境背景专题库、建筑设施和道路信息构建数据库；基于步骤S4中的数值模拟模型以及预测模型构建模型库；基于专家咨询、查阅文献、行业标准和设计规程构建知识库；基于模型库识别问题、确认问题和统筹数据库提出解决办法构建方法库。
10.一种基于多参量分析的矿区地表塌陷风险识别与预测系统，其特征在于：包括：
一个或多个处理器；以及
一个或多个存储器，其中存储有计算机可执行程序，当由所述处理器执行所述计算机可执行程序时，执行权利要求1-9中任一项所述的基于多参量分析的矿区地表塌陷风险识别与预测方法。
CN202210574837.1A 2022-05-25 2022-05-25 一种基于多参量分析的矿区地表塌陷风险识别与预测方法 Pending CN115115180A (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN202210574837.1A CN115115180A (zh) 2022-05-25 2022-05-25 一种基于多参量分析的矿区地表塌陷风险识别与预测方法

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN202210574837.1A CN115115180A (zh) 2022-05-25 2022-05-25 一种基于多参量分析的矿区地表塌陷风险识别与预测方法

Publications (1)

Publication Number Publication Date
CN115115180A true CN115115180A (zh) 2022-09-27

Family

ID=83326899

Family Applications (1)

Application Number Title Priority Date Filing Date
CN202210574837.1A Pending CN115115180A (zh) 2022-05-25 2022-05-25 一种基于多参量分析的矿区地表塌陷风险识别与预测方法

Country Status (1)

Country Link
CN (1) CN115115180A (zh)

Cited By (4)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN115392137A (zh) * 2022-10-27 2022-11-25 山东省地质矿产勘查开发局八〇一水文地质工程地质大队(山东省地矿工程勘察院) 一种基于岩溶塌陷水土耦合作用的三维模拟系统
CN116739310A (zh) * 2023-08-11 2023-09-12 威海市城市规划技术服务中心有限公司 一种基于地理信息的建筑场地勘测规划系统
CN116777075A (zh) * 2023-06-30 2023-09-19 江苏中煤地质工程研究院有限公司 一种采空区影响范围确定方法及系统
CN117111057A (zh) * 2023-08-23 2023-11-24 首都师范大学 一种煤矿采空区形变敏感性评价方法

Cited By (6)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN115392137A (zh) * 2022-10-27 2022-11-25 山东省地质矿产勘查开发局八〇一水文地质工程地质大队(山东省地矿工程勘察院) 一种基于岩溶塌陷水土耦合作用的三维模拟系统
CN116777075A (zh) * 2023-06-30 2023-09-19 江苏中煤地质工程研究院有限公司 一种采空区影响范围确定方法及系统
CN116777075B (zh) * 2023-06-30 2024-03-29 江苏中煤地质工程研究院有限公司 一种采空区影响范围确定方法及系统
CN116739310A (zh) * 2023-08-11 2023-09-12 威海市城市规划技术服务中心有限公司 一种基于地理信息的建筑场地勘测规划系统
CN116739310B (zh) * 2023-08-11 2023-10-24 威海市城市规划技术服务中心有限公司 一种基于地理信息的建筑场地勘测规划系统
CN117111057A (zh) * 2023-08-23 2023-11-24 首都师范大学 一种煤矿采空区形变敏感性评价方法

Similar Documents

Publication Publication Date Title
Chen et al. Machine learning-based classification of rock discontinuity trace: SMOTE oversampling integrated with GBT ensemble learning
CN115115180A (zh) 一种基于多参量分析的矿区地表塌陷风险识别与预测方法
Jiang et al. Observe the temporal evolution of deep tunnel's 3D deformation by 3D laser scanning in the Jinchuan No. 2 Mine
Pan et al. 3D scene and geological modeling using integrated multi-source spatial data: Methodology, challenges, and suggestions
Agliardi et al. Structurally-controlled instability, damage and slope failure in a porphyry rock mass
Adoko et al. Predicting tunnel convergence using multivariate adaptive regression spline and artificial neural network
Li et al. Risk assessment model for water and mud inrush in deep and long tunnels based on normal grey cloud clustering method
CN106777585B (zh) 一种区域浅表滑坡时空分异规律的esda分析法
CN108842821B (zh) 一种钻爆法修建海底隧道合理埋深的计算方法
Jiskani et al. Integrated 3D geological modeling of Sonda-Jherruck coal field, Pakistan
Wang et al. Automatic identification of rock discontinuity and stability analysis of tunnel rock blocks using terrestrial laser scanning
Xu et al. Research on monitoring and stability evaluation of ground subsidence in gypsum mine goaf
CN113376695B (zh) 一种适用于煤层底板复杂陷落柱的全波形反演方法
CN117745466A (zh) 一种尾矿库数智运维系统
Shi et al. Data-driven sequential development of geological cross-sections along tunnel trajectory
Wang et al. Multi-source data integration and multi-scale modeling framework for progressive prediction of complex geological interfaces in tunneling
CN113987095A (zh) 一种面向地下非法采矿识别的gis时空数据模型
RU2711178C1 (ru) Способ 3D прогнозирования свойств и строения геологических объектов на основе компьютерного анализа марковских свойств поверхностных геолого-геофизических полей
Tian et al. Forecasting high-probability zones of boulder occurrence for TBM tunnel engineering in tropical and subtropical granite regions
Makitaavola et al. InSAR as a practical tool to monitor and understand large-scale mining-induced ground deformations in a caving environment
Besimbayeva et al. Modern data analysis technologies used for geomechanical monitoring. Review
Hajian et al. Prior applications of neural networks in geophysics
Wu et al. Comparative study of five machine learning algorithms on prediction of the height of the water-conducting fractured zone in undersea mining
Liu Back analysis of rock mass properties in the regional fault zone under Lake Mälaren
Zhang Evaluation of mechanisms contributing to valley closure subsidence effects under irregular topographic conditions

Legal Events

Date Code Title Description
PB01 Publication
PB01 Publication
SE01 Entry into force of request for substantive examination
SE01 Entry into force of request for substantive examination