李靖宇，王 磊，朱尚軍，滕超群，江克貴

(安徽理工大學測繪學院，安徽 淮南 232001)

我國近70%的天然能源來自于煤炭，近年來煤炭資源在能源消費中的占比逐漸降低，但短期內(nèi)仍難以被取代。煤炭資源雖使用廣泛、成本較低，但煤炭開采會造成開采區(qū)巖層發(fā)生移動、變形和破壞，為提前實施有效的預(yù)防措施，必須準確地預(yù)計開采沉陷的影響范圍及變形值大小[1]。我國礦山開采沉陷預(yù)計的主要方法是概率積分法，由于模型簡單、結(jié)果精準等優(yōu)點，得到廣泛應(yīng)用[2]。提高開采沉陷預(yù)測精度[3]主要有兩個途徑：一是提高預(yù)測模型精度；二是提高預(yù)測參數(shù)的準確性與可靠性。在預(yù)測模型正確的情況下，概率積分法模型參數(shù)的精確性會直接影響到預(yù)計的結(jié)果[4]。因此，研究基于地表移動實測資料的概率積分參數(shù)穩(wěn)健反演方法對開采沉陷精準預(yù)計具有重要意義。

概率積分法模型參數(shù)反演歷程[5]中，經(jīng)歷了直接反演、實驗設(shè)計方法反演、優(yōu)化算法反演、智能算法反演的過程。智能算法是現(xiàn)階段概率積分法參數(shù)反演的主要方法，查劍鋒等[6]構(gòu)建了基于遺傳算法的概率積分法預(yù)計參數(shù)反演，抗干擾能力強，在求取參數(shù)方面具有一定的準確性、可靠性和穩(wěn)定性，可以避免計算發(fā)散、陷入局部最優(yōu)解；徐夢強等[7]提出采用粒子群算法反演概率積分法參數(shù)，搜索速度較快、效率高且算法簡單；蘇軍明等[8]利用模擬退火算法反演概率積分法參數(shù)，精度較高，并且計算過程簡單、通用；此后又有專家學者將果蠅算法[9]、量子遺傳算法[10]引入到概率積分法參數(shù)反演中，取得了良好的效果。上述智能算法具有求參精度高、計算速度快、可用于非線性問題等良好的尋優(yōu)特性。為獲得更好的尋優(yōu)特性，吳虎勝等[11]提出了狼群算法。狼群算法作為一種新智能算法，可以運用到復雜的、非線性的、多參數(shù)的尋優(yōu)過程中，求取參數(shù)時能夠保證較高準確性、可靠性和穩(wěn)定性，同時能夠以較大的速率找到最優(yōu)解，可以從多個點出發(fā)進行搜索，在避免點與點之間的影響的同時，防止陷入局部最優(yōu)解。目前，狼群算法被實際應(yīng)用在求解多維背包問題、水電站水庫優(yōu)化調(diào)度以及旅行商等問題中[12-15]，在概率積分法參數(shù)反演中還未得到應(yīng)用。

因此，本文擬將狼群算法引入到概率積分法參數(shù)反演中，構(gòu)建基于WPA的概率積分法模型參數(shù)反演方法，通過模擬實驗系統(tǒng)的討論方法的可靠性，最后開展工程應(yīng)用實驗，并且討論該算法在概率積分法參數(shù)反演中的可靠性和適用性。

1 基于WPA的概率積分法模型參數(shù)反演方法

1.1 概率積分法模型及其參數(shù)體系

如圖1所示，考慮單元開采引起變形的三維情況，整個工作面開采引起坐標為(x,y)的地表任意點A的下沉值的計算過程見式(1)～(4)。

(1)

W0=mqcosα

(2)

W0(x)=

(3)

W0(y)=

(4)

地表任意點A沿φ方向的水平移動值，其計算公式見式(5)和式(6)。

U(x,y,φ)=bri(x,y,φ)

(5)

(6)

式中：φ為x軸的正向逆時針到指定方向的角值；r為主要影響半徑；b為水平移動系數(shù)。

圖1 單元開采地表任意點下沉預(yù)計原理圖Fig.1 Schematic diagram of surface subsidenceprediction at any point of unit mining

綜上可知，概率積分法模型的參數(shù)體系為B=[q,tanβ,b,θ,S1,S2,S3,S4]，參數(shù)體系B的解算精度直接決定開采沉陷的預(yù)計精度。

1.2 基于WPA的概率積分法模型參數(shù)反演方法

1.2.1 WPA理論基礎(chǔ)

