王艷彬，黃寶柱，郭 兵，鄭慶學，劉義生 ，田秀國，李培現(xiàn)

(1.北京龍軟科技股份有限公司，北京 100190;2.開灤(集團)有限責任公司，河北 唐山 063000;3.中國礦業(yè)大學(北京) 地球科學與測繪工程學院，北京 100083)

目前，我國“三下” 壓煤約為137.9億t[1-2]，部分礦區(qū)“三下”壓煤量巨大，嚴重制約礦區(qū)的經濟社會發(fā)展，可靠的沉陷預測預報是降低礦區(qū)沉陷災害和土地損害的關鍵，對指導“三下”開采實踐有著重要的作用。波蘭學者李特威尼申于上世紀50年代將隨機介質理論引入巖層移動研究，在此基礎上，我國學者劉寶琛、廖國華等[3]發(fā)展成為概率積分法。目前概率積分法已經成為我國應用最為成熟、最為廣泛的預計方法。中國礦業(yè)大學、煤炭科學研究總院特采所等單位相繼開發(fā)了基于概率積分法的開采沉陷預計軟件。但目前大多數(shù)煤礦開采沉陷預計軟件需要和CAD或其他軟件平臺結合，CAD不具有GIS軟件的特點，對圖形屬性處理能力較弱，其他軟件平臺又和煤礦的日常生產工作聯(lián)系不緊密，因此，在前人研究的基礎上，依托于目前在煤礦企業(yè)應用廣泛的北京龍軟科技股份有限公司自主研發(fā)的LongRuan GIS平臺上進行沉陷預計與分析[4-6]，旨在進行地表移動變形預計與分析的自動化處理及結果的可視化輸出，降低沉陷預測工作難度，提高開采沉陷預測效率。

1 預計模型

概率積分法是一種以隨機介質理論為基礎的開采沉陷預計方法，概率積分法因其所用的移動和變形預計公式中含有概率積分(或真導數(shù))而得名。其發(fā)展已較為成熟，能較精確地描述地表沉陷的分布形態(tài)[7]。目前,已成為我國應用最為廣泛的傾斜煤層地表移動變形預計方法之一。

地表移動盆地內任意一點A(x,y)下沉值可以用式(1) 來計算：

(1)

式中,WA(x,y)為坐標A(x,y)的下沉值；Wcm為充分采動條件下地表最大下沉值；r為主要影響半徑；D為地下開采區(qū)域。

傾斜變形i，曲率變形K，水平移動U，水平變形ε都可由下沉表達式變換求導而得。

當煤層傾角大于55°時，煤層屬于急傾斜煤層。對于急傾斜煤層開采的地表下沉和變形預計，概率積分法存在明顯不足，此時多采用皮爾遜Ⅲ型公式法[8]，該方法是急傾斜煤層分階段開采地表移動變形預計的一種常用剖面函數(shù)法。

急傾斜煤層開采時走向剖面的移動和變形可按照公式(1)計算。而傾斜剖面不出現(xiàn)充分采動條件，在下沉盆地傾向主斷面上，以底板移動邊界為坐標原點、下山方向為軸的坐標系統(tǒng)上，傾向主斷面的下沉和水平移動按式(2)和式(3)計算[3]：

W(y)=a1WmaxZa2exp (-a3Z)

(2)

u(y)=q[B-P(y)]i(y)

(3)

其中,

變形公式為(4)～(6)：

(4)

(5)

ε(y)=qBb1(1-Z)b2exp[-b3(1-Z)]

(6)

據(jù)此，可以進行地表下沉盆地傾向主斷面上的全部移動變形計算。

對任意形狀的工作面采用沿走向或傾向把其切割為多個近似的小矩形工作面，計算每個小矩形開采引起的同名移動變形值進行疊加,得到整個開采點引起的移動變形值[9]。

2 軟件系統(tǒng)設計和實現(xiàn)

