耦合拐點偏移修正的Boltzmann函數預計模型及其參數反演方法

江克貴，王 磊，池深深，魏 濤，蔣 創

( 1. 安徽理工大學 測繪學院，安徽 淮南 232001；2. 中國礦業大學 環境與測繪學院，江蘇 徐州 221116 )

我國能源資源稟賦具有“富煤、貧油、少氣”的特點。據預測，即使到2030年煤炭資源仍占中國能源需求的50%以上[1]。大規模高強度的煤炭資源開發會引發一系列礦山地質環境災害問題，如含水層破壞、山體滑坡、地表塌陷和建構筑物損毀等。開采沉陷預計理論在評估潛在地質災害，分析礦區沉降機理方面具有至關重要的作用，因而開采沉陷預計及防治的相關技術與理論一直是我國礦山領域的研究熱點，基于實測資料求取精準、可靠的開采沉陷預計參數也是預計模型應用過程中的難 點。

常規概率積分法存在以下幾個問題：沉降變形預計曲線邊緣收斂過快，與實際變形情況有較大差別；不能較好地適用于復雜多樣地質采礦條件下的地表沉陷預計，不具有廣適應性；聯合時間函數的動態沉陷預計模型仍然存在一些問題。解決上述問題可以從2個途徑入手：① 構建精準廣適應的開采沉陷預計方法；② 提高預計參數的準確性和可靠性。對于途徑①，相關學者提出了新的開采沉陷預計方法或進行一定程度的改進。一類是針對不同地質采礦條件而提出的沉陷預計方法，戴華陽[2]等提出了地表移動的矢量預計法，為解決大傾向煤層開采沉陷問題提供了預計依據；郭文兵[3]提出了全采工作面疊加預計方法，提高了深部大采寬條帶開采預計精度；文獻[4-6]提出了形式各異的開采沉陷預計方法，在一定程度上解決了多種不同地質采礦條件下的開采沉陷預計問題。另一類是針對改進時間函數的動態開采沉陷預計方法，根據地表移動變形規律和采動程度，文獻[7]建立了一種新的Knothe時間函數求參模型，提高了地表動態下沉預測的準確性和可靠性；文獻[8]采用整體偏移修正和生長函數模型相結合的方法對時間函數進行優化，進一步提高了地表下沉動態預計的精度。近年來，一些學者基于Boltzmann函數提出了新的開采沉陷預計模型[9]，通過分析函數曲線形態和參數特點發現，與概率積分法相比，Boltzmann函數在曲線邊緣部分收斂緩慢，且其預計模型參數與概率積分法模型參數可以相互轉化；另外，研究表明，該模型能夠廣適應于不同開采條件下的地表沉陷預計[10]。鑒于此，在Boltzmann函數的基礎上，筆者借鑒概率積分法模型的推導經驗，通過理論分析，探討拐點偏移修正方法[11]，構建了一種耦合拐點偏移修正的Boltzmann函數預計模型。

針對途徑②，大量學者把智能優化算法引入 到開采沉陷預計參數求取中來，相對于傳統方 法[12-14]，避免了求參速度慢、對初值依賴性強、易陷入局部最優值、抗差能力差等缺陷，極大地提高了求取參數的準確性和可靠性。其中比較典型的有模矢法[15]、遺傳算法[16]、粒子群算法[17-18]等。筆者擬將FWA( Fireworks Algorithm )應用于反演預計模型參 數，提出了一種耦合拐點偏移修正的Boltzmann函數模型參數反演方法，并對該方法的抗差性能進行討論。

1 模型原理

玻爾茲曼( Boltzmann )函數模型在應用實踐中得到了廣泛的發展[19-20]。分析發現，Boltzmann函數可以生成一個S形曲線，跟半無限開采地表變形盆地主斷面的下沉曲線形態相似。對非充分開采、充分開采以及超充分開采的實測結果分析表明，Boltzmann函數可以較好地預計不同開采程度下的地表開采沉陷。基于Boltzmann函數的半無限地表變形盆地主斷面的下沉預計公式為

式中，W0=mq cos α，其中，m為煤層采厚，q為下沉系數，α為煤層傾角；R為新的主要影響半徑，與主要影響半徑r的物理意義相同，其主要影響下沉曲線的分布形態；S為拐點偏移距。

