分段Weibull時間函數開采沉陷動態預計研究

張永勝，胡海峰，蔡音飛，賀福帥

(太原理工大學 礦業工程學院，太原 030024)

地下開采是一個復雜的動態平衡過程，對其進行動態預測不僅可以實時掌握各個點的變形程度及位置，還可為構(建)筑物的保護提供可靠的數據支撐。目前我國使用的預測方法主要是概率積分法，但是概率積分法適用于靜態預計及地表移動變形穩定后的沉陷監測，無法動態化的表述變形過程。

動態預計的函數主要有：Knothe時間函數、Weibull時間函數、雙曲線時間函數、Gompertz時間函數、Logistic時間函數、正態分布時間函數等。我國主要采用Knothe時間函數，而研究發現其不太符合實際地表下沉[1-2]；常占強等[3]提出了分段Knothe時間函數，提高了預測精度；張兵等[4-5]對分段Knothe時間函數進行了優化，解決了關鍵點理論誤差以及終點無法收斂于1的問題，進一步提高了精度；劉玉成等[6-7]針對Knothe時間函數存在的不足，提出了冪指數Knothe時間函數模型；為了解決分段Knothe函數的時間影響系數c等參數為固定值的問題，郭旭煒等[8]對其進行改進，提出了一種新的分段Knothe時間函數；劉玉成等[9]將Weibull時間函數與概率積分法結合，證明Weibull時間函數可以較好地表述測點下沉過程；王小華等[10]將Weibull時間函數與負指數法進行結合，推導出了地表下沉、傾斜、曲率、水平移動、水平變形等公式，并通過實例進行驗證，結果表明兩者具有很好的一致性。在日常實際應用中，筆者發現Weibull時間函數在動態預計過程中存在一些不足：曲線初始期與實際下沉軌跡不一致；關鍵節點函數值與理論值不符。本文使用分段函數思想對該函數優化，構建了一種精度更高、更符合實際的分段Weibull時間函數，以山西某礦3214工作面實測數據為例對優化后的分段Weibull時間函數進行了精度驗證，并與Weibull時間函數、Knothe時間函數和優化后的分段Knothe時間函數[4]進行了對比分析。

1 Weibull時間函數及其分析

1.1 Weibull時間函數

地下開采引起的地表點移動變形是一個復雜的四維變化過程，在Knothe時間函數的基礎上，對t增加一個冪指數k能夠更好地描述此過程，也就是文獻[10]中提到的Weibull時間函數，函數模型如下：

Wt=Wm(1-e-ctk) .

(1)

式中：Wm為沉陷點穩定后的最大沉降量，mm；t為沉陷時間；c、k都是與采空區上覆巖層地質條件有關的參數。

對式(1)求導可得沉陷點下沉速度的函數模型為：

vt=Wmcktk-1e-ctk.

(2)

對式(2)求導可得沉陷點加速度的函數模型為：

at=Wmcke-ctk[(k-1)tk-2-ckt2(k-1)] .

(3)

令式(3)等于0，解出

(4)

結合式(1)-式(4)，做出Weibull時間函數模型的下沉、下沉速度以及下沉加速度的曲線圖，如圖1所示。

圖1 Weibull函數下沉、下沉速度和下沉加速度曲線圖Fig.1 Weibull function subsidence, subsidence velocity, and subsidence acceleration curve

由圖1可知，Weibull函數模型在t=0時刻，沉陷值、下沉速度、加速度都為0；在(0，t0)時間段，沉陷值和下沉速度持續增加，加速度為大于0的正數；在(t0，∞)時間段，沉陷值繼續增加，下沉速度減小至0，加速度為小于0的負數；在t=∞時刻，沉陷值達到最大，下沉速度、加速度都為0.沉陷點的下沉速度由0→vmax→0變化，加速度由0→+amax→0→-amax→0變化，與實際變化過程相符，且大致呈S型分布，證明能夠用于沉陷區的動態預測。

1.2 Weibull時間函數在動態預測中的局限性

利用Weibull時間函數對山西某礦3214工作面實測下沉曲線進行擬合，Weibull函數擬合出的沉陷點下沉過程與實測曲線相差較大，與實際不符，如圖2所示。表明Weibull時間函數模型無法準確表述地表下沉的變化情況。分析其原因如下：

1) 從圖2來看，這個函數的主要問題是初始期沒有斜率變大的階段，這就導致地表沉陷點的預測值提升過快，而實際上當地下開采產生的影響波及到地表點之后，地表點的沉陷值是緩慢增加的，在圖像上則表現為初始斜率較小且斜率逐漸增長，經過一段時間，地表點的下沉速度(圖像斜率)達到最大，再往后逐漸減小到0，直至地表趨于穩定狀態。

2) Weibull時間函數在下沉開始和下沉結束時刻的函數值等于0和1，符合實際情況，但在其他關鍵點(如下沉速度最大時刻)處的函數值與理論值不符。

