張小艷，許慧，姜水軍

(1.西安科技大學(xué) 計算機科學(xué)與技術(shù)學(xué)院，陜西 西安 710600; 2.中國神華神東煤炭集團，陜西 榆林 719315)

0 引言

近年來，隨著數(shù)字礦山概念的提出，三維地質(zhì)建模在行業(yè)領(lǐng)域中備受關(guān)注。因此，對綜采工作面煤炭賦存形態(tài)空間三維可視化成為煤炭生產(chǎn)管理過程中的一個重要環(huán)節(jié)。

目前三維地質(zhì)建模在行業(yè)領(lǐng)域中備受關(guān)注。范文瑤等[1]基于GOCAD軟件對研究區(qū)域進行三維建模與可視化操作；尚福華等[2]為了解決傳統(tǒng)三維地質(zhì)模型中不能有效地分析地層和斷層之間的形態(tài)和結(jié)構(gòu)，構(gòu)建出了一種TIN-Octree混合空間數(shù)據(jù)模型的精確三維地質(zhì)結(jié)構(gòu)方法；宋越等[3]將研究的重點放在了解決邊界不連續(xù)的問題層面上，提出了一種基于機器學(xué)習(xí)的擴展紋理樣本庫的方法。

空間數(shù)據(jù)插值是三維建模過程中的一個至關(guān)重要的環(huán)節(jié)。自20世紀(jì)初以來，就有大量的學(xué)者對此展開研究并取得了卓越的成效。Kriging算法被廣泛用于各類觀測的空間插值。本文提出了一種用自適應(yīng)PSO算法優(yōu)化Kriging插值算法(UPSO-Kriging)，將其用于工作面煤層高程值預(yù)測，為繪制工作面煤層賦存形態(tài)三維圖提供基礎(chǔ)數(shù)據(jù)，并以神東集團大柳塔煤礦604工作面的實際數(shù)據(jù)驗證了算法的可行性,進而基于規(guī)則網(wǎng)格法構(gòu)建DEM數(shù)字高程模型，以Three.js為平臺實現(xiàn)了瀏覽器端綜采工作面三維可視化，為煤炭企業(yè)開采提供科學(xué)依據(jù)。

1 Kriging插值法

Kriging插值法是一種基于協(xié)方差函數(shù)對空間隨機過程建模和預(yù)測的回歸算法[4]。一般公式可表示為：

(1)

式中:待測樣本點的值為Z*(x)，已知樣本點的值以及個數(shù)分別為Z(xi)、n，第i個已知樣本點相對于未知樣本點的權(quán)重為λi[5]。

當(dāng)Lagrange乘數(shù)法求取的最小估計方差滿足本征假設(shè)時，用變差函數(shù)表示Kriging方程組如式(2)所示：

(2)

式中：μ為Lagrange系數(shù)，γ(xi,xj)為兩個已知采樣點xi,xj的變差函數(shù)值。

2 UPSO-Kriging插值法

變差函數(shù)模型擬合中的參數(shù)確定是Kriging插值過程中的關(guān)鍵所在[6-7]。本文針對傳統(tǒng)Kriging插值法在變差函數(shù)模型擬合中存在的過擬合或欠擬合而導(dǎo)致算法整體精度下降的問題，提出一種自適應(yīng)PSO算法來優(yōu)化模型中的3個參數(shù)，提高算法整體預(yù)測精度。

2.1 自適應(yīng)PSO算法

2.1.1 PSO算法

粒子群優(yōu)化算法(PSO)是目前廣泛運用于各個領(lǐng)域的進化算法[8]，具有簡單、準(zhǔn)確、收斂速度快等特點。標(biāo)準(zhǔn)PSO算法更新公式如下[9]：

c2r2[pg(t)-xi(t)],

(3)

xi(t+1)=xi(t)+vi(t+1),

(4)

2.1.2 UPSO算法

針對PSO算法存在的搜索時間過長、易陷入局部極值的最優(yōu)化等問題[13]，本文采用張小艷等[14]提出的優(yōu)化方法。分別對算法尋優(yōu)過程中的粒子位置更新方式、迭代預(yù)選值處理以及基于正態(tài)分布的慣性調(diào)整策略進行優(yōu)化。

