一种兼顾风积沙流动量的煤矿开采地表沉陷预测方法

1.本发明涉及一种煤矿开采地表沉陷预测方法，尤其适用于一种兼顾风积沙流动量的煤矿开采地表沉陷预测方法，属于开采沉陷盆地的地表风积沙流动量测量领域。
技术背景
2.风积沙是地质作用形成的分布于地表及地表以下一定深度范围的并受风力作用改造的松散的沙质沉积物，在我国鄂尔多斯矿区广泛分布。由于风积沙压缩系数低、含水量低、流动性强，在采动影响下受到风的搬运作用较为明显，因此风积沙区开采的地表沉陷特征与一般地质采矿条件下开采时存在显著区别，具体表现为下沉值偏小，盆地边缘区抬升，盆地特征与流动量相关联等，这导致概率积分预测方法在风积沙矿区的适用性变差。并且由于煤矿开采和风积沙流动共同作用于地表土层的位移变化，因此直接对地表进行位移和高程的测量无法确定风积沙的流动方向和流动量，难以界定矿区风积沙运移对地表沉陷的影响程度。同时鄂尔多斯盆地煤炭和天然气赋存地点具有趋同性，上述问题的存在严重制约了气-煤交叉区天然气工程设施的建设，对采气井和埋地管线造成了安全威胁。因此，如何实现沉陷区风积沙流动量的精准测量、量化风积沙与沉陷特征之间的关联以及计算全盆地风积沙的流动量、实现矿区地表沉陷的精准预测是亟需解决的难题。
3.公开号为cn106593524a的现有技术：公开一种固体充填开采地表沉陷动态预计方法，将固体充填开采顶板随时间变化的下沉模型与knothe时间函数相结合，采用叠加原理将不同时刻的充填体压缩变形引起的地表下沉进行累加，建立了一套固体充填开采地表动态预计模型。采用本发明的方法，根据岩层时间影响系数、充填体压缩变形时间影响系数、采深、采厚、煤层倾角、预计工作面尺寸以及概率积分法预计参数等参数，可以预计出固体充填开采地表沉陷随时间变化的动态变形值。但是该方法计算复杂，同时并未考虑沙漠化环境中风沙对地表的影响，在我国西北等地的环境中使用预测精度低。


技术实现要素：