人工狼群狩獵區(qū)域是一個N×D的歐氏空間(N為人工狼總數(shù)，D為待求參數(shù)維度)，任意一匹人工狼的位置為Wn=(wn1,wn2,…,wnd,…,wnD)。在人工狼群中，每一匹人工狼都是隨機生成的，隨機生成人工狼群的公式見(7)和式(8)。

Wn=(wn1,wn2,…,wnd,…,wnD),

(1≤n≤N,1≤d≤D)

(7)

wnd=wmin+rand×(wmax-wmin)

(8)

式中：wnd為第n匹人工狼在第d維空間中的位置(即第n組隨機參數(shù)中的第d個參數(shù))；wmax、wmin分別為wnd中的最大值、最小值；rand為-1至1之間均勻分布的隨機數(shù)。

頭狼產(chǎn)生規(guī)則：在人工狼群中，頭狼感知氣味濃度最大，頭狼不斷變化且不參與游行、奔襲、圍攻行為。

游走行為：探狼是除頭狼外、感知氣味濃度較大的人工狼。探狼每次游走向不同方向前進，每次前進后返回原來位置，若探狼前進后感知氣味濃度超過頭狼感知氣味濃度，則更新頭狼并進行奔襲行為，若不是，則繼續(xù)游走直至更新頭狼或達到最大游走次數(shù)；游走結(jié)束后，若探狼感知最大氣味濃度大于當前位置感知氣味濃度，則更新探狼。第k次前進后，探狼的位置計算見式(9)。

(9)

式中：wnd為開始搜尋獵物前，探狼在第d維空間中的位置；stepad為探狼在第d維空間中的游走步長。

奔襲行為：猛狼是距頭狼較近的人工狼。猛狼在奔襲過程中，若感知氣味濃度超過頭狼感知氣味濃度，則更新頭狼且繼續(xù)奔襲行為，直至猛狼與頭狼之間的距離滿足要求。猛狼奔襲后的位置見式(10)。

(10)

圍攻行為：猛狼、探狼接收頭狼的信號后，開始圍攻獵物，若每次攻擊后感知氣味濃度超過原位置感知氣味濃度，則更新該人工狼。第k+1次攻擊后的位置見式(11)。

(11)

狼群更新規(guī)則：圍攻結(jié)束后，淘汰感知氣味濃度較小的人工狼，同時隨機生成相同數(shù)量的人工狼。

1.2.2 參數(shù)反演方法構(gòu)建

設(shè)開采區(qū)任意一點(x,y)的實測下沉值為Wxy實、水平移動值為Uxy實，將狼群算法反演參數(shù)帶入概率積分法模型，得到的該點下沉值為Wxy預(yù)、水平移動值為Uxy預(yù)，以實測值與預(yù)計值之差平方和最小為準則，構(gòu)建適應(yīng)度函數(shù)為式(12)。

F=∑[(Wxy實-Wxy預(yù))2+(Uxy實-Uxy預(yù))2]

(12)

令概率積分法參數(shù)B=[q,tanβ,b,θ,S1,S2,S3,S4]為人工狼狩獵區(qū)域，以適應(yīng)度函數(shù)值作為判斷人工狼感知氣味濃度大小的依據(jù)。在每一次迭代過程中，人工狼群在狩獵區(qū)域中感知氣味濃度并對獵物(概率積分法參數(shù))進行捕獲，達到終止條件后，輸出最優(yōu)概率積分法參數(shù)。 求參主要流程如下所述。

1) 設(shè)置狼群數(shù)量N、最大迭代次數(shù)Dmax、最大游走次數(shù)Tmax等參數(shù)；輸入觀測點實測下沉值Wxy實、水平移動值Uxy實及概率積分法參數(shù)的波動范圍ΔB=[Δq,Δtanβ,Δb,Δθ,ΔS1,ΔS2,ΔS3,ΔS4]，根據(jù)式(7)和式(8)生成人工狼群Bi=[qi,tanβi,bi,θi,S1i,S2i,S3i,S4i](i=1…N)。

2) 以適應(yīng)度函數(shù)值為判斷人工狼感知氣味濃度Yn大小的依據(jù)，適應(yīng)度函數(shù)值越小則人工狼感知氣味濃度越大，選擇感知氣味濃度最大的人工狼為頭狼，記其位置、感知氣味濃度分別為wlead、Ylead。

3) 選擇除頭狼外、感知氣味濃度較大的人工狼為探狼，探狼根據(jù)式(9)進行游走行為，直至達到最大游走次數(shù)Tmax或更新頭狼。