統計分析發現[14]，R≈r/4.13，其中，r=H/tan β，H為煤層采深，tan β為主要影響角正切。并且當x→∞時，有W( x )→W0；當x=s時，W( x )=W0/2；這些結果都符合開采沉陷一般規律。

對式( 1 )求微分，可得Boltzmann函數預計模型的單元影響函數，即

根據概率積分法的推導經驗和彈性力學知識，可得單元水平移動函數為

對單元水平移動函數進行積分，可得地表沿x方向的水平移動為

令b'=B/R，統計分析可得b'=4.13b，則

式中，B為常數；b為水平移動系數。

煤壁附近采空區上方頂板會出現懸頂或過度破壞垮落，其作用效果相當于工作面采空區兩端煤壁向內或向外分別平移了S的距離。文獻[11]表 明，概率積分法將拐點偏移距的影響全部加在了終采線和傾向上山處；由此，根據實測數據擬合求取的開采沉陷預計參數會因為擬合度低而出現偏差較大的現象。針對此問題，在Boltzmann函數預計模型的基礎上對拐點偏移進行修正，走向和傾向主斷面有限開采下沉和水平移動計算公式分別為

式中，

其中，D3和D1分別為工作面走向和傾向開采長度； S1，S2，S3，S4分別為工作面傾向下山、傾向上山、開切眼和終采線拐點偏移距；θ為開采影響傳播角。

井下開采尺寸為D3×D1的工作面引起地面任意一點( x，y )沿φ方向( x軸的正向逆時針到指定方向的角值 )的下沉和水平移動分別為

綜上所述，模型預計參數為

地質采礦參數為

式中，?為工作面傾向方位角，直接參與工作面坐標轉換的過程。

構建模型的概化如式( 9 )所示，工作面地質采礦參數和模型預計參數作為輸入變量控制地表相應預計點的形變；相反，當地表監測點的形變獲取后，在工作面地質采礦參數已知的情況下，求取精確可靠的模型預計參數顯得尤為重要；同時構建模型屬于高度非線性模型。因此，本文引入一種FWA。

2 模型參數反演方法

2.1 FWA

FWA是受到夜空中煙花爆炸的啟發而提出的一種具有局部搜索能力和全局搜索能力自調節機制的群體智能算法[21]。FWA的性能取決于爆炸算子( 即爆炸產生火花操作 )、變異算子( 高斯變異產生火花操作 )、映射規則( 映射到指定范圍操作 )和選擇策略( 選擇下一代火花操作 )[22-23]。

( 1 ) 爆炸算子

爆炸算子是FWA的核心組成部分，它由爆炸數目、爆炸半徑、位移操作3部分組成。假設煙花種群大小為N，對于煙花xi，其爆炸數目和爆炸半徑的計算公式分別為

式中，Di為第i個煙花爆炸產生的火花數目；D為常 數，用來調整產生的爆炸火花數目；f( xi)為第i個煙花的適應度值；ymax為當前種群中最大適應度值，ymax=max( f( xi) )；ξ為一個機器最小量，用來避免除零操作；Ai為第i個煙花的爆炸半徑；A為常數，用來調整爆炸半徑的大小；ymin為當前種群中最小適應度值，ymin=min( f( xi) )。

為了避免適應度較好的煙花位置產生過多的爆炸火花，同時適應度較差的煙花位置產生過少的爆炸火花，對爆炸的火花數目進行如下限制：

式中，a，b為常數；round( )為根據四舍五入原則的取整函數。

位移操作是對煙花的每一維進行位移，即

式中，rand( 0，Ai)表示在幅度Ai內生成的均勻隨機 數；xik為第i個煙花在第k個維度上的位置信息。

( 2 ) 變異算子

為了進一步增加種群的多樣性，FWA在隨機選擇的煙花中執行高斯變異操作，即

式中，Gaussian( 1，1 )表示高斯分布。

( 3 ) 映射規則

爆炸算子和變異算子產生的火花可能會超出邊界，利用映射規則將超出域值的火花映射到一個新的位置。

( 4 ) 選擇策略

在候選者集合中( 包括煙花、爆炸火花和高斯變異火花 )，適應度最好的個體會首先被確定性地選擇到下一代作為煙花，剩下的后代煙花采用輪盤賭法進行選擇。候選者xi被選擇的概率計算公式 為