圖2 Weibull函數預測下沉與實測下沉對比圖Fig.2 Comparison diagram of Weibull function’s predicted subsidence and measured subsidence

2 分段Weibull時間函數

2.1 動態預計分段Weibull時間函數的建立

針對上述Weibull時間函數的局限性，筆者利用常占強等[3]提出的分段函數思想，將沉陷點的下沉過程分為加速下沉和減速下沉兩個階段，假設其在時序上和數值上是相等的，以沉陷點的最大下沉速度時刻τ(以d為單位)為分界點，對這兩個階段分別建模，并結合偏差改正、生長函數模型等手段消除理論誤差，提出了一種新的、精度更高的分段Weibull時間函數：

(5)

式中：t為沉陷時間；τ為沉陷點最大下沉速度時刻；T為沉陷總時間；c為時間影響系數；k為冪指數參數，與上覆巖層地質條件有關。

對式(5)求一階導可得沉陷點下沉速度的函數模型為：

(6)

化簡得：

(7)

對式(5)求二階導可得沉陷點加速度的函數模型為：

(8)

化簡得：

(9)

根據上述函數模型，令Wm=1，繪制出在τ=200 d時分段Weibull時間函數下沉曲線、下沉速度曲線、下沉加速度曲線如圖3、4所示。圖3中下沉曲線表明，曲線初始期斜率逐漸增加，衰退期斜率逐漸減小，沉陷值處于0～1之間，斜率最大處即最大下沉速度時刻函數值為0.5Wm，形態與理想時間函數基本一致；下沉速度曲線圖中，下沉速度經歷了0→vmax→0變化過程，且在活躍期速度增加迅速；從圖4下沉加速度曲線可以看出，下沉加速度趨勢符合實際情況，在t=0到t=τ時間段，從0增加到+amax，在t=τ到t=T時間段，從-amax減小到0，但分布情況與實際不符，初始期和衰退期加速度為0的過程較長，加速度變化較快，縮減了變形時間，表明其地表點移動具有突發性。綜上所述，分段Weibull時間函數模型符合理想時間函數具有的特征，能夠真實客觀地反映地形切割強烈的黃土丘陵地區地表點的下沉過程。

圖3 修正后下沉、下沉速度曲線Fig.3 Corrected subsidence and subsidence velocity curves

圖4 修正后下沉加速度曲線Fig.4 Corrected subsidence acceleration curve

2.2 最大下沉速度時刻τ

最大下沉速度時刻τ是指以首次觀測時刻或工作面開始回采時刻為起點，沉陷范圍內地表點達到最大下沉速度時所經歷的時間差。

以每一期沉陷值與其最大值的比值(即下沉系數)為縱坐標，沉陷時間為橫坐標，繪制下沉系數關于沉陷時間的函數曲線。依據分段思想，沉陷值在τ時刻即下沉速度最大時刻沉陷值為0.5Wm，由此可確定函數圖像上縱坐標為0.5所對應的橫坐標值即為τ.圖5中，以山西某礦33K2工作面A18號點為例，依據上述方法求取的τ值為113 d，但由于該礦區屬于山西丘陵黃土地區，實際上下沉速度最大時刻τ為118 d，兩者略有差別，使用兩個τ值分別進行預測，計算它們的相對誤差，發現只有0.45%的偏差，對預測結果影響較小，表明使用上述方法求取τ值符合沉陷區動態預計的精度要求。

圖5 τ的確定Fig.5 Determination of τ

求出本文研究礦區3214，3215等5個工作面各個地表點對應的τ值，分析發現其與距開切眼距離x(即采動程度)存在多項式關系，對τ和x進行回歸分析，得到相關性較好的回歸方程：

τ=3.022 7×10-4×x2+0.178 7x+110.628 9 .

(10)

利用上述回歸方程，即可預測出地表任意點最大下沉速度時刻τ值。

圖6 τ關于x的回歸方程Fig.6 Regression equation of τ with respect to x

2.3 參數確定

研究表明，在深厚比較大(一般大于30)的情況下，地下開采引起的地表動態移動變形具有明顯的連續性，且開采沉陷預計模型參數與地質采礦條件有著密不可分的關系[11-14]，具有各向異性，因此在時間影響系數c和冪指數參數k兩個參數調控下，動態預計模型的精度會更高。

2.3.1時間影響系數c

上覆巖塊硬度越大，深厚比越小，巖塊的碎漲系數越小，相應的時間影響系數c也就越??；反之越大。開采動態預計模型具有相似推廣性，因此在相同地質條件下，時間影響系數c可用經驗類比法確定，也可用以下公式計算[15-16]：

(11)

(12)

式中：P為覆巖參數；H為開采深度，m；h為松散層厚度，m；γ為覆巖平均容重，kg/m3；E為覆巖平均彈性模量，MPa；σt為覆巖平均抗拉強度，MPa；υ為覆巖平均泊松比。