4) 選擇距離頭狼較近的人工狼為猛狼，猛狼根據(jù)式(10)向頭狼奔襲，若奔襲過程中猛狼感知氣味濃度Yn大于Ylead，則更新頭狼，繼續(xù)召喚行為；若Yn小于Ylead，則猛狼繼續(xù)奔襲，直至猛狼與頭狼間距離Dn小于要求距離D。

5) 根據(jù)式(11)，探狼聯(lián)合猛狼對獵物進行圍攻。

6) 對頭狼進行更新，淘汰部分感知氣味濃度較小的人工狼，并根據(jù)式(7)和式(8)隨機生成等數(shù)量人工狼。

7) 判斷是否達到精度要求或最大迭代次數(shù)Dmax，若達到，則輸出頭狼位置，即概率積分法預(yù)計參數(shù)，否則轉(zhuǎn)至第二步驟。

基于WPA的概率積分模型參數(shù)反演方法程序流程圖如圖2所示。

圖2 WPA參數(shù)反演程序流程圖Fig.2 Flow chart of WPA parameter inversion program

2 模擬實驗

2.1 模擬實驗工作面概況

以淮南礦區(qū)煤系地層為背景，設(shè)計模擬實驗工作面參數(shù)為：煤層采厚m=3.0 m，煤層傾角α=5°，走向長度D3=800 m，傾向長度D1=500 m，平均采深H=400 m，工作面采用全部垮落法管理頂板。假設(shè)模擬工作面的概率積分法模型參數(shù)為：下沉系數(shù)q=0.8，主要影響正切角tanβ=2.5，水平移動系數(shù)b=0.25，開采影響傳播角θ=85°，上拐點、下拐點、左拐點、右拐點偏移距分別為S1=60 m、S2=60 m、S3=60 m、S4=60 m。本次模擬實驗在移動盆地內(nèi)設(shè)計了沿走向、傾向的兩條主斷面的移動和變形的觀測線，相鄰兩個監(jiān)測點之間相距20 m，走向觀測線1 180 m，共60個觀測點，傾向觀測線880 m，共45個觀測點。模擬工作面設(shè)計圖如圖3所示。

2.2 狼群算法參數(shù)反演準確性

使用WPA進行概率積分參數(shù)反演，并進行模擬實驗參數(shù)設(shè)計值與模擬實驗參數(shù)反演值的對比，通過參數(shù)平均值、反演參數(shù)相對誤差以及反演參數(shù)中誤差來體現(xiàn)狼群算法反演參數(shù)的準確性。智能算法在一定范圍內(nèi)具有隨機性，為避免偶然誤差的影響，在相同的實驗條件下連續(xù)進行了10次實驗。設(shè)置人工狼數(shù)量為50，最大迭代次數(shù)為50，實驗結(jié)果見表1。

由表1可知,在參數(shù)相對誤差方面，q、tanβ、b、θ的參數(shù)相對誤差均不超過0.3%，S1、S2、S3、S4對求參敏感性較低,但參數(shù)相對誤差均小于2.1%；WPA反演概率比積分法模型參數(shù)具有較高的準確度。在參數(shù)中誤差方面，q、tanβ、b、θ的參數(shù)中誤差均不超過1.03，S1、S2、S3、S4的參數(shù)中誤差均小于8.9 m；WPA反演概率比積分法模型參數(shù)具有較高的精度。綜上所述，WPA反演概率積分法模型參數(shù)具有較好的穩(wěn)健性。

圖3 模擬工作面設(shè)計圖Fig.3 Simulated working face design drawing

表1 WPA反演參數(shù)值與設(shè)計值對比Table 1 Comparison of WPA inversion parameter valueswith design values

為更加直觀地看出下沉與水平移動的擬合情況，得到10次模擬實驗反演參數(shù)均值計算出的下沉值曲線和水平移動值曲線見圖4和圖5。

由圖4和圖5可以看出，下沉值和水平移動值的實測曲線與擬合曲線基本一致，點位實測值與擬合值間的絕對值誤差均不超過5 mm，達到預(yù)期精度。

2.3 基于WPA的概率積分法模型參數(shù)反演方法的可靠性研究

為進一步驗證基于WPA算法的概率積分法模型參數(shù)反演方法的可靠性，對模擬實驗10次反演結(jié)果中的概率積分法模型參數(shù)分別進行波動性分析，對比WPA反演參數(shù)與設(shè)計參數(shù)，結(jié)果見圖6。

圖4 模擬實驗WPA擬合下沉值與實測下沉值對比圖Fig.4 The comparison diagram of the simulated WPAfitted subsidence value and the measuredsubsidence value