為驗證UPSO算法的有效性，本文分別選取4個函數(shù)作為標(biāo)準(zhǔn)測試函數(shù)，與PSO算法在python2.7環(huán)境下進行仿真測試。其中包括兩個單峰函數(shù)：Sphere函數(shù)與Schwefel函數(shù)；兩個多峰函數(shù)：Ackley函數(shù)與Rastrigrin函數(shù)。單峰函數(shù)用來測試算法的收斂速度以檢驗其局部搜索能力，多峰函數(shù)用來檢驗算法的全局搜索能力從而衡量算法綜合性能。

圖1 Sphere函數(shù)收斂測試Fig.1 Sphere function convergence test diagram

圖2 Schwefel函數(shù)收斂測試Fig.2 Schwefel function convergence test diagram

圖3 Ackley函數(shù)收斂測試Fig.3 Ackley function convergence test diagram

圖4 Rastrigrin函數(shù)收斂測試Fig.4 Rastrigrin function convergence test diagram

由上圖可知，不管是在單峰測試函數(shù)下的局部搜索能力還是多峰函數(shù)下的全局搜索能力，UPSO算法相比PSO算法來說,得到最優(yōu)解的收斂速度更快，陷入局部極值的次數(shù)明顯少于原始算法。這均表明UPSO算法整體性能得到了提升。

2.2 UPSO-Kriging插值法

UPSO算法優(yōu)化Kriging插值法，主要通過UPSO算法迭代尋優(yōu)得到最優(yōu)的模型參數(shù)進行變異函數(shù)模型的擬合，進而得到預(yù)測結(jié)果。算法步驟如下:

1)確定待求解問題的參數(shù)與目標(biāo)函數(shù)。以變差函數(shù)模型中的3個模型參數(shù)作為求解參數(shù)，即塊金值c0，偏基臺值c，變程a，目標(biāo)函數(shù)為不同距離下的實驗與理論變差函數(shù)差值平方和[15],即

(5)

2)初始化粒子群。設(shè)群體規(guī)模為M，將3個模型參數(shù)作為一個粒子個體，即xi=(c0,c,a),i=1,2,…,M。種群中每個個體的參數(shù)初始化取值如下：

(6)

其中，Rmax是實驗變差函數(shù)的最大值，Amax為最大滯后距[16]。

3)確定UPSO算法的有關(guān)運行參數(shù)。其中包括慣性權(quán)重最大最小取值vmax,vmin;學(xué)習(xí)因子最大最小取值c1s,c1f,c2s,c2f；迭代次數(shù)Iter；速度vmax,vmin；位置xmax,xmin。

4)UPSO優(yōu)化運算。優(yōu)化過程中，通過目標(biāo)函數(shù)來評價群體中每個個體的適合度值是否達(dá)到標(biāo)準(zhǔn)，通過更新迭代過程中的速度和位置進行迭代尋優(yōu)運算。若當(dāng)前迭代中的最優(yōu)解不滿足終止條件則繼續(xù);若滿足，則停止。最終結(jié)果即最優(yōu)變差函數(shù)參數(shù)，由此可擬合出相應(yīng)的變差函數(shù)表達(dá)式。

5)Kriging插值預(yù)測。將最終得到的變差函數(shù)代入式(2)即可得到權(quán)重系數(shù)，再通過式(1)計算出待預(yù)測點的屬性值。

3 實驗驗證

3.1 數(shù)據(jù)準(zhǔn)備

本文運用神東集團大柳塔煤礦604綜采工作面，寬度300 m，深度3 400 m，煤層厚度在6.24～40.17 m的69組煤層煤樣數(shù)據(jù)作為樣本數(shù)據(jù)集，驗證UPSO-Kriging算法應(yīng)用于綜采工作面煤層高程屬性預(yù)測的可行性。部分樣本數(shù)據(jù)如表1所示。

表1 采樣點數(shù)據(jù)