2.3.2冪指數參數k

假定Wm為單位值，令τ=200 d，繪制不同k對應的下沉曲線見圖7.由圖7可知，在時間影響系數c取定值的情況下，參數k控制預計點的變形路徑，與起止大小無關，顯然，這樣更能客觀地描述開采動態過程。依據經驗類比或上述公式確定好時間影響系數c后，用最小二乘法擬合冪指數參數k，如果存在足夠的實測數據，也可通過類比確定k值。

圖7 參數k對下沉曲線的影響Fig.7 Influence of parameter k on the subsidence curve

3 實例驗證

山西某礦3214工作面地面位于丘陵地帶，地形復雜，切割強烈，大部分被黃土覆蓋，且黃土沖溝發育，基巖在較大的黃土沖溝中出露，測區范圍內北高南低最大高差達120 m，東西方向地勢平緩，3214工作面走向長1 214 m，傾向長156 m，煤層平均厚度1.98 m，平均傾角5.5 °，平均開采深度535.5 m，開采速度為2.5 m/d，開采方法為走向長壁后退式采煤方法。

3.1 單點實例驗證

為驗證優化后Weibull時間函數模型的適用性，以3214工作面走向觀測線上的最大下沉點C39號點為研究對象，利用上述4種時間函數對其隨時間變化情況分別進行動態預計，預計時通過回歸方程求得C39號點對應的最大下沉速度時刻τ=230 d，模型參數的選取采用最小二乘擬合法，確定出分段Weibull時間函數的參數c=0.000 4，k=1.67，Weibull時間函數的參數c=0.000 4，k=1.67，Knothe時間函數的參數c=0.02.預計結果如圖8所示。

圖8 最大下沉點C39號點預測對比Fig.8 Prediction comparison of the maximum subsidence point C39

從圖8可以看出，分段Weibull時間函數預測的下沉曲線與實測下沉曲線最為接近，優化后的分段Knothe時間函數次之，未優化的Weibull時間函數和Knothe時間函數擬合度最差，為了定量分析4種函數模型預測的準確性，利用式(13)、(14)計算均方誤差m和相對誤差f，將其作為精度評定的指標，見表1.

(13)

(14)

式中：d為地面點預測沉陷值與實測沉陷值之差；n為觀測次數；Wm為沉陷點穩定后的最大沉降量，mm.

表1 不同時間函數精度對比Table 1 Accuracy comparison of different time functions

通過表1可知，采用分段Knothe時間函數、Weibull、Knothe時間函數進行預測的相對誤差分別為4.95%、42.57%、52.51%，而采用分段Weibull時間函數進行預測，相對誤差為1.48%，由此可見，使用分段Weibull時間函數進行動態預計具有更高的適用性。

3.2 主斷面實例驗證

為了進一步驗證優化后Weibull時間函數模型的可靠性，采用分段Weibull時間函數對3214工作面走向主斷面進行預計，對比分析實際測量值與預計結果，如圖9所示。為確定預計精度，分別計算11期測量值和預測值之間的均方誤差和相對誤差，結果見表2.對走向主斷面而言，最大均方誤差為7.61 mm，最大相對誤差為5.25%，表明預計曲線與實測曲線高度一致，且最終預計值與實測最大下沉值相差很小，更直接地證實了本文所提出的分段Weibull時間函數模型的可靠性。

圖9 走向主斷面動態預計對比Fig.9 Dynamic comparison of strike main section

表2 走向主斷面精度分析Table 2 Accuracy analysis of strike main section

4 結論

本文在Weibull時間函數的基礎上，建立了分段Weibull時間函數模型，結合工程實例驗證了本模型的可靠性，并得出以下結論：

1) Weibull時間函數預測開采沉陷，其下沉值、下沉速度、下沉加速度等物理量的變化比較符合實際動態下沉過程，但是存在一定的局限性：下沉曲線在初始期斜率變化太快，在最大下沉速度時刻處的函數值不為0.5.

2) 提出分段Weibull時間函數，優化后函數的下沉、下沉速度和下沉加速度曲線形態符合理想時間函數特征，利用圖解法求取最大下沉速度時刻τ，建立了τ與距開切眼距離x的回歸方程，并對巖性參數c和k的確定方法進行了探討。

3) 利用優化后的Weibull時間函數對山西某礦3214工作面的最大下沉點進行預測，將預測結果與Weibull、Knothe、優化后分段Knothe時間函數的預測值及實測值進行對比，發現優化后Weibull時間函數的預測精度分別提高了41.09%、51.03%和3.47%；對比分析了3214工作面走向主斷面的11次預測結果和實測值，最大相對誤差為5.25%，表明本文函數模型能夠很好地預測地表動態下沉，進一步證明了分段Weibull時間函數模型的適用性和可靠性。