4.针对上述技术中的问题，提供一种步骤简单、精准度高且兼顾采动影响和风积沙流动特点的地表风积沙厚度变化测量及计算方法。
5.为实现上述技术目的，本发明的一种兼顾风积沙流动量的煤矿开采地表沉陷预测方法，其步骤如下：
6.首先在采动工作面地表不同位置均匀间隔布设附有测钉的风积沙厚度观测站，即沿走向和倾向固定距离铺设透水隔沙聚酯密目网，并且网格四周用长测钉固定；
7.之后采用实时动态rtk测量技术和钢尺量距方法进行长周期观测，获取不同采动程度下地表沉陷盆地不同位置地表移动变形值和地表风积沙厚度变化值δh
沙
，地表移动变形值包括下沉w
测
、倾斜i
测
、曲率k
测
、水平移动u
测
、水平变形ε
测
；
8.然后将地表移动变形值和地表风积沙厚度变化值δh
沙
进行相关分析和多元线性回归，得到采动影响下风积沙厚度变化的预测公式；
9.最后依据移动观测分析报告中记载的本矿区的概率积分预计参数，利用概率积分法预计地表下沉值w
预
，根据已建立的风积沙厚度变化预测公式计算出风积沙流动造成的地表沉陷值w
沙
，并且能通过二者相减得出单一采动因素条件下矿区地表实际沉陷值w
实
。
10.进一步，沿规划开采工作面走向和倾向，间隔固定距离布设附有测钉的风积沙厚度观测站，风积沙厚度观测站起始点应处于不受采动干扰的区域；所述观测站采取地表铺设2m
×
2m透水隔沙聚酯密目网均匀覆盖0.2m松软表土，透水隔沙聚酯密目网四周用1m长钢质测钉固定，钢质测钉埋深》0.5m。
11.进一步，随着工作面开采的推进，定期观测风积沙厚度观测站的测钉三维坐标cx
,y,z
及风积沙厚度变化值δh
沙
；
12.对每个观测站的四个测钉均进行rtk坐标点位测量，获得所有观测站测钉的三维坐标c
x,y,z
的数据集，将观测站的四个测钉的x、y、z轴坐标平均值作为该观测站的点位坐标值以代表该网格区域的移动变形特征值；
13.通过第i个观测站不同周期三维坐标观测结果的分析和计算，得出每个观测站在不同采动程度下的地表移动变形值：下沉w
测
、倾斜i
测
、曲率k
测
、水平移动u
测
、水平变形ε
测
，计算公式如下：
14.3.1.m次观测时第i个观测站的下沉：
[0015][0016]
式中，wn是第i个观测站的下沉值；h
im
分别为首次和m次观测时第i个观测站的高程；
[0017]
3.2.相邻两个观测站的倾斜：
[0018]
x轴(倾向)方向倾斜分量：
[0019]
y轴(走向)方向倾斜分量：
[0020]
总倾斜量：
[0021]
式中，δx
i-i+1
为第i个观测站和第i+1个观测站的x坐标差；δy
i-i+1
为第i个观测站和第i+1个观测站的y坐标差；w
i+1
,wi分别表示第i+1个观测站和第i个观测站的下沉量；
[0022]
3.3第i个观测站附近的曲率，即第i-1个观测站与第i+1个观测站之间的曲率：
[0023][0024]
式中，i
i-1-i
,i
i-i-1
分别表示第i+1个观测站至第i个观测站和第i个观测站至第i-1个观测站的倾斜；δx
i+1-i
,δy
i+1-i
,δy
i-i-1
分别表示第i+1个观测站至第i个观测站和第i个观测站至第i-1个观测站的x坐标差和y坐标差；
[0025]
3.4.第i个观测站的水平移动：
[0026]
x轴方向水平移动量：u
ix
＝x
im-x
i0
[0027]
y轴方向水平移动量：u
iy
＝y
im-y
i0
[0028]
总水平移动量：
[0029]
式中，x
im
,x
i0
,y
im
,y
i0
分别表示m次观测和首次观测第i个观测站的横、纵坐标；
[0030]
3.5.第i个观测站至i+1个观测站间的水平变形：
[0031]
x轴方向水平变形量：
[0032]
y轴方向水平变形量：
[0033]
总水平变形量：
[0034]
式中，(x
i+1～i
)m,(x
i+1～i
)0,(y
i+1～i
)m,(y
i+1～i
)0分别表示第i+1个观测站和第i个观测站在m次观测时和首次观测时的横纵坐标。
[0035]
通过毫米钢尺量取地表至聚酯密目网之间的距离减去初始覆土厚度值获得风积沙厚度变化值δh
沙
，在每个透水隔沙聚酯密目网内选取4～8处起伏变化处测量高度，结果取平均值作为第i个观测站的风积沙厚度变化值考虑季节性风速和降水量的影响，应持续观测1年/12期及以上；
[0036]
进一步，通过在不同采动程度条件下不同位置的地表移动变形值w
测
、i
测
、k
测
、u
测
、ε
测
和风积沙厚度变化值δh
沙
的获取，对二者进行相关分析和线性回归，得出沉陷全盆地任意位置关于地表移动变形值的风积沙流动量预测方法；
[0037]
其中采动影响下风积沙流动量的回归公式具体为：
[0038]
δh
沙
＝v1w
测
+v2i
测
+v3k
测
+v4u
测
+v5ε
测
+b
[0039]
式中，w为地表下沉值；i为地表倾斜值；k为地表曲率变形值；u为地表水平移动值；ε为地表水平变形值，v1、v2、v3、v4、v5为相关系数，由已有测量结果反演得出；b为常数项。
[0040]
进一步，结合矿区采掘规划、采掘工程平面图以及钻孔柱状图资料，收集整理规划开采工作面的相关信息，包括工作面面宽a、面长b、采深h、采厚m、倾角α，以及预测参数的信息，包括下沉系数q、主要影响角正切tanβ，拐点偏移距s、开采影响传播角θ等信息，采用概率积分法进行工作面开采地表移动变形预测，得到不同位置的地表移动及变形预测值w
测
、i
测
、k
测
、u
测
、ε
测
；
[0041]
概率积分法预测公式为：
[0042]
5.1.任意点a(x，y)的下沉值w(x,y)
[0043]
w(x，y)＝w
cmc′
xc′y[0044]
其中：
[0045][0046][0047]
式中：w
cm
表示充分采动条件下地表最大下沉值，w
cm
＝mq cosα；m表示采出煤层厚度；q表示地表下沉系数；α表示煤层倾角；c'
x
，c'y表示为待求点在走向和倾向主断面上投影
点处的下沉分布系数；l，l表示为采区拐点平移后走向长度及倾斜方向在地表的计算开采宽度；r，r1，r2表示分别为走向、下山、上山方向的主要影响半经；x、y表示待求点坐标；
[0048]
5.2.地表任意点a(x,y)沿方向倾斜变形值为：
[0049][0050][0051][0052][0053]
式中：t(x,y)m表示待求点的最大倾斜值表示最大倾斜值方向与ox轴的夹角，即沿逆时针方向旋转的角度，t
x
，ty表示分别为待求点沿走向和倾向主断面上投影点处迭加后的倾斜变形值，
[0054]
5.3.地表任意点a(x，y)沿φ方向的曲率变形为：
[0055][0056][0057][0058][0059][0060]
式中：k(x,y)
max
，k(x,y)
min
表示分别为待求点最大、最小曲率变形值；k
x
，ky表示分别为待求点沿走向及倾向在主断面投影处迭加后的曲率值；
[0061]
5.4.地表任意点a(x，y)沿φ方向的水平移动值为：
[0062][0063][0064]
u(x,y)
cm
＝u
xc′ycosφu+u
yc′
x
sinφu[0065]
φu＝arctg(u
yc′
x
/u
xc′y)u
x
，
[0066]
式中：φu表示为最大水平移动方向与ox轴的夹角，uy表示分别为待求点沿走向和倾向在主断面投影点处的水平移动值，对于倾斜方向需加cy·wcm
·
cotθ。
[0067]
5.5.地表任意点a(x，y)沿方向的水平变形值为：
[0068][0069][0070]
[0071][0072][0073]
式中：ε(x,y)
max
，ε(x,y)
min
表示为待求点最大、最小水平变形值；ε
x
，εy表示为待求点沿走向及倾向在主断面投影处迭加的水平变形值。
[0074]
进一步，根据地表移动变形的概率积分预测值，结合已建立的风积沙厚度变化公式，可以推导出风积沙流动量对地表下沉的影响大小：
[0075]w沙
＝v1w
预
+v2i
预
+v3k
预
+v4u
预
+v5ε
预
+b
[0076]
并且可以计算出剔除风积沙流动量影响的煤矿开采地表实际下沉值w
实
：
[0077]w实
＝w
预-w
沙
。
[0078]
有益效果：本方法在考虑自然因素和采煤因素的基础上实现了采动影响的风积沙流动量精准测量，构建了地表移动变形指标和风积沙流动量之间的模型关系，在矿区常用概率积分法预测的基础上能实现风积沙流动量的精准计算，可为界定风积沙对采煤沉陷的影响程度、实施气-煤交叉区建(构)筑物精准保护提供技术支撑。
附图说明
[0079]
图1为本发明实施例中风积沙流动对沉陷范围及程度的影响示意图；
[0080]
图2为本发明实施例中煤矿开采地表沉陷盆地不同位置风积沙流动量测量装置示意图；具体为透水隔沙聚酯密目网的铺设俯视图和斜视图；
[0081]
图3为本发明的采煤沉陷盆地风积沙流动量测量及计算方法示意图。
具体实施方式
[0082]
下面将结合图和具体实施过程对本发明做进一步详细说明：
[0083]
如图1和图2所示，本发明的一种兼顾风积沙流动量的煤矿开采地表沉陷预测方法，首先在采动工作面地表不同位置均匀间隔布设附有测钉的风积沙厚度观测站，即沿走向和倾向固定距离铺设透水隔沙聚酯密目网，并且网格四周用长测钉固定；之后采用rtk测量技术和钢尺量距方法进行长周期观测，获取忽略风积沙流动影响不同位置的地表移动变形值(下沉w
测
、倾斜i
测
、曲率k
测
、水平移动u
测
、水平变形ε
测
)和地表风积沙厚度变化值δh
沙
；然后将地表移动变形值和地表风积沙厚度变化值δh
沙
进行相关分析和多元线性回归，得到采动影响下风积沙流动量的计算公式；最后根据建立的公式，利用概率积分法计算相应的地表沉陷预测值w
预
，可以根据已建立的风积沙厚度变化预测公式计算出风积沙流动造成的地表沉陷值w
沙
，并能计算出剔除风积沙流动量δh
沙
影响的煤矿开采地表实际沉陷量w
实
。
[0084]
如图3所示具体步骤为：
[0085]
步骤1：首先沿规划开采工作面走向和倾向，间隔固定距离(50m)布设附有测钉的风积沙厚度观测站，观测站起始点应处于不受采动干扰的区域。观测站采取地表铺设2m
×
2m透水隔沙聚酯密目网并均匀覆盖0.2m松软表土，网格四周用1m长钢质测钉固定(埋深》0.5m)；
[0086]
步骤2：伴随工作面开采的推进，周期(30d)观测测钉点位坐标、高程及风积沙厚度变化值。测钉的三维坐标c
x,y,z
采取rtk测量方法，取四个测钉三维坐标测量值的平均值作为
第i个观测点坐标值风积沙厚度变化值δh
沙
采取毫米钢尺量取地表至聚酯密目网之间的距离，在格网内选取4～8处起伏变化处测量高度，结果取平均值作为第i个观测站的风积沙厚度变化值δhi。考虑季节性风速和降水量的影响，应持续观测1年(12期)及以上；
[0087]
步骤3：收集地表风积沙厚度变化δh
沙
观测的成果，通过与地表变形实测值(下沉w
测
、倾斜i
测
、曲率k
测
、水平移动u
测
、水平变形ε
测
)进行相关分析，根据相关分析进行多元线性回归得到采动影响下的地表风积沙流动量δh
沙
的计算公式：
[0088]
δh
沙
＝v1w
测
+v2i
测
+v3k
测
+v4u
测
+v5ε
测
+b
[0089]
多元线性回归公式形式如文中所示，计算过程采用现有软件如matlab进行计算，过程不予给出。
[0090]
步骤4：结合矿区采掘规划、采掘工程平面图以及钻孔柱状图资料，收集整理规划开采工作面的相关信息，包括工作面面宽a、面长b、采深h、采厚m、倾角α，以及预测参数的信息，包括下沉系数q、主要影响角正切tanβ，拐点偏移距s、开采影响传播角θ等信息。通过对矿区类似地质条件下规划开采工作面进行地表移动变形预测，可以得到规划工作面的不同位置的地表移动及变形预测值(下沉w
预
、倾斜i
预
、曲率k
预
、水平移动u
预
、水平变形ε
预
)，然后将规划工作面的预测值带入步骤3得到的风积沙流动量计算公式，就可以计算出规划开采工作面采后地表任意位置由风积沙流动引起的下沉值w
沙
，并可以计算出剔除风积沙流动量δh
沙
影响的煤矿开采地表实际沉陷量w
实
。其中：
[0091]w沙
＝v1w
预
+v2i
预
+v3k
预
+v4u
预
+v5ε
预
+b
[0092]w实
＝w
预-w
沙
式中：v1,v2,v3,v4,v5,b为已知量，由步骤3反演得到。

