煤矿地下开采引起的地表移动变形复杂多变,沉陷预计为采空区上方建筑物、铁路、高速公路、大坝等建构筑物采后保护和治理提供了理论基础与决策依据。静态预计提供的是地表移动变形最终结果,往往并不能满足实践需求,地面点动态变化规律在研究地表建构筑物的保护中也尤为重要。国内外学者对此进行了大量地探索和研究,提出了不同的时间函数和地表动态沉陷预计理论[1]。
常用于动态开采沉陷预计的时间函数主要有以下几种:波兰学者KNOTHE于1952年提出了Knothe时间函数[2-5];1982年SROKA和 SCHOBER在Knothe 时间函数的基础上提出了Sroka-Schober时间函数也称双参时间函数[6];1999年KOWALSKI[7]提出了广义时间函数;2007年GONZALEZ 等[2]研究后认为可以使用正态分布函数作为时间函数等。在我国由采矿引起的地面沉降变形的动态预计大多数使用的是Knothe时间函数,国内很多学者对Knothe时间函数进行了深入地研究。崔希民等[8]提出了时间函数影响系数的确定方法,结合概率积分法,提出基于Knothe时间函数的概率积分动态预计方法,推导了动态预计公式;刘玉成等[9]提出了幂指数Knothe时间函数,定义了岩性参数,并证明了该时间函数符合地表下沉的动态过程;针对工作面匀速开采过程中及开采后,胡海峰等[10]推导了基于Knothe动态沉降假设和Litwiniszyn随机介质理论的动态沉降预测模型;针对Knothe时间函数时空不完备性,常占强等[11]提出了分段Knothe时间函数,提高了动态预计精度;张兵等[12]以分段Knothe时间函数为基础进行了改进,拓宽了该时间函数的适用范围,并且提高了预计精度;廉旭刚等[13]提出了通过倾向微元化达到顾及主要影响半径随采深变化的方法,使得缓倾斜煤层使用Knothe模型得以实现;郭博婷等[14]结合山区滑移影响函数,提出了基于Knothe模型的山区开采地表动态下沉预计方法。2007年GONZALEZ 等[2]将正态分布时间函数应用到阿斯图里亚斯煤矿区的地表动态下沉预测,然后和Knothe时间函数以及双曲时间函数的预测结果进行对比分析,说明了该时间函数具有更高的精度;张兵等[15]给出了正态分布时间函数参数确定方法;李春意等[1]讨论了正态分布时间函数的完备性,给出了相应的沉陷预计模型。笔者深入分析和研究了正态分布时间函数,通过对比分析得出该时间函数存在理论缺陷:关键节点处函数值存在理论误差;时间参数c的选取具有局限性,只适用于c>2的情况。通过分析理论缺陷产生的原因,采用整体偏差修正和生长函数模型优化的方法,构建了一个理论上更加完善,适用范围更广,同时精度有所提高的时间函数。以大同矿区某矿工作面开采为例,对优化后时间函数进行了对比和验证。
动态预计时间函数是一个比例函数,其函数值代表地表某一点某一时刻的下沉量占该点最大下沉量的比例,设时间函数为F(t)∈[0,1],有
W(t)=W0F(t)
(1)
式中:W0为地表最大下沉量。
大量实测资料表明煤矿开采所引起地表下沉的实际过程:工作面开始推进,经过一定的时间后,上覆岩层的冒落、断裂和弯曲才会传递到地表,引起地表的下沉和水平移动等移动变形,地表的下沉经历0→W0(W0为地表最大下沉量,下文中V、a代表速度与加速度),下沉速度经历0→Vmax→0,下沉加速度经历0→+amax→0→-amax→0这样一个过程[16-17],这个过程可以近似认为是对称的,且地表点下沉速度最大时,其下沉量大致等于该点最大下沉量的一半[18],即:下沉总时间为T,当t=T/2时,
由图1可以看出,正态分布原函数及正态分布密度函数与正态分布一阶导数在曲线特征上符合地表点下沉、下沉速度与下沉加速度的动态变化趋势,并且具有很好的时空完备性。
μ—样本均值;σ—中误差
图1 正态分布函数
Fig.1 Normal distribution funclion
一般情况下,地表下沉都具有一定的延迟时间,即工作面开始推进至地表出现移动变形的时间,用tp表示;地表下沉具有一定的持续时间,即地表开始出现移动变形至地表移动变形结束的时间,用tc表示[1]。由正态分布函数的性质与图1可知,直线x=μ为正态分布函数的对称轴,此时地表下沉达到最大下沉量的1/2,下沉时间为总时间的1/2,即
μ=tp+tc/2
(2)
引入时间参数c(c的取值由实际开采情况决定),取t=μ-cσ为地表下沉的开始时间,取t=μ+cσ为地表下沉的结束时间,即
(3)
由式(2)和式(3)建立的均值μ和标准差σ与延迟时间和持续时间的关系,代入正态分布原函数中即可得出正态分布时间函数为
(4)
化简得
(5)
通过正态分布时间函数与理想时间函数对比分析可知该函数存在适用范围局限和预计精度随时间参数c的减小而降低等问题。
1)对于正态分布密度函数来言,其积分域只有是(-,+)时,原函数的值域才会是[0,1],也就意味着使用该函数作为时间函数的原型时,只有下沉时间无限长才会达到最大下沉,与实际情况相悖。
图2 正态分布时间函数与理想时间函数对比
Fig.2 Comparison of normal distribution time
function and ideal time function
图3 正态分布时间函数与时间参数关系
Fig.3 Relationships of normal distribution time
function and time parameters
2)由图2可知,正态分布时间函数只满足开始下沉时刻函数值等于0的函数要求,其他各关键节点处的函数值均不等于理论值。
3)由图3可知,时间参数c取值越小,时间曲线变化越平缓,终点函数值越远离于理论值1。导致该时间函数的时间参数的取值范围受到限制,只有c>2时为其合理取值。
为了探究上节提到问题的根源所在,对正态分布时间函数的形成过程进行分析。对于该时间函数的3个参数(tp、tc、c)来说,延迟时间tp与持续时间tc的取值是通过实测得到的,而时间参数c的取值直接决定正态分布时间函数的质量高低。由图4可知,正态分布时间函数是在正态分布密度函数的基础上选取了有限积分区域对其积分得到的,虽然函数曲线保持了原有形态,但是发生了整体的偏移,导致该时间函数最终无法收敛于1,各个关键节点处的函数值不等于理论值。
图4 正态分布密度函数
Fig.4 Normal distribution density function
时间参数的值越小,时间函数的理论误差越大,是因为正态分布时间函数的函数值是正态分布密度曲线在积分域为[tp,tp+tc]的积分值,即两时刻垂线(图4中虚线)、时间横轴和密度曲线所围成的面积,两时刻垂线两侧的面积即为时间函数的理论误差。由式3可知,持续时间tc是确定的,时间参数c与密度函数标准差σ成反比,c取值越大,σ取值越小,密度函数越向中间靠拢,有效积分域越大,函数理论误差越小;相反函数理论误差越大。通过计算得出正态分布时间函数的起点误差为
(6)
终点误差为
(7)
为了使时间函数的起点值为0,所以起点误差和终点误差都叠加表现在了时间函数内部及终点处,且终点误差最大。
通过上节中对正态分布时间函数模型的分析得出误差的大小是由有效积分域的大小决定的,通过计算得出理论误差值为2ΔF(t);函数曲线如图5所示,本次模型修正的思想是:
1)使用起点误差ΔF(t)作为改正数,对时间函数采用整体偏差修正,使该函数恢复到理论模型的位置。该方法可以有效地使理论误差重新对称分配到函数中点的两侧,改变了原函数理论误差从t=tp开始累积至函数终点导致函数后段误差过大的缺陷,有利于进一步对该时间函数的修正。
图5 修正函数示意
Fig.5 Schematic diagram of revised function
令作为改正项,将原时间函数进行整体偏差修正,得
(8)
化简得
(9)
2)改正后正态分布时间函数的理论误差是对称分布的,对称点为(tp+tc/2,1/2),对称点处理论误差为0,其他位置理论误差不均匀分布,且随至对称点距离的增大而增加,随时间参数取值减小而增加,因此需构建1个与时间参数呈反比、与距对称点距离呈正比的误差修正函数模型,减小修正后正态分布时间函数模型的理论误差。构建生长函数模型为
(10)
将式(10)化简为简化生长函数模型为
(11)
将式(11)代入改正后函数F1(t),得出优化后正态分布时间函数模型为
(tp≤t≤tp+tc)
(12)
图6 优化后正态分布时间函数与时间参数关系
Fig.6 Relationships between time function and time
parameter of normal distribution after optimization
对比图3与图6可知,正态分布时间函数经上述优化,在关键时间节点t=tp,F(t)=0;t=tp+tc/2,F(t)=1/2;t=tp+tc,F(t)=1,符合理论要求,且没有改变原本函数的时空完备性。由表1知,优化后的时间函数解决了原本时间参数只能取c>2的问题,并消除了原时间函数在关键节点处的理论误差。
表1 原时间函数与改进时间函数理论误差对比
Table 1 Comparison of theoretical error between original time function and improved time function
时间参数原正态分布时间函数函数终值理论偏差相对理论误差改进正态分布时间函数函数终值理论偏差相对理论误差相对优化率/%c=10.6830.31731.7%10046.41c=20.9550.0454.5%1004.71c=30.9970.0030.30%1000.30c=40.9990.0010.10%1000.10
为了验证优化后正态分布时间函数的适用性与可靠性,选取大同甲煤矿8KK0工作面A观测线A25号点和乙煤矿8KK3工作面2号观测线D17号点作为实例验证对象,分别使用正态分布时间函数与优化后正态分布时间函数进行动态预计。根据实测数据可知甲煤矿A25点于2017年12月1日开始沉降,2018年7月16日沉降趋于稳定,乙煤矿D17点于2015年12月11号开始沉降,2016年6月3日沉降趋于稳定,时间参数使用最小二乘最佳逼近的方式计算求得,反算间隔为0.1,通过计算得出A25点的预计参数为:tp=0,tc=227,c=1.7,D17点的预计参数为:tp=0,tc=176,c=2。由图7 可以看出,优化后时间函数模型相比原时间函数模型更加吻合实际地表下沉情况,实测数据均匀的分布在优化后时间函数两侧。根据表2中预计结果可以看出,优化后时间函数模型相比于原时间函数在预计精度方面有所提高,本实例中A25点标准误差M提高了55.4 mm,相对误差F提高了5.11%,相对优化率最高达到11.3%;D17点标准误差提高了17.3 mm,相对误差提高了3.23%,相对优化率最高达到9.9%。
图7 优化前后时间函数模型预计结果与实测数据对比
Fig.7 Comparison of predicted and measured data of time
function model before and after optimization
表2 A25和D17号点优化前后时间函数模型预计结果对比
Table 2 Comparison of predicted results of time function models before and after optimization of Point A25 and D17
点号观测时间(年-月-日)下沉量/mm实测优化前预计优化后预计原函数预计误差/mm原函数预计精度优化后预计误差/mm优化后预计精度优化率/%A252017-12-0100002019-01-16-165.8-115.229-133.90450.5712018-03-03-367.0-340.342-377.69426.6582018-04-08-640.0-559.808-611.77680.1922018-05-10-823.0-736.590-801.54986.4102018-06-12-990.0-867.069-945.427122.931M=82.433 mmF=7.97%031.896-10.69428.22421.45144.573M=27.004 mmF=2.87%011.304.408.107.907.90D172015-12-1100002016-01-03-31.5-24.723-27.8386.7772016-01-26-68.0-77.457-83.688-9.4572016-02-13-129.0-142.147-150.816-13.1472016-03-02-228.0-222.586-233.6945.4142016-03-23-348.0-320.773-334.72627.2272016-04-08-412.0-386.315-402.43525.6852016-04-21-462.0-428.63-446.51133.372016-05-03-480.0-457.945-477.45122.0552016-05-19-516.0-483.751-505.42532.2492016-06-03-521.0-497.294-521.0023.706M=27.786 mmF=5.33%03.662-15.688-21.816-5.69413.2749.56515.4892.54910.5750M=10.475 mmF=2.11%09.90-9.20-6.70-0.104.003.903.904.104.204.60
1)通过对比和分析总结出正态分布时间函数在理论上和实际应用中存在的问题:①在各关键节点处函数值不等于理论值,函数模型存在理论误差;②正态分布时间函数对时间参数的选取具有局限性,只适用于c>2的情况,导致一些情况下无法使用正态分布时间函数进行动态预计。
2)通过对正态分布时间函数分析,总结出该时间函数存在问题的原因是时间参数变换影响该时间函数的形态,导致时间函数有效积分域亏损随时间参数减小而增加,进而使该时间函数的理论误差变大。
3)结合问题产生的原因,确立了优化思想。首先通过对函数模型以整体偏差修正的方式将函数模型恢复到理论位置;再构建生长函数对修正后模型进行进一步优化,最终优化后函数模型在各关键节点处的函数值都等于理论值,更加符合实际情况,各节点间预计误差也得到了一定的控制,从而提高了函数整体的预计精度。并且优化后时间参数的选取不受限制,拓宽了正态分布时间函数的应用范围。
4)使用大同甲煤矿8KK0工作面A观测线A25点乙煤矿8KK3工作面2观测线D17点的实测数据验证优化后正态分布时间函数模型可靠性。通过预计结果与实测数据对比,使用优化后正态分布时间函数进行动态预计比原时间函数模型的预计结果精度有了一定程度的提高。对于矿区地表重要建筑物损坏的实时预防和维护具有实际的指导意义。
[1] 李春意,高永格,崔希民.基于正态分布时间函数地表动态沉陷预测研究[J].岩土力学,2016,37(S1):108-116.
LI Chunyi,GAO Yongge,CUI Ximin.Progressive subsidence prediction of ground surface based on the normal distribution time function[J].Rock and Soil Mechanics,2016,37(S1):108-116.
[2] GONZALEZ N C,ALVAREZ F M I,MENENDEZ D A,et al.The influence of time on subsidence in the Central Asturian Coalfield[J].Bulletin of Engineering Geology & the Environment,2007,66(3):319-329.
[3] 许国胜,李德海,侯得峰,等.厚松散层下开采地表动态移动变形规律实测及预测研究[J].岩土力学,2016,37(7):2056-2062.
XU Guosheng,LI Dehai,HOU Defeng,et al.Measurement and prediction of the transient surface movement and deformation[J].Rock and Soil Mechanics,2016,37(7):2056-2062.
[4] NIE L,WANG H,XU Y,et al.A new prediction model for mining subsidence deformation:the arc tangent function model[J].Natural Hazards,2015,75(3):2185-2198.
[5] KNOTHE S.Time influence on a formation of a subsidence surface[J].Archiwum Gómictwai Hutnictwa,Kraków,1952,1(1):1
[6] JAROSZ A,KARMIS M,SROKA A .Subsidence development with time :experiences from longwall operations in the Appalachian coalfield[J].International Journal of Mining & Geological Engineering,1990,8(3):261-273.
[7] KOWALAKI A.Surface subsidence and rate of its increments based on measurements and theory[J].Archives of Mining Science,2001,46(4):391-406.
[8] 崔希民,缪协兴,赵英利,等.论地表移动过程的时间函数[J].煤炭学报,1999,24(5):453-456.
CUI Ximin,MIAO Xiexing,ZHAO Yingli,et al.Discussion on the time function of time dependent surface movement[J].Journal of China Coal Society,1999,24(5):453-455.
[9] 刘玉成,曹树刚,刘延保.改进的Konthe地表沉陷时间函数模型[J].测绘科学,2009,34(5):16-17,31.
LIU Yucheng,CAO Shugang,LIU Yanbao.The improved Knothe time function for surface subsidence[J].Science of Surveying and Mapping,2009,34(5):16-17,31.
[10] HU H F,LIAN X G,Li Y .Application study on the dynamic prediction model for determining the mining subsidence[J].Journal of Engineering Science & Technology Review,2006,9(2):80-85.
[11] 常占强,王金庄.关于地表点下沉时间函数的研究:改进的克诺特时间函数[J].岩石力学与工程学报,2003,22(9):1496-1496.
CHANG Zhanqiang,WANG Jinzhuang.Study on time function of subsidence:the improved Knothe time function[J].Chinese Journal of Rock Mechanics and Engineering,2003,22(9):1496-1499.
[12] 张 兵,崔希民.开采沉陷动态预计的分段Knothe时间函数模型优化[J].岩土力学,2017,38(2):541-548.
ZHANG Bing,CUI Ximin.Optimization of segmented Knothe time unction model for dynamic prediction of mining subsidence[J].Rock and Soil Mechanics,2017,38(2):541-548.
[13] 廉旭刚,戴华阳,胡海峰.煤矿开采动态地表下沉预计算法[J].辽宁工程技术大学学报:自然科学版,2016,35(10):1032-1036.
LIAN Xugang,DAI Huayang,HU Haifeng.Prediction algorithm of dynamic surface subsidence induced by underground coal mine[J].Journal of Liaoning Technical University:Natural Science 2016,35(10):1032-1036
[14] 郭博婷,胡海峰,廉旭刚.基于Knothe模型的山区开采地表动态下沉预计方法[J].煤矿安全,2016,47(3):190-193.
GUO Boting,HU Haifeng,LIAN Xugang.Prediction method of surface subsidence in mountain area based on Knothe Model[J].Safety in Coal Mines,2016,47(3):190-193.
[15] 张 兵,崔希民,胡青峰.开采沉陷动态预计的正态分布时间函数模型研究[J].煤炭科学技术,2016,44(4):140-145.
ZHANG Bing,CUI Ximin,HU Qingfeng.Study on normal distributed time function model to dynamically predict mining subsidence[J].Coal Science and Technology,2016,44(4):140-145.
[16] 崔希民,陈立武.沉陷大变形动态监测与力学分析[M].北京:煤炭工业出版社,2004.
[17] 王军保,刘新荣,李 鹏,等.MMF模型在采空区地表沉降预测中的应用[J].煤炭学报,2012,37(3):411-415.
WANG Junbao,LIU Xinrong,LI Peng,et al.Study on prediction of surface subsidence in mined-out region with the MMF model[J].Journal of China Coal Society,2012,37(3):411-415.
[18] 何国清,杨 伦,凌赓娣,等.矿山开采沉陷学[M].徐州:中国矿业大学出版社,1995.