式中，R( xi)為當前個體xi與其他個體的距離之和； d( xi，xj)表示任意2個個體xi和xj之間的歐氏距離；xj∈K表示第 j個位置屬于集合K。

2.2 方法構建

綜合地表移動變形、預計模型以及FWA原理，提出了一種耦合拐點偏移修正的Boltzmann函數模型參數反演方法( BFWA-MSPP )，具體步驟如下：

( 1 ) 數據準備。首先，利用地表移動觀測站觀測目標區域地表進入穩定期后的下沉和水平移動量，獲取工作面地質采礦條件參數；然后，選取適當的FWA參數，確定開采沉陷預計參數取值空間，進而生成煙花初始種群。

( 2 ) 適應度評價。煙花初始種群用作迭代種 群；基于構建的Boltzmann函數預計模型，結合獲取的地質采礦條件，第i組開采沉陷參數預計目標區域任意點 j的下沉和水平移動分別為 Wpi,j和 Upi,j；假設地表點 j的實測下沉和水平移動值分別為 Wmj和Umj，以預計值和實測值之差絕對值的和為準則進行適應度評價，適應度的構造如式( 17 )所示；然后判斷是否滿足迭代終止條件，如果滿足，輸出最優解，否則執行步驟( 3 )。

( 3 ) 生成火花。根據煙花迭代種群的適應度，結合FWA參數，計算每個煙花的爆炸半徑與爆炸火花數目；利用式( 13 )對各個煙花進行位移偏移，生成多個爆炸火花；另外，在煙花種群中，利用變異算子生成多個變異火花。

( 4 ) 映射及后代選擇。利用映射規則把超出預計參數范圍的火花映射到可行域范圍內。根據選擇策略，在候選者集合中選擇下一代。

( 5 ) 循環執行步驟( 2 )～( 4 )，直到滿足終止條件，輸出最優開采沉陷預計參數，模型參數反演方法技術 路線如圖1所示。

圖1 BFWA-MSPP技術路線 Fig. 1 Technology roadmap of BFWA-MSPP

3 仿真試驗

3.1 數據模擬

以淮南礦區煤系地層為背景，模擬工作面上覆巖層巖性為中硬，工作面地質采礦條件為：煤層采高m=3 m，煤層采深H=400 m，煤層傾角α=5°，工作面傾向方位角?=0°，工作面開采尺寸D3×D1=800 m×200 m。分別沿工作面走向和傾向主斷面 模擬布設了2條監測線E和N，測點間距為30 m，共 布設88個監測點，模擬的工作面和布設的監測線 如圖2所示。基于構建模型和模擬的開采沉陷預計參數( 表1 )預計開采沉陷穩定后的下沉和水平移 動。

3.2 反演試驗

根據模擬的地質采礦條件和預計的形變數據，利用本文構建的BFWA-MSPP進行開采沉陷預計參數反演試驗，求參結果見表1。

從表1可以看出，反演參數q，tan β，b，θ的相對誤差最大不超過1.2%，中誤差最大不超過1.2；反演拐點偏移距S的相對誤差和中誤差分別控制在6.2%和6.5 m以內；此外，BFWA-MSPP擬合下沉和水平移動誤差為-11.6～5.7 mm，中誤差為±6.1 mm。試驗結果表明，構建的方法可以準確可靠地求取開采沉陷預計參數。

圖2 模擬工作面和監測線示意 Fig. 2 Schematic diagram of simulated the working face and the observation lines

表1 仿真試驗BFWA-MSPP求參結果 Table 1 Solved results of BFWA-MSPP for simulation experiment

3.3 抗差試驗

在實際測量中，觀測值會不可避免地摻入一些隨機誤差和粗差[24-25]，隨機誤差和粗差的存在勢必影響求參結果的質量。為了探討BFWA-MSPP抵抗觀測誤差的能力，基于3.1節的模擬數據，設計試驗如下：在所有監測點的下沉值和水平移動值中分別加入均值為0、標準差為30 mm和10 mm的高斯噪 聲；在此基礎上，將所有監測值再加入-200 ～200 mm( 服從均勻分布 )的粗差；然后，以包含測量誤差的數據為觀測值，利用BFWA-MSPP進行求參，求參結果見表2。