技术特征：
1.一种兼顾风积沙流动量的煤矿开采地表沉陷预测方法，其特征在于步骤如下：首先在采动工作面地表不同位置均匀间隔布设附有测钉的风积沙厚度观测站，即沿走向和倾向固定距离铺设透水隔沙聚酯密目网，并且网格四周用长测钉固定；之后采用实时动态rtk测量技术和钢尺量距方法进行长周期观测，获取不同采动程度下地表沉陷盆地不同位置地表移动变形值和地表风积沙厚度变化值δh
沙
，地表移动变形值包括下沉w
测
、倾斜i
测
、曲率k
测
、水平移动u
测
、水平变形ε
测
；然后将地表移动变形值和地表风积沙厚度变化值δh
沙
进行相关分析和多元线性回归，得到采动影响下风积沙厚度变化的预测公式；最后依据移动观测分析报告中记载的本矿区的概率积分预计参数，利用概率积分法预计地表下沉值w
预
，根据已建立的风积沙厚度变化预测公式计算出风积沙流动造成的地表沉陷值w
沙
，并且能通过二者相减得出单一采动因素条件下矿区地表实际沉陷值w
实
。2.根据权利要求1所述的一种兼顾风积沙流动量的煤矿开采地表沉陷预测方法，其特征在于：沿规划开采工作面走向和倾向，间隔固定距离布设附有测钉的风积沙厚度观测站，风积沙厚度观测站起始点应处于不受采动干扰的区域；所述观测站采取地表铺设2m
×
2m透水隔沙聚酯密目网均匀覆盖0.2m松软表土，透水隔沙聚酯密目网四周用1m长钢质测钉固定，钢质测钉埋深>0.5m。3.根据权利要求1所述的一种兼顾风积沙流动量的煤矿开采地表沉陷预测方法，其特征在于：随着工作面开采的推进，定期观测风积沙厚度观测站的测钉三维坐标c
x,y,z
及风积沙厚度变化值δh
沙
；对每个观测站的四个测钉均进行rtk坐标点位测量，获得所有观测站测钉的三维坐标c
x,y,z
的数据集，将观测站的四个测钉的x、y、z轴坐标平均值作为该观测站的点位坐标值以代表该网格区域的移动变形特征值；通过第i个观测站不同周期三维坐标观测结果的分析和计算，得出每个观测站在不同采动程度下的地表移动变形值：下沉w
测
、倾斜i
测
、曲率k
测
、水平移动u
测
、水平变形ε
测
，计算公式如下：3.1.m次观测时第i个观测站的下沉：式中，w
n
是第i个观测站的下沉值；h
im
分别为首次和m次观测时第i个观测站的高程；3.2.相邻两个观测站的倾斜：x轴(倾向)方向倾斜分量：y轴(走向)方向倾斜分量：总倾斜量：式中，δx
i-i+1
为第i个观测站和第i+1个观测站的x坐标差；δy
i-i+1
为第i个观测站和第i+1个观测站的y坐标差；w
i+1
,w
i
分别表示第i+1个观测站和第i个观测站的下沉量；3.3第i个观测站附近的曲率，即第i-1个观测站与第i+1个观测站之间的曲率：
式中，i
i-1-i
,i
i-i-1
分别表示第i+1个观测站至第i个观测站和第i个观测站至第i-1个观测站的倾斜；δx
i+1-i
,δy
i+1-i
,δy
i-i-1
分别表示第i+1个观测站至第i个观测站和第i个观测站至第i-1个观测站的x坐标差和y坐标差；3.4.第i个观测站的水平移动：x轴方向水平移动量：u
ix
＝x
im-x
i0
y轴方向水平移动量：u
iy
＝y
im-y
i0
总水平移动量：式中，x
im
,x
i0
,y
im
,y
i0
分别表示m次观测和首次观测第i个观测站的横、纵坐标；3.5.第i个观测站至i+1个观测站间的水平变形：x轴方向水平变形量：y轴方向水平变形量：总水平变形量：式中，(x
i+1～i
)
m
,(x
i+1～i
)0,(y
i+1～i
)
m
,(y
i+1～i
)0分别表示第i+1个观测站和第i个观测站在m次观测时和首次观测时的横纵坐标。通过毫米钢尺量取地表至聚酯密目网之间的距离减去初始覆土厚度值获得风积沙厚度变化值δh
沙
，在每个透水隔沙聚酯密目网内选取4～8处起伏变化处测量高度，结果取平均值作为第i个观测站的风积沙厚度变化值考虑季节性风速和降水量的影响，应持续观测1年/12期及以上。4.根据权利要求1所述的一种兼顾风积沙流动量的煤矿开采地表沉陷预测方法，其特征在于：通过在不同采动程度条件下不同位置的地表移动变形值w
测
、i
测
、k
测
、u
测
、ε
测
和风积沙厚度变化值δh
沙
的获取，对二者进行相关分析和线性回归，得出沉陷全盆地任意位置关于地表移动变形值的风积沙流动量预测方法；其中采动影响下风积沙流动量的回归公式具体为：δh
沙
＝v1w
测
+v2i
测
+v3k
测
+v4u
测
+v5ε
测
+b式中，w为地表下沉值；i为地表倾斜值；k为地表曲率变形值；u为地表水平移动值；ε为地表水平变形值，v1、v2、v3、v4、v5为相关系数，由已有测量结果反演得出；b为常数项。5.根据权利要求1所述的一种兼顾风积沙流动量的煤矿开采地表沉陷预测方法，其特征在于：结合矿区采掘规划、采掘工程平面图以及钻孔柱状图资料，收集整理规划开采工作面的相关信息，包括工作面面宽a、面长b、采深h、采厚m、倾角α，以及预测参数的信息，包括下沉系数q、主要影响角正切tanβ，拐点偏移距s、开采影响传播角θ等信息，采用概率积分法进行工作面开采地表移动变形预测，得到不同位置的地表移动及变形预测值w
测
、i
测
、k
测
、u
测
、ε
测
；概率积分法预测公式为：
5.1.任意点a(x，y)的下沉值w(x,y)w(x,y)＝w
cm
c
′
x
c
′
y
其中：其中：式中：w
cm
表示充分采动条件下地表最大下沉值，w
cm
＝mqcosα；m表示采出煤层厚度；q表示地表下沉系数；α表示煤层倾角；c'
x
，c'
y
表示为待求点在走向和倾向主断面上投影点处的下沉分布系数；l，l表示为采区拐点平移后走向长度及倾斜方向在地表的计算开采宽度；r，r1，r2表示分别为走向、下山、上山方向的主要影响半经；x、y表示待求点坐标；5.2.地表任意点a(x,y)沿方向倾斜变形值为：为：为：为：式中：t(x,y)
m
表示待求点的最大倾斜值表示最大倾斜值方向与ox轴的夹角，即沿逆时针方向旋转的角度，t
x
，t
y
表示分别为待求点沿走向和倾向主断面上投影点处迭加后的倾斜变形值，5.3.地表任意点a(x，y)沿φ方向的曲率变形为：为：为：为：为：式中：k(x,y)
max
，k(x,y)
min
表示分别为待求点最大、最小曲率变形值；k
x
，k
y
表示分别为待求点沿走向及倾向在主断面投影处迭加后的曲率值；5.4.地表任意点a(x，y)沿φ方向的水平移动值为：
u(x,y)
cm
＝u
x
c
′
y
cosφ
u
+u
y
c
′
x
sinφ
u
φ
u
＝arctg(u
y
c
′
x
/u
x
c
′
y
)u
x
，式中：φ
u
表示为最大水平移动方向与ox轴的夹角，u
y
表示分别为待求点沿走向和倾向在主断面投影点处的水平移动值，对于倾斜方向需加c
y
·
w
cm
·
cotθ。5.5.地表任意点a(x，y)沿方向的水平变形值为：为：为：为：为：式中：ε(x,y)
max
，ε(x,y)
min
表示为待求点最大、最小水平变形值；ε
x
，ε
y
表示为待求点沿走向及倾向在主断面投影处迭加的水平变形值。6.根据权利要求1所述的一种兼顾风积沙流动量的煤矿开采地表沉陷预测方法，其特征在于：根据地表移动变形的概率积分预测值，结合已建立的风积沙厚度变化公式，可以推导出风积沙流动量对地表下沉的影响大小：w
沙
＝v1w
预
+v2i
预
+v3k
预
+v4u
预
+v5ε
预
+b并且可以计算出剔除风积沙流动量影响的煤矿开采地表实际下沉值w
实
：w
实
＝w
预-w
沙
。

技术总结
本发明涉及一种兼顾风积沙流动量的煤矿开采地表沉陷预测方法，属于开采沉陷盆地的地表风积沙流动量测量领域。首先在采动工作面地表不同位置均匀间隔布设附有测钉的透水隔沙聚酯密目网作为风积沙厚度观测站，采用RTK技术和钢尺量距方法观测，获取不同采动程度下地表沉陷盆地不同位置地表移动变形值和地表风积沙厚度变化值；然后将测量参数进行相关分析和多元线性回归，得到采动影响下风积沙流动量的计算公式；最后利用概率积分法计算相应的地表沉陷预测值，计算出风积沙流动造成的地表沉陷值，计算出剔除风积沙流动量影响的煤矿开采地表实际沉陷量。其实施方便，计算步骤简易，能有效界定鄂尔多斯矿区风积沙运移对地表开采沉陷的影响程度。沉陷的影响程度。沉陷的影响程度。


技术研发人员：郑辉 孟凡贞 李怀展 张连贵 郭广礼 陈延康 查剑锋 宫亚强 李天雨 唐路 元亚菲 霍文奇
受保护的技术使用者：中国矿业大学
技术研发日：2023.06.25
技术公布日：2023/9/23