3.2 基于UPSO-Kriging的高程值預(yù)測

基于樣本數(shù)據(jù)集合中同一煤層下的數(shù)據(jù)，分別優(yōu)化UPSO-Kriging插值中的球面、指數(shù)以及高斯模型，并與PSO-Kriging以及Kriging的插值預(yù)測效果進行比較，其中PSO算法相關(guān)參數(shù)設(shè)置如表2所示。變差函數(shù)模型參數(shù)的邊界區(qū)間如表3所示。

表2 PSO算法的參數(shù)設(shè)置

表3 變差參數(shù)的邊界設(shè)置

將實際采樣點高程與預(yù)測點之間的均方根誤差(RMSE)和平均絕對百分誤差(MAE)作為算法評價指標(biāo)進行分析，其中MAE和RMSE公式分別為：

(7)

(8)

由表4可知，在UPSO-Kriging插值方法下將指數(shù)模型作為變差函數(shù)理論模型時，效果最佳。結(jié)合69組初始樣本點數(shù)據(jù)繪制高程值分布折線A；然后分別根據(jù)UPSO-Kriging算法、PSO-Kriging算法計算出待預(yù)測點的高程值，繪制出相應(yīng)高程預(yù)測折線B及折線C，最終實驗對比圖如圖5所示。

由表4和圖5可以看出， UPSO-Kriging算法在各種評價指標(biāo)下插值的結(jié)果都優(yōu)于其他算法。因此，基于UPSO-Kriging算法的工作面煤層高程值預(yù)測完全可行。

表4 估算精度比較

圖5 采樣點高程分布對比實驗Fig.5 Broken line diagram of elevation distribution of sampling points

4 工作面煤層三維模型

4.1 基于規(guī)則網(wǎng)格的DEM模型建立

數(shù)字高程模型是利用有限的地形數(shù)據(jù)對地形進行數(shù)值模擬，通常分為規(guī)則網(wǎng)格DEM模型、TIN模型、等高線模型等[17]。本文結(jié)合煤層實際地形特點，基于matlab仿真實現(xiàn)的單層DEM高程模型以及多層DEM高程模型分別如圖6,7所示。

圖6 單層DEM模型Fig.6 Single layer DEM model

4.2 基于Three.js的三維模型建立

Three.js 是一種基于javascript編寫的第三方WEBGL庫[18]。由于其速度快、易用性強、易于交互等特點，近年來被廣泛應(yīng)用于三維建?？梢暬I(lǐng)域。本文運用Three.js平臺實現(xiàn)瀏覽器端三維模型建立，并實現(xiàn)人機交互操作。

人機交互功能主要采用Three.js 中的THREE. OrbitControls控件實現(xiàn)，在這種實現(xiàn)模式下，用戶可以通過鼠標(biāo)對界面圖形進行放大、縮小、旋轉(zhuǎn)、平移等相關(guān)操作，從不同的角度觀察工作面煤層分布情況。如圖8所示，其中黃色為矸石，深色為煤層。

圖7 多層DEM模型Fig.7 Multi-layer DEM model

圖8 人機交互Fig.8 Human-computer interaction

5 結(jié)論

論文以神東集團大柳塔煤礦604工作面的實際數(shù)據(jù)，運用本文提出的UPSO-Kriging插值法對工作面分層結(jié)構(gòu)的高程進行預(yù)測，進而基于規(guī)則網(wǎng)格法構(gòu)建DEM數(shù)字高程模型，以Three.js為平臺實現(xiàn)了瀏覽器端綜采工作面三維可視化，并與傳統(tǒng)的Kriging算法和PSO-Kriging算法進行對比，結(jié)果表明，利用UPSO-Kriging方法求得的空間插值數(shù)據(jù)更能反映實際的地質(zhì)特性，能夠比較準(zhǔn)確地反映工作面地質(zhì)結(jié)構(gòu)。從實用意義上說,工作面煤層賦存狀態(tài)三維可視化有著十分重要的意義，可以為煤炭企業(yè)的透明開采、智能開采、分質(zhì)開采提供依據(jù)。