石曉宇, 張燕海, 楊可明, 姚樹一, 王劍
(1.中國礦業(yè)大學(北京) 地球科學與測繪工程學院, 北京 100083;2.淮北礦業(yè)(集團)股份有限公司 通防地測部, 安徽 淮北 235000)
煤炭是目前我國主要的供給能源[1]。地下煤炭資源的大面積開采使得采區(qū)周圍巖層原始應(yīng)力平衡被打破[2],引起巖層移動和變形,不僅直接造成礦區(qū)地表沉陷,還會帶來一系列潛在的災(zāi)害問題。對礦區(qū)地表沉陷規(guī)律的研究是采礦工程中一個重大課題,其中如何更高效、更系統(tǒng)、更安全地對礦區(qū)地表進行沉陷監(jiān)測與預(yù)計備受關(guān)注。
地下采煤引起的巖層移動和變形是一個復(fù)雜的時空過程,傳統(tǒng)的地表沉陷預(yù)計方法(如概率積分法、剖面函數(shù)法等)僅能獲得地表沉陷最終穩(wěn)定后的空間變化,不能較好地反映開采沉陷隨時間變化的動態(tài)過程[3]。眾多學者綜合考慮因巖體本身復(fù)雜性而難以全面研究采煤工作面上覆巖層移動及地表沉陷動態(tài)過程的情況,以局部觀測點為切入點,建立沉陷區(qū)單點沉陷量時間序列模型,從而研究礦區(qū)地表動態(tài)沉陷規(guī)律[4-6]。GM(1,1)與灰色Verhulst模型為灰色系統(tǒng)理論中2個基礎(chǔ)的時間序列模型,可在小樣本、貧信息、不確定性等劣勢情況下,鑒別系統(tǒng)因素間發(fā)展趨勢的相異程度,尋找系統(tǒng)變化規(guī)律,并對事物發(fā)展趨勢作出模糊性的短期、中長期預(yù)測[7-8]。GM(1,1)是用一階微分方程對1個變量建立的模型[9],主要解決具有指數(shù)變化規(guī)律的生成序列,只能描述單調(diào)變化過程;灰色Verhulst模型是由德國生物學家Verhulst在研究物種繁殖趨勢時提出的一階非線性動態(tài)模型[10],其基本思想是生物繁殖數(shù)量受到周圍環(huán)境的影響,繁殖曲線會以“S”形呈現(xiàn),最終達到一個平衡穩(wěn)定的狀態(tài)。將這2種灰色模型(Grey Model,GM)應(yīng)用于礦區(qū)地表沉陷預(yù)計,建立沉陷量時間序列模型,可較好地揭示地表單點連續(xù)動態(tài)變化的過程,實現(xiàn)地表沉陷動態(tài)預(yù)計。然而,現(xiàn)有文獻只針對2種GM在礦區(qū)沉陷預(yù)計方面的應(yīng)用進行研究,并沒有結(jié)合模型自身特點,分析和比較二者在礦區(qū)地表沉陷預(yù)計中的適用性。此外,目前礦區(qū)地表單點動態(tài)沉陷預(yù)計方法主要基于傳統(tǒng)的水準測量數(shù)據(jù)[11],監(jiān)測方法單一,成本高,且觀測點易破壞,不能保證地表形變信息的實時性。合成孔徑雷達差分干涉測量(Differential Interferometry Synthetic Aperture Radar,D-InSAR)技術(shù)通過雷達干涉圖的差分獲得地表形變信息,以高精度、全天候、連續(xù)空間覆蓋、實時性等優(yōu)勢成為新興的形變監(jiān)測手段,被廣泛應(yīng)用于礦區(qū)地表沉陷監(jiān)測中[12-15]。將D-InSAR技術(shù)與GM相結(jié)合,避免監(jiān)測與預(yù)計獨立進行,可較好地形成系統(tǒng)化及模塊化流程,實現(xiàn)礦區(qū)地表沉陷監(jiān)測與動態(tài)預(yù)計一體化。
本文以淮北礦業(yè)(集團)股份有限公司袁店二礦7221工作面為試驗區(qū)域,融合D-InSAR技術(shù)和GM(1,1)與灰色Verhulst模型,通過建立沉陷量與時間的GM進行地表沉陷動態(tài)預(yù)計,實現(xiàn)監(jiān)測與預(yù)計一體化,并針對2種GM的特點,結(jié)合預(yù)計結(jié)果分析和比較其在礦區(qū)地表沉陷預(yù)計中的適用性。
袁店二礦7221工作面上方地表主要覆蓋物為村莊和農(nóng)田,四鄰無其他采掘活動。該工作面為俯斜開采,工作面寬約110 m,走向長約544 m,煤層傾角為10~20°,基巖厚度為132~162 m,采高為4.5 m,于2017-12-06開始回采,2018-06-03停采。為了控制地表采動損害在建筑物等安全保護等級范圍之內(nèi),分別在地表公路或建筑體上設(shè)置水準測量與變形監(jiān)測點,以監(jiān)測地表及重要建筑物變形情況。但該區(qū)域內(nèi)人工溝渠縱橫分布,監(jiān)測點稀少、空間配置不合理且常被人為破壞,無法準確對重要建筑物進行損害評估。因此,有必要在工作面開采過程中對該區(qū)域進行地表沉陷D-InSAR監(jiān)測。
為了盡可能減少D-InSAR技術(shù)相位失相干對地表沉陷監(jiān)測結(jié)果精度的影響,選取冬季和初春時段(該時段植被蕭疏且地物波譜特征無明顯變化)的9景哨兵衛(wèi)星(Sentinel-1A)雷達影像作為D-InSAR監(jiān)測的數(shù)據(jù)源。影像時間為2017-11-04—2018-02-08,相鄰時間間隔(12 d)的影像可以形成8個干涉對,見表1,其中VV為同向垂直極化方式,IW為干涉測量寬幅模式。外部數(shù)字高程模型選用分辨率為90 m的SRTM(Shuttle Radar Topography Mission,航天飛機雷達地形測繪使命)數(shù)據(jù)。
表1 Sentinel-1A雷達影像干涉對
根據(jù)表1數(shù)據(jù),利用D-InSAR技術(shù)二軌法獲取8組7221工作面時間序列沉陷分布。為方便分析,以2017-11-04雷達影像為參考對象,將各組沉陷量依次累加,即可得到各組時間序列沉陷分布。由于2017-11-04—28地形無明顯變化,所以僅給出后6組時間序列沉陷分布,如圖1所示??煽闯鲈诨夭沙跗?2017-12-06—10),沉陷量為-10~0 mm。隨著工作面推進,進入沉陷發(fā)展階段(圖1(b)—圖1(f)),在2017-12-22開始形成下沉盆地,與7221工作面推進方向一致。下沉盆地自西南向東北方向擴張,沉陷量和平面范圍隨之增大。至2018-02-08,最大沉陷量為-85 mm,位于7221工作面切眼附近,離注漿孔(注3、注4、注5,位于工作面東北部)較近的村莊未發(fā)生明顯沉陷,說明注漿充填開采工藝有很好的局域減沉效果。
(a) 2017-11-04—2017-12-10 (b) 2017-11-04—2017-12-22 (c) 2017-11-04—2018-01-03
(d) 2017-11-04—2018-01-15 (e) 2017-11-04—2018-01-27 (f) 2017-11-04—2018-02-08
圖1 7221工作面時間序列沉陷分布
Fig.1 Subsidence distribution of 7221 working face in time series
為更好地對D-InSAR監(jiān)測結(jié)果進行定量研究,提取下沉盆地中沿7221工作面傾向l1和走向l2(l1,l2位置如圖2所示)的沉陷剖面線,如圖3所示。l1長807 m,方向為自西向東,設(shè)其最西端為零點;l2長1 089 m,方向為自南向北,設(shè)其最南端為零點。從圖3(a)可看出,2017-11-04—2018-02-08,l1最大沉陷量為-81 mm,2017-11-28—2018-01-03沉陷剖面線變陡,沉陷量顯著增大,后期進入沉陷衰減階段,開始微小沉陷。從圖3(b)可看出,2017-11-04—2018-02-08,l2最大沉陷量為-84 mm,2017-11-28—2018-01-03沉陷速度較快,動態(tài)最大沉陷量不斷增加,2018-01-03—02-08沉陷速度開始衰減,但平面影響范圍不斷增大。對比圖3(a)和圖3(b)可看出,開采對傾向沉陷影響范圍小于走向,且傾向沉陷階段發(fā)育早于走向,下沉盆地近似出現(xiàn)平底狀態(tài)。
基于D-InSAR監(jiān)測獲取的2017-11-16—2018-02-08 7221工作面地表沉陷量,從各期沉陷分布(以表1中干涉對組合為參考, 2017-11-16—28為第1期觀測數(shù)據(jù),2017-11-16—12-10為第2期觀測數(shù)據(jù),以此類推至2018-02-08,共7期數(shù)據(jù))中提取43個水準點的像元信息,并繪制沉陷量曲線,如圖4所示。選取H2,H3,E2,G5,F(xiàn)4為觀測點,對7期數(shù)據(jù)建立GM(1,1)與灰色Verhulst模型方程,獲得前7期沉陷量擬合值,并預(yù)計第8,9期(2018-02-20,2018-03-04)沉陷量。1—7期各觀測點沉陷實際值與2種GM擬合值對比如圖5所示。
圖2 l1,l2位置
(a) 傾向l1
(b) 走向l2
GM精度檢驗一般分為事中檢驗和事后檢驗[16]。事中檢驗包含殘差檢驗、關(guān)聯(lián)度檢驗和后驗差檢驗。殘差檢驗屬于算術(shù)檢驗,定義為模型值與實際值之間的差,是一種較為客觀的逐點檢驗法;關(guān)聯(lián)度檢驗屬于幾何檢驗,表征2個事物之間的關(guān)聯(lián)程度,關(guān)聯(lián)值越大(越接近1)則相關(guān)性越好,模型等級越高;后驗差檢驗屬于統(tǒng)計檢驗,是按照均方差比值C和小誤差概率P來評估灰色模型優(yōu)劣,C越小越好,相反P越大越好。事后檢驗即預(yù)測檢驗,通常采用殘差檢驗和相對誤差檢驗來了解預(yù)測誤差。本文對于GM預(yù)計結(jié)果的精度檢驗分為2步:采用殘差檢驗、關(guān)聯(lián)度檢驗和后驗差檢驗對GM(1,1)和灰色Verhulst模型進行擬合精度檢驗,結(jié)果見表2;對預(yù)計的第8,9期沉陷值進行殘差和相對誤差檢驗,結(jié)果見表3和表4。
(a) 水準點位置
(b) 各水準點沉陷量曲線
(a) H2
(b) H3
(c) E2
(d) G5
(e) F4
表2 GM(1,1)和灰色Verhulst模型擬合精度檢驗Table 2 Fitting precision checking of GM(1,1) and grey Verhulst model
表3 第8期沉陷量預(yù)計精度檢驗Table 3 Prediction precision checking of subsidence value in phase 8
(1) 從圖5(a)、圖5(b)可看出,實際沉陷量曲線后期呈平底飽和狀態(tài),此時GM(1,1)和灰色Verhulst模型前期(1,2,3期)擬合結(jié)果差別不大,但后期(4—7期)灰色Verhulst模型較GM(1,1)收斂快,擬合結(jié)果與實際沉陷量吻合較好。圖5(c)中實際沉陷量曲線后期包含非單調(diào)擺動過程,使用GM(1,1)擬合誤差較大。圖5(d)、圖5(e)中實際沉陷量曲線呈近似單峰型,GM(1,1)擬合曲線更接近實際值,且后期GM(1,1)和灰色Verhulst模型擬合結(jié)果相差不大。
表4 第9期沉陷量預(yù)計精度檢驗Table 4 Prediction precision checking of subsidence value in phase 9
(2) 從表2可看出,觀測點H2,H3使用灰色Verhulst模型擬合結(jié)果明顯優(yōu)于GM(1,1),平均殘差分別減小了約26.5%,36.7%。根據(jù)后驗差精度等級(表5),觀測點H2的GM(1,1)擬合結(jié)果只能達到合格,而其灰色Verhulst模型擬合結(jié)果精度等級達到良好;觀測點G5采用GM(1,1)擬合時平均殘差較灰色Verhulst模型減小了1.249 mm,關(guān)聯(lián)度提高了約3.2%;觀測點E2,F(xiàn)4采用2種GM擬合的效果相近,精度等級均達到良好,且平均殘差之差均小于1 mm。
表5 后驗差精度等級Table 5 Precision level of posterior error
(3) 從表3、表4可看出,采用GM(1,1)和灰色Verhulst模型預(yù)計第8,9期5個觀測點的沉陷量時精度與表2中擬合精度較一致,證明可根據(jù)各點擬合精度(殘差、關(guān)聯(lián)度和后驗差)選擇沉陷預(yù)計適用模型。GM(1,1)對觀測點H3的預(yù)計誤差最大,最大殘差絕對值為56 mm,顯然該點沉陷量曲線更適宜選用灰色Verhulst模型來描述?;疑玍erhulst模型對第8期觀測點G5的預(yù)計精度最低,最大相對誤差為21%,而使用GM(1,1)預(yù)計時殘差減少了6 mm,且相對誤差減小了5%,預(yù)計效果更佳。
可見對于單點沉陷曲線后期呈平底飽和狀態(tài)的觀測點(如H2,H3),地表進入沉陷衰退階段時選用灰色Verhulst模型進行沉陷動態(tài)預(yù)計,精度明顯優(yōu)于GM(1,1),最大相對誤差減少了55%;地表進入沉陷發(fā)展階段后,對于單點沉陷量曲線呈近似單峰型的觀測點(如G5,F(xiàn)4),選用GM(1,1)預(yù)計沉陷量更合適,最大相對誤差可減少35%。
(1) 井下采煤地表移動是一個先加速后減速的過程,即隨著工作面推進和開采時間推移,地表點會經(jīng)歷開始移動、劇烈移動、緩慢移動及停止移動4個階段,其移動軌跡是一個隨機變化的過程。
(2) 礦山開采沉陷開始至活躍前期,地表單點沉陷量曲線呈近似單峰型,采用GM(1,1)進行短期預(yù)計效果更佳;當?shù)V區(qū)地表單點沉陷量曲線呈平底飽和狀態(tài),沉陷進入衰退階段,此時采用灰色Verhulst模型預(yù)計的結(jié)果與實際值誤差較小,適用于中長期預(yù)計。
(3) GM(1,1)和灰色Verhulst模型較傳統(tǒng)地表沉陷預(yù)計方法而言,是假設(shè)地表沉陷歷程符合已知時間序列所描述的動態(tài)過程,通過先驗信息來研究或預(yù)測未知領(lǐng)域的時間函數(shù)模型,而開采引起的地表沉陷是一個復(fù)雜的時空變化過程,僅依靠時間序列模型并不能實現(xiàn)準確預(yù)計,后續(xù)會將時間序列模型與某種常規(guī)靜態(tài)方法(如概率積分法)相結(jié)合,建立更精確的地表任意點動態(tài)沉陷預(yù)計模型。