開采沉陷預計分析系統(tǒng)基于LongRuan GIS平臺開發(fā)完成，僅使用LongRuan GIS平臺即可完成煤層工作面數(shù)據(jù)預處理、沉陷預計計算和預計結果展示分析等功能。

本沉陷預計分析系統(tǒng)分為五大模塊，即錄入工作面數(shù)據(jù)模塊、沉陷預計計算模塊、繪制沉陷預計等值線模塊、繪制沉陷變形剖面圖模塊、繪制損害等級云圖模塊。本系統(tǒng)中傾斜煤層開采沉陷預計采用概率積分法模型，急傾斜煤層開采傾向方向沉陷預計采用皮爾遜Ⅲ型公式法模型。

2.1 傾斜和急傾斜煤層工作面空間數(shù)據(jù)庫的建立

為了實現(xiàn)壓煤村莊搬遷與塌陷地治理補償業(yè)務數(shù)據(jù)的標準化、規(guī)范化和存儲共享化，建立了傾斜和急傾斜煤層工作面空間數(shù)據(jù)庫。圖1所示為傾斜煤層工作面數(shù)據(jù)管理界面，圖2所示為急傾斜煤層工作面數(shù)據(jù)管理界面。以傾斜煤層工作面數(shù)據(jù)錄入為例，本系統(tǒng)規(guī)定煤層工作面數(shù)據(jù)包括煤層參數(shù)、工作面參數(shù)和預計參數(shù)，煤層參數(shù)包括煤層傾角、煤層編號、煤層厚度，工作面參數(shù)包括指向下山方位角和任意形狀工作面各個角點的(x,y)坐標和每個角點的采深，預計參數(shù)包括下沉系數(shù)、影響傳播系數(shù)、水平移動系數(shù)、主要影響角正切和拐點偏移距。其中工作面的角點坐標可以從LongRuan GIS地圖上拾取，也可以從Ascii文件中導入。傾斜和急傾斜煤層工作面空間數(shù)據(jù)庫的建立為開采沉陷預測、煤層底板等高線修正提供有效的基礎數(shù)據(jù)支持。

圖1 傾斜煤層工作面數(shù)據(jù)管理界面

圖2 急傾斜煤層工作面數(shù)據(jù)管理界面

本系統(tǒng)中錄入工作面數(shù)據(jù)模塊用于對空間數(shù)據(jù)庫中傾斜、急傾斜煤層工作面數(shù)據(jù)統(tǒng)一管理，包括對工作面數(shù)據(jù)進行增、刪、改、查操作，為了方便管理，在空間數(shù)據(jù)庫中設置了第一等級項目名稱和第二等級工作面名稱兩個等級，項目名稱等級包含工作面名稱等級，項目名稱等級中包含具有同一屬性的回采工作面數(shù)據(jù)，方便對同一項目下回采工作面數(shù)據(jù)預計地表的移動變形值，例如同一煤層的回采工作面數(shù)據(jù)或者同一年份的回采工作面數(shù)據(jù)。工作面數(shù)據(jù)預處理工作包括在礦井采掘工程平面圖中拾取各個回采工作面邊界坐標，通過圖上鉆孔煤層數(shù)據(jù)計算回采工作面的角點深度，獲取回采工作面對應的煤層編號、煤層厚度和煤層傾角，以及煤礦測量人員提供的回采工作面對應的預計參數(shù)，以上信息組合成回采工作面數(shù)據(jù)，錄入到空間數(shù)據(jù)庫中。

2.2 不規(guī)則三角網模型存儲

系統(tǒng)使用不規(guī)則三角網[10]模型存儲沉陷預計結果數(shù)據(jù)，程序中自動預計了傾向和走向兩個方向的地表移動變形值，可以預計指定方向的地表移動變形值。利用固定間距格網點上的移動變形值信息生成DTin數(shù)據(jù)并存儲在指定的LongRuan GIS平臺地圖圖層中。DTin數(shù)據(jù)是沉陷預計可視化圖形繪制和預計分析的基礎數(shù)據(jù)。