表2 BFWA-MSPP抗差求參結果 Table 2 Anti-error result of BFWA-MSPP

從表2可以看出，在所有觀測點加入觀測誤差后，反演q，tan β，b，θ的相對誤差不超過13%，反演 拐點偏移距S的相對誤差控制在32%以內；就擬合效果來看( 圖3 )，在含有觀測誤差干擾的情況下，擬合曲線和模擬觀測曲線較吻合。需要說明的是，本 文為了極限測試算法的穩定性加入了較大的誤差，但在工程實踐中很少有如此大的觀測誤差( 通過精密測量儀器和數據處理方法，觀測誤差可以得到控制 )。試驗結果表明，BFWA-MSPP可以有效抵抗觀測誤差的影響，在監測實施困難區域，可適當降低地表形變監測精度要求。

圖3 BFWA-MSPP抗差擬合下沉和水平移動效果 Fig. 3 Effect diagram of fitting the subsidence and the horizontal movement by anti-error of BFWA-MSPP

4 工程應用

4.1 礦區概況

以淮南礦區顧橋煤礦1414( 1 )工作面為例，該工作面采用綜合機械化采煤，一次采全高，全部垮落法管理頂板。工作面地質采礦條件為：平均采高3 m，工作面傾角5°，平均采深735 m，工作面開采尺寸為2 115.433 m×241.210 m，走向為超充分采動，傾向為非充分采動，總體為非充分采動。地表移動觀測站共布設6個控制點和145個監測點，點間距 為30/60 m，傾向觀測線長度為1 500 m，走向觀測 線長度為3 480 m，工作面和觀測線的布設如圖4所 示。

圖4 工作面和觀測線實地布設 Fig. 4 Field layout of the working face and the observation lines

4.2 試驗設計及結果分析

分別利用BFWA-MSPP和PIM-MSPP對1414( 1 )工作面地表移動觀測站最后一期觀測數據( 地表進入穩定期 )進行開采沉陷預計參數反演，反演結果見表3。

表3 工程應用BFWA-MSPP和PIM-MSPP求參結果比較 Table 3 Comparison between the engineering applications of BFWA-MSPP and PIM-MSPP

從表3中可以看出，2種方法反演參數的差值主要集中在下沉系數q和拐點偏移距S處，就差值所占總體比例來說，下沉系數的差值百分比超過了 10%，拐點偏移距的差值百分比平均約為50%；此 外，開切眼和終采線處的拐點偏移距過大，主要由于工作面走向兩側形成較大砌體梁結構所致。從預計擬合效果來看( 圖5 )，無論是曲線邊緣部分，還是曲線拐點附近，BFWA-MSPP擬合值更加靠近實測值，2種方法預計擬合下沉和水平移動中誤差分別為±83 mm和±108 mm。

5 結 論

( 1 ) 基于Boltzmann函數，融入拐點偏移修正方法，構建了耦合拐點偏移修正的Boltzmann函數預計模型，由于該模型屬于高度非線性模型，通過引入FWA，建立了耦合拐點偏移修正的Boltzmann函數模型參數反演方法( BFWA-MSPP )。

圖5 BFWA-MSPP和PIM-MSPP預計擬合下沉和水平移動比較 Fig. 5 Comparison of BFWA-MSPP and PIM-MSPP fitting the subsidence and the horizontal movement

( 2 ) 仿真試驗結果表明，BFWA-MSPP求參的平均相對誤差約為2.1%，平均中誤差為2.3，BFWAMSPP預計擬合下沉和水平移動中誤差為±6.1 mm。在觀測值包含一定的隨機誤差和粗差后，反演參數的平均相對誤差為13.8%，BFWA-MSPP可以有效抵抗觀測誤差的影響。

( 3 ) 將BFWA-MSPP和PIM-MSPP應用在顧橋 煤礦1414( 1 )工作面，對地表開采穩定后工作面的 開采沉陷預計參數進行了求取，預計擬合下沉和 水平移動中誤差分別為±83 mm 和±108 mm；BFWA-MSPP求參結果為：q=1.07，tan β=1.91，b=0.41，θ=88.32°，S1=-18.89 m，S2=14.74 m，S3=48.89 m，S4=64.42 m。