圖5 模擬實驗WPA擬合水平移動值與實測水平移動值對比圖Fig.5 The comparison diagram of the simulated WPAfitted horizontal movement value and the measuredhorizontal movement value

圖6 模擬實驗10次反演參數(shù)的波動情況Fig.6 The fluctuation of inversion parameters was simulated for 10 times

從圖6可以看出，WPA反演參數(shù)波動范圍較小，反演參數(shù)較為接近于真值，基于WPA算法的概率積分法模型參數(shù)反演方法具有良好的可靠性。

3 工程實例

3.1 礦區(qū)概況

淮南顧橋礦南二采區(qū)1414(1)工作面是顧橋煤礦南區(qū)的第一個采面，工作面從礦床邊界開始向巷道往回開采，采用機械化掘進，垮落法管理頂板，一次采全高。工作面長約2 120 m，寬約251 m，走向為充分采動，傾向為非充分采動，總體為非充分采動。平均采高為3.0 m，煤層傾角在平均為5°，工作面平均深度為735 m，在距離切眼1 144 m處以及停采線976 m處布置傾向觀測線，傾向線長度為1 500 m，相鄰點之間點間距為30 m，共布設(shè)了3個控制點和50個監(jiān)測點。 在下山方向偏離工作面中心線的距離39 m處設(shè)置走向觀測線，長度為3 480 m，相鄰點之間點間距為30 m/60 m，共布設(shè)了3個控制點和95個監(jiān)測點。

3.2 實驗結(jié)果及分析

以淮南顧橋礦南二采區(qū)1414(1)工作面的實測數(shù)據(jù)為基礎(chǔ)，將WPA反演參數(shù)均值計算所得的下沉值、水平移動值與監(jiān)測點的下沉值、水平移動值進行對比,為避免實驗的偶然性，在相同條件下進行10次實驗，然后分別計算反演參數(shù)平均值、中誤差，結(jié)果見表2。

由表2可知，在參數(shù)中誤差方面，q、tanβ、b、θ的參數(shù)中誤差均不超過0.89，S1、S2、S3、S4的參數(shù)中誤差均不超過4.5 m；WPA反演顧橋南礦1414(1)工作面的概率積分法參數(shù)為：q=0.935 7，tanβ=1.988 1，b=0.422 8，θ=84.138 3°，S1=-9.136 6 m，S2=-21.283 9 m，S3=55.954 4 m，S4=34.111 3 m。

將10次工程實例反演參數(shù)均值帶入概率積分法模型，計算出的下沉值曲線和水平移動值曲線見圖7和圖8。

表2 WPA反演參數(shù)工程實例結(jié)果Table 2 Engineering example results of WPAinversion parameters

圖7 工程實例WPA擬合下沉值與實測下沉值對比圖Fig.7 The comparison diagram of WPA fitted subsidencevalue and the measured subsidence value

圖8 工程實例WPA擬合水平移動值與實測水平移動值對比圖Fig.8 The comparison diagram of WPA fitted horizontalmoving value and measured horizontal moving value

從圖7和圖8可以看出，下沉值和水平移動值的實測曲線與擬合曲線基本一致，點位實測值與擬合值間的絕對值誤差均不超過300 mm。WPA反演概率積分法參數(shù)取均值后，其下沉與水平移動值的擬合中誤差為114.81 mm，其中，下沉值擬合中誤差約為最大下沉值的2.7%，水平移動值擬合中誤差約為最大水平移動值的8.6%，滿足工程應(yīng)用精度[18]。

4 結(jié) 論

1) 模擬實驗結(jié)果表明，WPA反演參數(shù)q、tanβ、b、θ的參數(shù)相對誤差均不超過0.3%，參數(shù)中誤差均不超過1.03；S1、S2、S3、S4對求參敏感性較低，但參數(shù)相對誤差均小于2.1%，S1、S2、S3、S4的參數(shù)中誤差均不超過8.9 m。WPA反演參數(shù)波動范圍較小，反演參數(shù)較為接近于真值。以上證明了WPA在概率積分法模型參數(shù)反演中，具有良好的穩(wěn)健性，且反演參數(shù)具有較高的準確性、精度。

2) 利用本文提出的基于WPA的概率積分法模型參數(shù)反演方法對淮南顧橋礦南二采區(qū)1414(1)工作面概率積分法參數(shù)反演，將由反演結(jié)果計算的下沉值、水平移動擬合值與真實值對比，獲得的下沉值與水平移動擬合中誤差為114.81 mm，下沉值擬合中誤差約為最大下沉值的2.7%，水平移動值擬合中誤差約為最大水平移動值的8.6%，滿足工程應(yīng)用精度。