DTin數(shù)據(jù)中包括節(jié)點數(shù)據(jù)DTinNode和三角面數(shù)據(jù)DTinFace，數(shù)據(jù)結構分別如下：

struct DTinNode

{

int no;//>--節(jié)點編號

int marker;//>--節(jié)點類型，-1為默認【無類型】

double pos[3];//>--節(jié)點坐標X,Y,w, 【w為下沉值】

double value[12];//>--節(jié)點屬性value[0-2]走向、傾向、指定方向傾斜值；value[3-5]走向、傾向、指定方向曲率值；value[6-8]走向、傾向、指定方向水平移動值；value[9-12]走向、傾向、指定方向的水平變形值。

}

struct DTinFace

{

int no;//>--三角面片編號

int node_no[3];//>--頂點編號

int nbour_no[3];//>--關聯(lián)三角面片編號

int marker;//>--標識：外部布爾運算應用，默認為-1

}

DTinNode數(shù)據(jù)結構[11]中保存了三個方向的傾斜值、曲率值、水平移動值、水平變形值，共計有12個屬性值，生成一次三角網通過屬性的選取可以獲取最多三個方向的傾斜值、曲率值、水平移動值、水平變形值，減少了數(shù)據(jù)存儲量，提高了計算速度。

3 軟件可視化交互

開采沉陷預計分析系統(tǒng)操作流程為：首先通過錄入工作面數(shù)據(jù)模塊在空間數(shù)據(jù)庫中錄入傾斜、急傾斜煤層工作面數(shù)據(jù)，通過沉陷預計計算模塊選擇工作面名稱進行沉陷預計計算，并把預計結果保存到Longruan GIS指定圖層下，再通過獲取Longruan GIS指定圖層下的預計結果繪制沉陷預計等值線、沉陷變形剖面圖、損害等級云圖、統(tǒng)計沉陷預計信息。系統(tǒng)通過可視化的交互界面，降低了開采地表變形分析基礎數(shù)據(jù)獲取的難度，提高了開采沉陷預測的工作效率。

3.1 沉陷預計計算

建立傾斜和急傾斜煤層工作面空間數(shù)據(jù)庫后，利用沉陷預計計算模塊，如圖3所示，預計地表移動變形。沉陷預計是整個系統(tǒng)實現(xiàn)中最重要的部分，通過傾斜或急傾斜按鈕讀取空間數(shù)據(jù)庫中的工作面數(shù)據(jù)，用戶可選擇工程名稱以及工程名稱對應的工作面名稱進行組合沉陷預計，其中預計范圍可通過程序自適應生成，也可手動拾取坐標的方式設置。程序中自動計算了傾向和走向兩個方向的變形，還可以選擇計算指定方向的變形。利用預計結果信息生成不規(guī)則三角網模型數(shù)據(jù)(以下簡稱DTin數(shù)據(jù))并存儲在LongRuan GIS平臺指定圖層下。

圖3 沉陷預計計算交互界面

系統(tǒng)利用C++編程語言編寫了基于概率積分法的預計模型算法和基于皮爾遜Ⅲ型公式法的預計模型算法，通過傾斜和急傾斜煤層工作面數(shù)據(jù)標識自動識別對應的模型算法，計算地表移動變形值。

3.2 沉陷預計等值線繪制

繪制沉陷預計等值線模塊如圖4所示，用于獲取LongRuan GIS平臺指定圖層下的預計結果DTin數(shù)據(jù)，繪制下沉、傾斜、曲率、水平移動、水平變形等值線專題圖及注記。通過LongRuan GIS平臺中自帶的“修改”菜單可以對繪制的等值線實體進行編輯、裁剪、延長、打斷等操作。

圖4 生成沉陷預計等值線交互界面

