王艷彬,黃寶柱,郭 兵,鄭慶學(xué),劉義生 ,田秀國(guó),李培現(xiàn)
(1.北京龍軟科技股份有限公司,北京 100190;2.開灤(集團(tuán))有限責(zé)任公司,河北 唐山 063000;3.中國(guó)礦業(yè)大學(xué)(北京) 地球科學(xué)與測(cè)繪工程學(xué)院,北京 100083)
目前,我國(guó)“三下” 壓煤約為137.9億t[1-2],部分礦區(qū)“三下”壓煤量巨大,嚴(yán)重制約礦區(qū)的經(jīng)濟(jì)社會(huì)發(fā)展,可靠的沉陷預(yù)測(cè)預(yù)報(bào)是降低礦區(qū)沉陷災(zāi)害和土地?fù)p害的關(guān)鍵,對(duì)指導(dǎo)“三下”開采實(shí)踐有著重要的作用。波蘭學(xué)者李特威尼申于上世紀(jì)50年代將隨機(jī)介質(zhì)理論引入巖層移動(dòng)研究,在此基礎(chǔ)上,我國(guó)學(xué)者劉寶琛、廖國(guó)華等[3]發(fā)展成為概率積分法。目前概率積分法已經(jīng)成為我國(guó)應(yīng)用最為成熟、最為廣泛的預(yù)計(jì)方法。中國(guó)礦業(yè)大學(xué)、煤炭科學(xué)研究總院特采所等單位相繼開發(fā)了基于概率積分法的開采沉陷預(yù)計(jì)軟件。但目前大多數(shù)煤礦開采沉陷預(yù)計(jì)軟件需要和CAD或其他軟件平臺(tái)結(jié)合,CAD不具有GIS軟件的特點(diǎn),對(duì)圖形屬性處理能力較弱,其他軟件平臺(tái)又和煤礦的日常生產(chǎn)工作聯(lián)系不緊密,因此,在前人研究的基礎(chǔ)上,依托于目前在煤礦企業(yè)應(yīng)用廣泛的北京龍軟科技股份有限公司自主研發(fā)的LongRuan GIS平臺(tái)上進(jìn)行沉陷預(yù)計(jì)與分析[4-6],旨在進(jìn)行地表移動(dòng)變形預(yù)計(jì)與分析的自動(dòng)化處理及結(jié)果的可視化輸出,降低沉陷預(yù)測(cè)工作難度,提高開采沉陷預(yù)測(cè)效率。
概率積分法是一種以隨機(jī)介質(zhì)理論為基礎(chǔ)的開采沉陷預(yù)計(jì)方法,概率積分法因其所用的移動(dòng)和變形預(yù)計(jì)公式中含有概率積分(或真導(dǎo)數(shù))而得名。其發(fā)展已較為成熟,能較精確地描述地表沉陷的分布形態(tài)[7]。目前,已成為我國(guó)應(yīng)用最為廣泛的傾斜煤層地表移動(dòng)變形預(yù)計(jì)方法之一。
地表移動(dòng)盆地內(nèi)任意一點(diǎn)A(x,y)下沉值可以用式(1) 來計(jì)算:
(1)
式中,WA(x,y)為坐標(biāo)A(x,y)的下沉值;Wcm為充分采動(dòng)條件下地表最大下沉值;r為主要影響半徑;D為地下開采區(qū)域。
傾斜變形i,曲率變形K,水平移動(dòng)U,水平變形ε都可由下沉表達(dá)式變換求導(dǎo)而得。
當(dāng)煤層傾角大于55°時(shí),煤層屬于急傾斜煤層。對(duì)于急傾斜煤層開采的地表下沉和變形預(yù)計(jì),概率積分法存在明顯不足,此時(shí)多采用皮爾遜Ⅲ型公式法[8],該方法是急傾斜煤層分階段開采地表移動(dòng)變形預(yù)計(jì)的一種常用剖面函數(shù)法。
急傾斜煤層開采時(shí)走向剖面的移動(dòng)和變形可按照公式(1)計(jì)算。而傾斜剖面不出現(xiàn)充分采動(dòng)條件,在下沉盆地傾向主斷面上,以底板移動(dòng)邊界為坐標(biāo)原點(diǎn)、下山方向?yàn)檩S的坐標(biāo)系統(tǒng)上,傾向主斷面的下沉和水平移動(dòng)按式(2)和式(3)計(jì)算[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ù)此,可以進(jìn)行地表下沉盆地傾向主斷面上的全部移動(dòng)變形計(jì)算。
對(duì)任意形狀的工作面采用沿走向或傾向把其切割為多個(gè)近似的小矩形工作面,計(jì)算每個(gè)小矩形開采引起的同名移動(dòng)變形值進(jìn)行疊加,得到整個(gè)開采點(diǎn)引起的移動(dòng)變形值[9]。
開采沉陷預(yù)計(jì)分析系統(tǒng)基于LongRuan GIS平臺(tái)開發(fā)完成,僅使用LongRuan GIS平臺(tái)即可完成煤層工作面數(shù)據(jù)預(yù)處理、沉陷預(yù)計(jì)計(jì)算和預(yù)計(jì)結(jié)果展示分析等功能。
本沉陷預(yù)計(jì)分析系統(tǒng)分為五大模塊,即錄入工作面數(shù)據(jù)模塊、沉陷預(yù)計(jì)計(jì)算模塊、繪制沉陷預(yù)計(jì)等值線模塊、繪制沉陷變形剖面圖模塊、繪制損害等級(jí)云圖模塊。本系統(tǒng)中傾斜煤層開采沉陷預(yù)計(jì)采用概率積分法模型,急傾斜煤層開采傾向方向沉陷預(yù)計(jì)采用皮爾遜Ⅲ型公式法模型。
為了實(shí)現(xiàn)壓煤村莊搬遷與塌陷地治理補(bǔ)償業(yè)務(wù)數(shù)據(jù)的標(biāo)準(zhǔn)化、規(guī)范化和存儲(chǔ)共享化,建立了傾斜和急傾斜煤層工作面空間數(shù)據(jù)庫(kù)。圖1所示為傾斜煤層工作面數(shù)據(jù)管理界面,圖2所示為急傾斜煤層工作面數(shù)據(jù)管理界面。以傾斜煤層工作面數(shù)據(jù)錄入為例,本系統(tǒng)規(guī)定煤層工作面數(shù)據(jù)包括煤層參數(shù)、工作面參數(shù)和預(yù)計(jì)參數(shù),煤層參數(shù)包括煤層傾角、煤層編號(hào)、煤層厚度,工作面參數(shù)包括指向下山方位角和任意形狀工作面各個(gè)角點(diǎn)的(x,y)坐標(biāo)和每個(gè)角點(diǎn)的采深,預(yù)計(jì)參數(shù)包括下沉系數(shù)、影響傳播系數(shù)、水平移動(dòng)系數(shù)、主要影響角正切和拐點(diǎn)偏移距。其中工作面的角點(diǎn)坐標(biāo)可以從LongRuan GIS地圖上拾取,也可以從Ascii文件中導(dǎo)入。傾斜和急傾斜煤層工作面空間數(shù)據(jù)庫(kù)的建立為開采沉陷預(yù)測(cè)、煤層底板等高線修正提供有效的基礎(chǔ)數(shù)據(jù)支持。
圖1 傾斜煤層工作面數(shù)據(jù)管理界面
圖2 急傾斜煤層工作面數(shù)據(jù)管理界面
本系統(tǒng)中錄入工作面數(shù)據(jù)模塊用于對(duì)空間數(shù)據(jù)庫(kù)中傾斜、急傾斜煤層工作面數(shù)據(jù)統(tǒng)一管理,包括對(duì)工作面數(shù)據(jù)進(jìn)行增、刪、改、查操作,為了方便管理,在空間數(shù)據(jù)庫(kù)中設(shè)置了第一等級(jí)項(xiàng)目名稱和第二等級(jí)工作面名稱兩個(gè)等級(jí),項(xiàng)目名稱等級(jí)包含工作面名稱等級(jí),項(xiàng)目名稱等級(jí)中包含具有同一屬性的回采工作面數(shù)據(jù),方便對(duì)同一項(xiàng)目下回采工作面數(shù)據(jù)預(yù)計(jì)地表的移動(dòng)變形值,例如同一煤層的回采工作面數(shù)據(jù)或者同一年份的回采工作面數(shù)據(jù)。工作面數(shù)據(jù)預(yù)處理工作包括在礦井采掘工程平面圖中拾取各個(gè)回采工作面邊界坐標(biāo),通過圖上鉆孔煤層數(shù)據(jù)計(jì)算回采工作面的角點(diǎn)深度,獲取回采工作面對(duì)應(yīng)的煤層編號(hào)、煤層厚度和煤層傾角,以及煤礦測(cè)量人員提供的回采工作面對(duì)應(yīng)的預(yù)計(jì)參數(shù),以上信息組合成回采工作面數(shù)據(jù),錄入到空間數(shù)據(jù)庫(kù)中。
系統(tǒng)使用不規(guī)則三角網(wǎng)[10]模型存儲(chǔ)沉陷預(yù)計(jì)結(jié)果數(shù)據(jù),程序中自動(dòng)預(yù)計(jì)了傾向和走向兩個(gè)方向的地表移動(dòng)變形值,可以預(yù)計(jì)指定方向的地表移動(dòng)變形值。利用固定間距格網(wǎng)點(diǎn)上的移動(dòng)變形值信息生成DTin數(shù)據(jù)并存儲(chǔ)在指定的LongRuan GIS平臺(tái)地圖圖層中。DTin數(shù)據(jù)是沉陷預(yù)計(jì)可視化圖形繪制和預(yù)計(jì)分析的基礎(chǔ)數(shù)據(jù)。
DTin數(shù)據(jù)中包括節(jié)點(diǎn)數(shù)據(jù)DTinNode和三角面數(shù)據(jù)DTinFace,數(shù)據(jù)結(jié)構(gòu)分別如下:
struct DTinNode
{
int no;//>--節(jié)點(diǎn)編號(hào)
int marker;//>--節(jié)點(diǎn)類型,-1為默認(rèn)【無(wú)類型】
double pos[3];//>--節(jié)點(diǎn)坐標(biāo)X,Y,w, 【w為下沉值】
double value[12];//>--節(jié)點(diǎn)屬性value[0-2]走向、傾向、指定方向傾斜值;value[3-5]走向、傾向、指定方向曲率值;value[6-8]走向、傾向、指定方向水平移動(dòng)值;value[9-12]走向、傾向、指定方向的水平變形值。
}
struct DTinFace
{
int no;//>--三角面片編號(hào)
int node_no[3];//>--頂點(diǎn)編號(hào)
int nbour_no[3];//>--關(guān)聯(lián)三角面片編號(hào)
int marker;//>--標(biāo)識(shí):外部布爾運(yùn)算應(yīng)用,默認(rèn)為-1
}
DTinNode數(shù)據(jù)結(jié)構(gòu)[11]中保存了三個(gè)方向的傾斜值、曲率值、水平移動(dòng)值、水平變形值,共計(jì)有12個(gè)屬性值,生成一次三角網(wǎng)通過屬性的選取可以獲取最多三個(gè)方向的傾斜值、曲率值、水平移動(dòng)值、水平變形值,減少了數(shù)據(jù)存儲(chǔ)量,提高了計(jì)算速度。
開采沉陷預(yù)計(jì)分析系統(tǒng)操作流程為:首先通過錄入工作面數(shù)據(jù)模塊在空間數(shù)據(jù)庫(kù)中錄入傾斜、急傾斜煤層工作面數(shù)據(jù),通過沉陷預(yù)計(jì)計(jì)算模塊選擇工作面名稱進(jìn)行沉陷預(yù)計(jì)計(jì)算,并把預(yù)計(jì)結(jié)果保存到Longruan GIS指定圖層下,再通過獲取Longruan GIS指定圖層下的預(yù)計(jì)結(jié)果繪制沉陷預(yù)計(jì)等值線、沉陷變形剖面圖、損害等級(jí)云圖、統(tǒng)計(jì)沉陷預(yù)計(jì)信息。系統(tǒng)通過可視化的交互界面,降低了開采地表變形分析基礎(chǔ)數(shù)據(jù)獲取的難度,提高了開采沉陷預(yù)測(cè)的工作效率。
建立傾斜和急傾斜煤層工作面空間數(shù)據(jù)庫(kù)后,利用沉陷預(yù)計(jì)計(jì)算模塊,如圖3所示,預(yù)計(jì)地表移動(dòng)變形。沉陷預(yù)計(jì)是整個(gè)系統(tǒng)實(shí)現(xiàn)中最重要的部分,通過傾斜或急傾斜按鈕讀取空間數(shù)據(jù)庫(kù)中的工作面數(shù)據(jù),用戶可選擇工程名稱以及工程名稱對(duì)應(yīng)的工作面名稱進(jìn)行組合沉陷預(yù)計(jì),其中預(yù)計(jì)范圍可通過程序自適應(yīng)生成,也可手動(dòng)拾取坐標(biāo)的方式設(shè)置。程序中自動(dòng)計(jì)算了傾向和走向兩個(gè)方向的變形,還可以選擇計(jì)算指定方向的變形。利用預(yù)計(jì)結(jié)果信息生成不規(guī)則三角網(wǎng)模型數(shù)據(jù)(以下簡(jiǎn)稱DTin數(shù)據(jù))并存儲(chǔ)在LongRuan GIS平臺(tái)指定圖層下。
圖3 沉陷預(yù)計(jì)計(jì)算交互界面
系統(tǒng)利用C++編程語(yǔ)言編寫了基于概率積分法的預(yù)計(jì)模型算法和基于皮爾遜Ⅲ型公式法的預(yù)計(jì)模型算法,通過傾斜和急傾斜煤層工作面數(shù)據(jù)標(biāo)識(shí)自動(dòng)識(shí)別對(duì)應(yīng)的模型算法,計(jì)算地表移動(dòng)變形值。
繪制沉陷預(yù)計(jì)等值線模塊如圖4所示,用于獲取LongRuan GIS平臺(tái)指定圖層下的預(yù)計(jì)結(jié)果DTin數(shù)據(jù),繪制下沉、傾斜、曲率、水平移動(dòng)、水平變形等值線專題圖及注記。通過LongRuan GIS平臺(tái)中自帶的“修改”菜單可以對(duì)繪制的等值線實(shí)體進(jìn)行編輯、裁剪、延長(zhǎng)、打斷等操作。
圖4 生成沉陷預(yù)計(jì)等值線交互界面
本系統(tǒng)根據(jù)間隔距離生成格網(wǎng)點(diǎn)矩陣,每個(gè)格網(wǎng)點(diǎn)(x,y)都記錄了在DTinNode數(shù)據(jù)中的點(diǎn)號(hào),通過選擇移動(dòng)變形參數(shù)、方向類型,即可獲取對(duì)應(yīng)的移動(dòng)變形值做為格網(wǎng)點(diǎn)(x,y)的Z值,根據(jù)得到的(x,y,z)坐標(biāo)值,便可通過編寫的繪制等值線函數(shù)生成對(duì)應(yīng)的等值線[12]專題圖,如圖5所示,代碼為:
圖5 急傾斜煤層開采下沉等值線專題圖
GetPntDataByDTinNode(vSurPnt, srcDtin) ; //>-- 通過DTin數(shù)據(jù)獲取格網(wǎng)點(diǎn)坐標(biāo)(x,y)和對(duì)應(yīng)的移動(dòng)變形參數(shù)
GetIsolinesWithGridMethod_NoInter(vvIsolines, vSurPnt, iRowNum, iColNum, strIsolines, 0, 0, 0, m_dbConstPrecision, CallBack_Param)) ; //>-- 利用格網(wǎng)法追蹤等值線-不插值提高速度
LongRuan GIS平臺(tái)可以查詢繪制的等值線專題圖中指定點(diǎn)的坐標(biāo)、等值線的二三維長(zhǎng)度、封閉等值線區(qū)域的面積,通過插值下沉等值線邊界線面積即可以計(jì)算地表沉陷的影響范圍。
系統(tǒng)繪制沉陷變形剖面圖模塊用于獲取指定圖層的DTin數(shù)據(jù)繪制走向和傾向的主斷面圖,如圖6所示,以及選擇剖面線繪制指定位置的斷面圖,方便查看建筑物、鐵路、公路、河流等對(duì)應(yīng)位置的下沉情況[13]。繪制損害等級(jí)云圖模塊可按礦區(qū)建筑物損害等級(jí)劃分標(biāo)準(zhǔn),在LongRuan GIS地圖中利用傾斜、曲率、水平變形DTin數(shù)據(jù)生成臨界值等值線,再通過編寫的最小閉合圈算法繪制損害等級(jí)云圖,由此可判斷建筑物是否受到開采影響以及所受采動(dòng)影響的程度[14]。
圖6 走向和傾向主斷面圖
開灤集團(tuán)錢家營(yíng)煤礦開采7煤,地面標(biāo)高為20~23 m,采用綜放開采,煤層采深為310 ~420 m,工作面位于礦井的東北部,煤層走向方位角80°~95°,煤層傾角α=7°~15°,煤層厚度4.5~5.6 m。根據(jù)鄰近采區(qū)實(shí)測(cè)參數(shù)及開采工作面地質(zhì)條件,選用預(yù)計(jì)參數(shù)為:下沉系數(shù)0.85,水平移動(dòng)系數(shù)0.3,主要影響角正切2.5,開采影響角系數(shù)k=0.6,開采影響角90-k×α。采用本系統(tǒng)對(duì)開灤集團(tuán)錢家營(yíng)礦7煤多個(gè)傾斜煤層工作面進(jìn)行沉陷預(yù)計(jì)計(jì)算,預(yù)計(jì)結(jié)果如圖7所示,預(yù)計(jì)最大下沉值4.37 m,沉陷區(qū)域影響面積為6 588 870 m2,沉陷區(qū)域體積為10 590 572 m3,繪制了地表下沉等值線,確定了地面建筑物損害等級(jí),并計(jì)算了各個(gè)等級(jí)區(qū)域的面積,與現(xiàn)場(chǎng)調(diào)查結(jié)果比對(duì)后認(rèn)為,預(yù)計(jì)結(jié)果可以滿足礦區(qū)生產(chǎn)規(guī)劃要求。
圖7 多工作面的下沉等值線和損害等級(jí)云圖
(1)基于LongRuan GIS平臺(tái)的數(shù)據(jù)處理和空間數(shù)據(jù)分析能力,開展沉陷預(yù)計(jì)分析系統(tǒng)可視化研究。通過動(dòng)態(tài)庫(kù)的方式進(jìn)行加載,完成煤層工作面數(shù)據(jù)預(yù)處理、沉陷預(yù)計(jì)計(jì)算和預(yù)計(jì)結(jié)果展示等功能。
(2)采用可視化的參數(shù)界面對(duì)空間數(shù)據(jù)庫(kù)中的煤層工作面數(shù)據(jù)統(tǒng)一管理,建立煤層工作面預(yù)處理、煤層工作面數(shù)據(jù)入庫(kù)、沉陷預(yù)計(jì)計(jì)算、預(yù)計(jì)結(jié)果展示和分析等功能。
(3)結(jié)合工程實(shí)例表明,該系統(tǒng)用戶界面友好, 操作簡(jiǎn)單方便、計(jì)算結(jié)果可靠,可為礦區(qū)開采沉陷治理、沉陷區(qū)居民搬遷、土地?fù)p毀修復(fù)提供有力支持。