喻成林 ,張宏貞 ,范洪冬
(1.徐州市自然資源和規(guī)劃局,江蘇 徐州 221008;2.中國礦業(yè)大學(xué) 環(huán)境與測繪學(xué)院,江蘇 徐州 221116)
地下資源開采后,破壞其原有平衡狀態(tài),經(jīng)過覆巖傳遞到地表,是1 個(gè)復(fù)雜的動(dòng)態(tài)移動(dòng)過程,涉及地質(zhì)、采礦、空間等多方面因素,許多學(xué)者在開采沉陷動(dòng)態(tài)變形原理、模型、精度等方面進(jìn)行了研究,并取得了大量成果。開采沉陷動(dòng)態(tài)預(yù)測模型主要有動(dòng)態(tài)力學(xué)模型和時(shí)間函數(shù)模型2 類,目前主要采用時(shí)間函數(shù)構(gòu)建開采沉陷的動(dòng)態(tài)變形過程;常用的時(shí)間函數(shù)有Knothe 時(shí)間函數(shù)模型、雙曲線時(shí)間函數(shù)模型、Gompertz 時(shí)間函數(shù)模型、logistic 曲線時(shí)間函數(shù)模型及 Weibull 曲線時(shí)間函數(shù)模型等。
王悅漢等[1]將動(dòng)態(tài)移動(dòng)過程分為頂板初次垮落前、頂板垮落未充滿采空區(qū)、頂板垮落已充滿采空區(qū)、表土移動(dòng)4 個(gè)階段,采動(dòng)巖體分為連續(xù)介質(zhì)、擬連續(xù)介質(zhì)、非連續(xù)介質(zhì)3 種介質(zhì),分階段建立了動(dòng)態(tài)力學(xué)預(yù)測模型;吳侃等[2]提出了基于時(shí)間序列分析的動(dòng)態(tài)預(yù)測模型;崔希民等[3]建立基于Knothe 時(shí)間函數(shù)的概率積分預(yù)計(jì)方法,根據(jù)臨界開采尺寸確定時(shí)間影響系數(shù);張兵等[4-5]建立了1 種新的優(yōu)化分段Knothe 時(shí)間函數(shù),并提出了基于實(shí)測數(shù)據(jù)的反算時(shí)間函數(shù)對比求參法、基于概率積分參數(shù)的直接計(jì)算法的2 種優(yōu)化分段Knothe 時(shí)間函數(shù)參數(shù)求取方法;劉玉成等[6-7]分析了下沉動(dòng)態(tài)過程的Knothe 時(shí)間函數(shù)模型的不足,提出了冪指數(shù)的Knothe 時(shí)間函數(shù),分析了常見地表沉陷時(shí)間函數(shù)模型的下沉、下沉速度、下沉加速度與時(shí)間關(guān)系,認(rèn)為Weibull 曲線時(shí)間函數(shù)模型能完整地描述地表沉陷的動(dòng)態(tài)過程;孫闖等[8]建立了考慮松散層下沉特點(diǎn)的雙因素Knothe 時(shí)間函數(shù);王軍保等[9]借鑒巖石流變力學(xué)中非定常流變模型,構(gòu)建了變時(shí)間影響系數(shù)的Knothe 的地表動(dòng)態(tài)沉陷預(yù)測模型;許國勝等[10]回歸分析得到地表動(dòng)態(tài)移動(dòng)變形參數(shù)與地質(zhì)及開采技術(shù)參數(shù)之間的關(guān)系;李春意等[11-12]進(jìn)行了基于Logistic 時(shí)間函數(shù)、正態(tài)分布時(shí)間函數(shù)的動(dòng)態(tài)沉陷預(yù)測;張永勝等[13]提出了以地表點(diǎn)最大下沉速度時(shí)刻為分界點(diǎn),結(jié)合偏差改正、生長函數(shù)模型的1 種分段Weibull時(shí)間函數(shù);劉東海[14]根據(jù)常村觀測資料,指出Weibull 時(shí)間函數(shù)參數(shù)不固定,兩參數(shù)變化規(guī)律為距邊界變化類偏態(tài)分布;JAROSZ 等[15]考慮開采產(chǎn)生的空間收斂,建立了雙參數(shù)Sroka-Schober 時(shí)間函數(shù);KOWALSKI[16]在Knothe 和Sroka-Schober雙參數(shù)時(shí)間函數(shù)基礎(chǔ)上,提出了廣義時(shí)間函數(shù);González Nicieza 等[17]在Knothe 時(shí)間函數(shù)基礎(chǔ)上,引入正態(tài)分布函數(shù),建立了正態(tài)分布時(shí)間函數(shù)。
綜上,總體認(rèn)為地表動(dòng)態(tài)下沉量與時(shí)間的關(guān)系曲線呈“S”形曲線,下沉速度為0→Vmax→0 過程,下沉加速度為0→+amax→-amax→0 過程。
工作面上方的覆巖及地表在井下資源開采后,經(jīng)歷彎曲、垮落、再次壓縮等過程,但在工作面不同位置的巖層完整厚度、位移、受力等經(jīng)歷過程不同,根據(jù)工作面回采、覆巖受力與垮落情況,將工作面回采結(jié)束后的地表分為3 個(gè)部分,工作面動(dòng)態(tài)分析分區(qū)圖如圖1。
圖1 工作面動(dòng)態(tài)分析分區(qū)圖Fig.1 Partition diagram of dynamic analysis of working face
1)豎向整體受壓區(qū)(Ⅰ區(qū))。處于工作面外側(cè)煤柱區(qū)域,工作面開采邊界至開采影響邊界之間。此部分區(qū)域在整個(gè)回采期間,覆巖經(jīng)歷了彎曲,但沒有垮落、再次壓縮等過程,在豎向上是完整的。在整個(gè)開采期間豎向方向整體是受壓狀態(tài)。
2)中間區(qū)(Ⅱ區(qū))。處于工作面內(nèi)測靠近煤柱區(qū)域,工作面開采邊界至覆巖動(dòng)態(tài)完全區(qū)邊界(覆巖垮落角確定邊界與工作面開采邊界之間)。此部分區(qū)域在整個(gè)回采期間,覆巖經(jīng)歷了彎曲、垮落、非完全再次壓縮過程等部分或全部過程,主要是由于巖梁的存在導(dǎo)致垮落非充分及垮落巖石受力非充分,在豎向上部分破裂,且離工作面邊界(煤柱)越遠(yuǎn),破裂高度越大,即完整覆巖厚度減小。在整個(gè)回采期間主要經(jīng)歷如下幾個(gè)階段:①受到影響到工作面推到下方前,在豎直方向上整體是受到壓縮變形,豎向受力增大;②推過后至覆巖垮落前,在豎直方向上整體是受到拉伸變形、豎向受力減?。虎劭迓浜笾烈苿?dòng)穩(wěn)定,在豎直方向上整體是受到壓縮變形,豎向受力增大至穩(wěn)定。垮落巖石部分再次受到壓縮變形。
3)動(dòng)態(tài)完全區(qū)(Ⅲ區(qū))。處于工作面內(nèi)測遠(yuǎn)離煤柱區(qū)域,覆巖動(dòng)態(tài)完全區(qū)(覆巖垮落角確定邊界至采空區(qū)中心)。此部分區(qū)域在整個(gè)回采期間,覆巖經(jīng)歷了彎曲、垮落、完整再次壓縮過程,在豎向上部分破裂,破裂高度一致。在整個(gè)回采期間主要經(jīng)歷如下幾個(gè)階段:①受到影響到工作面推到下方(垮落前),在豎直方向上整體是受到壓縮變形,豎向受力增大;②到再次受到壓縮前,在豎直方向上整體是受到拉伸變形,豎向受力減?。虎鄣揭苿?dòng)穩(wěn)定,在豎直方向上整體是受到壓縮變形,豎向受力增大至穩(wěn)定。
根據(jù)開采沉陷巖梁相關(guān)理論,地表下沉可簡化為各種類型巖梁彎曲,其彎曲下沉值與剛度(彈性模量與慣性矩之積)成反比,其中慣性矩與巖梁高度的三次方成正比。結(jié)合圖1,在同種情況下,上覆巖層剛度、慣性矩變化規(guī)律為:Ⅰ區(qū)(整體一致,最大)>Ⅱ區(qū)(與離開切眼距離成反比)>Ⅲ區(qū)(整體一致,最小)。
目前以時(shí)間函數(shù)為基礎(chǔ)的各種地表動(dòng)態(tài)沉陷預(yù)測模型參數(shù)主要包括:最大變形值、時(shí)間函數(shù)(時(shí)間函數(shù)不同,其參數(shù)不同),各類參數(shù)與覆巖、開采等因素有關(guān),體現(xiàn)其差異,但目前時(shí)間函數(shù)都采用相同參數(shù),即在工作面不同位置(煤柱內(nèi)外)的時(shí)間函數(shù)參數(shù)都一樣,這與處在工作面不同位置的受力過程、變形過程、巖層破壞情況及巖石力學(xué)性質(zhì)等明顯不符。
Weibull 時(shí)間函數(shù)f(t)數(shù)學(xué)表達(dá)式為:
式中:k、m為參數(shù),正數(shù);t為時(shí)間。
據(jù)此建立開采沉陷過程中地表下沉?xí)r間函數(shù)、下沉速度時(shí)間函數(shù)、下沉加速度時(shí)間函數(shù),開采沉陷動(dòng)態(tài)下沉過程曲線如圖2。
圖2 Weibull 時(shí)間函數(shù)下沉、速度、加速度變化圖Fig.2 Weibull time function subsidence, velocity and acceleration variation diagram
式中:W0為該點(diǎn)下沉的最終值。
下沉速度時(shí)間函數(shù)v(t):
根據(jù)前面分析可知:工作面不同位置在回采期間,其覆巖移動(dòng)、破壞、受力、巖梁剛度等不同,則其動(dòng)態(tài)時(shí)間函數(shù)參數(shù)應(yīng)不同,m、k決定了Weibull 時(shí)間函數(shù)的形態(tài)變化過程,因此m、k在工作面不同位置其存在差異,以垮落區(qū)邊界為原點(diǎn)、指向煤柱方向?yàn)檎鶕?jù)前述巖梁彎曲特征簡述及相關(guān)巖石力學(xué)理論,建立其與地面點(diǎn)距離垮落區(qū)域邊界距離S的函數(shù),工作面不同位置m、k變化情況圖如圖3。
圖3 工作面不同位置m、k 變化情況圖Fig.3 Changes of m and k at different positions of working face
在煤柱側(cè)區(qū)域(Ⅰ區(qū))m、k為:
式中:a1、b1、c1、a2、b2為系數(shù);S0為煤柱邊界距離;H為工作面采高。
從m、k計(jì)算公式及圖3 中可以看出:
1)在垮落區(qū)外側(cè)為負(fù)指數(shù)關(guān)系,以工作面邊界為界,此區(qū)域分為2 部分;Ⅰ區(qū)(煤柱范圍)m變化較小,基本為常數(shù);Ⅱ區(qū)(中間區(qū))m變化較大,在垮落區(qū)內(nèi)側(cè),m為一常數(shù)。
2)以工作面邊界、垮落區(qū)邊界為界,此區(qū)域分為3 部分;Ⅰ區(qū)(煤柱范圍)k為常數(shù);Ⅱ區(qū)(中間區(qū))m為線性變化;在垮落區(qū)內(nèi)側(cè),k為一常數(shù)。
由于不同位置上覆巖層經(jīng)歷的受力、破壞、變形過程不同,根據(jù)開采沉陷相關(guān)研究,可采用如下方法確定垮落范圍d。① 方法1:采用覆巖力學(xué)理論,根據(jù)巖梁的相關(guān)破裂理論確定垮落邊界;② 方法2:采用概率積分法拐點(diǎn)偏移距;③ 方法3:采用頂板覆巖垮落角φ確定,d=H/tan φ。
某工作面采用長壁一次性采全高采煤法回采。走向?yàn)?74 m,傾斜寬為220.7 m。工作面開采深度為470~515 m。根據(jù)走向觀測數(shù)據(jù),采用上面的計(jì)算模型,獲得的反演模型系數(shù)見表1。
表1 反演模型系數(shù)Table 1 Coefficients of inversion model
根據(jù)工作面開采情況,以工作面內(nèi)側(cè)90 m 為垮落邊界,以此確定各點(diǎn)距離邊界的S值,根據(jù)表1 反演動(dòng)態(tài)時(shí)間函數(shù) φ(t)如下:
根據(jù)工作面開采情況及參數(shù)確定情況,分別以開采結(jié)束后實(shí)測結(jié)果作為最終值、概率積分法預(yù)測值作為最終值,進(jìn)行開采第523 d 的動(dòng)態(tài)預(yù)測,實(shí)測結(jié)果、概率積分法預(yù)測值作為最終值的動(dòng)態(tài)預(yù)測與實(shí)測對比如圖4(圖中距離以工作面邊界為原點(diǎn),采空區(qū)方向?yàn)檎?。
圖4 擬合與實(shí)測對比圖(523 d)Fig.4 Comparison between fitting and actual measurement (523 d)
以實(shí)測結(jié)果作為最終值:動(dòng)態(tài)預(yù)測中誤差為110 mm,為最大下沉值的3.5%;以概率積分法預(yù)測結(jié)果作為最終值:動(dòng)態(tài)預(yù)測中誤差為130 mm,為最大下沉值的4.2%。
1)分析了工作面上方不同位置的覆巖在開采過程中在豎向上的受力、變形、破壞等變化過程,根據(jù)煤柱、垮落情況將覆巖與地表劃分為豎向整體受壓區(qū)、中間區(qū)、動(dòng)態(tài)完全區(qū)3 個(gè)區(qū)域。
2)以Weibull 時(shí)間函數(shù)動(dòng)態(tài)預(yù)測模型為例,分析了m、k2 個(gè)參數(shù)的空間分布規(guī)律及特地,構(gòu)建了2 個(gè)參數(shù)的分區(qū)數(shù)學(xué)模型及選取準(zhǔn)則。
3)分別以實(shí)測結(jié)果、概率積分法預(yù)測結(jié)果作為工作面開采結(jié)束后的最終值,依據(jù)本文動(dòng)態(tài)預(yù)測模型及參數(shù)確定方法,動(dòng)態(tài)預(yù)測中誤差分別為110、130 mm,為最大下沉值的3.5%、4.2%,證明本模型的可靠性。