本系統(tǒng)根據(jù)間隔距離生成格網點矩陣，每個格網點(x,y)都記錄了在DTinNode數(shù)據(jù)中的點號，通過選擇移動變形參數(shù)、方向類型，即可獲取對應的移動變形值做為格網點(x,y)的Z值，根據(jù)得到的(x,y,z)坐標值，便可通過編寫的繪制等值線函數(shù)生成對應的等值線[12]專題圖，如圖5所示，代碼為：

圖5 急傾斜煤層開采下沉等值線專題圖

GetPntDataByDTinNode(vSurPnt, srcDtin) ; //>-- 通過DTin數(shù)據(jù)獲取格網點坐標(x,y)和對應的移動變形參數(shù)

GetIsolinesWithGridMethod_NoInter(vvIsolines, vSurPnt, iRowNum, iColNum, strIsolines, 0, 0, 0, m_dbConstPrecision, CallBack_Param)) ; //>-- 利用格網法追蹤等值線-不插值提高速度

3.3 沉陷預計結果分析

LongRuan GIS平臺可以查詢繪制的等值線專題圖中指定點的坐標、等值線的二三維長度、封閉等值線區(qū)域的面積，通過插值下沉等值線邊界線面積即可以計算地表沉陷的影響范圍。

系統(tǒng)繪制沉陷變形剖面圖模塊用于獲取指定圖層的DTin數(shù)據(jù)繪制走向和傾向的主斷面圖，如圖6所示，以及選擇剖面線繪制指定位置的斷面圖，方便查看建筑物、鐵路、公路、河流等對應位置的下沉情況[13]。繪制損害等級云圖模塊可按礦區(qū)建筑物損害等級劃分標準，在LongRuan GIS地圖中利用傾斜、曲率、水平變形DTin數(shù)據(jù)生成臨界值等值線，再通過編寫的最小閉合圈算法繪制損害等級云圖，由此可判斷建筑物是否受到開采影響以及所受采動影響的程度[14]。

圖6 走向和傾向主斷面圖

4 工程應用

開灤集團錢家營煤礦開采7煤，地面標高為20～23 m，采用綜放開采，煤層采深為310 ～420 m，工作面位于礦井的東北部，煤層走向方位角80°～95°，煤層傾角α=7°～15°，煤層厚度4.5～5.6 m。根據(jù)鄰近采區(qū)實測參數(shù)及開采工作面地質條件，選用預計參數(shù)為：下沉系數(shù)0.85，水平移動系數(shù)0.3，主要影響角正切2.5，開采影響角系數(shù)k=0.6，開采影響角90-k×α。采用本系統(tǒng)對開灤集團錢家營礦7煤多個傾斜煤層工作面進行沉陷預計計算，預計結果如圖7所示，預計最大下沉值4.37 m，沉陷區(qū)域影響面積為6 588 870 m2，沉陷區(qū)域體積為10 590 572 m3，繪制了地表下沉等值線，確定了地面建筑物損害等級，并計算了各個等級區(qū)域的面積，與現(xiàn)場調查結果比對后認為，預計結果可以滿足礦區(qū)生產規(guī)劃要求。

圖7 多工作面的下沉等值線和損害等級云圖

5 結 論

(1)基于LongRuan GIS平臺的數(shù)據(jù)處理和空間數(shù)據(jù)分析能力，開展沉陷預計分析系統(tǒng)可視化研究。通過動態(tài)庫的方式進行加載，完成煤層工作面數(shù)據(jù)預處理、沉陷預計計算和預計結果展示等功能。

(2)采用可視化的參數(shù)界面對空間數(shù)據(jù)庫中的煤層工作面數(shù)據(jù)統(tǒng)一管理，建立煤層工作面預處理、煤層工作面數(shù)據(jù)入庫、沉陷預計計算、預計結果展示和分析等功能。

(3)結合工程實例表明，該系統(tǒng)用戶界面友好, 操作簡單方便、計算結果可靠，可為礦區(qū)開采沉陷治理、沉陷區(qū)居民搬遷、土地損毀修復提供有力支持。