夏嘉南，李根生，卞正富，雷少剛，金慧霞,3

(1.中國礦業(yè)大學(xué) 公共管理學(xué)院，江蘇 徐州 221116；2.中國礦業(yè)大學(xué) 礦山生態(tài)修復(fù)教育部工程研究中心，江蘇 徐州 221116；3.溫州職業(yè)技術(shù)學(xué)院，浙江 溫州 325035)

露天開采作為我國草原礦區(qū)的重要開采方式，其內(nèi)排土場占地面積隨開采工作面的不斷推進(jìn)而增大。在此過程中，傳統(tǒng)地貌重塑方法往往將原始地貌轉(zhuǎn)變?yōu)椤芭_階式”或“錐狀堆積式”的人工堆疊地貌，使原有的地貌及水文特征消失，導(dǎo)致礦區(qū)重塑地貌與周邊自然地貌割裂。雖然在后續(xù)復(fù)墾設(shè)計中，內(nèi)排土場通過輸水管網(wǎng)布設(shè)以滿足復(fù)墾區(qū)重建植被的生長需求，但未能實現(xiàn)礦區(qū)內(nèi)外地貌水系的“無痕”拼接，從而在高強(qiáng)度降雨條件下，易引發(fā)邊界處嚴(yán)重表土水蝕等問題。相對人工規(guī)則地貌，自然穩(wěn)定地貌作為長期演變的結(jié)果，具有較好的適宜性與系統(tǒng)性。因此，可嘗試在采排復(fù)過程中保留自然地貌及其水文特征以破解上述困境。

當(dāng)前，GeoFluv模型作為主要的近自然地貌重塑方法，被廣泛應(yīng)用于采損區(qū)地貌重塑，且塑形效果在美國墨西哥州和懷俄明州、西班牙、澳大利亞等礦區(qū)均得到了良好的檢驗，并在國內(nèi)得到了初步的運(yùn)用。然而，受模型框架影響，重塑區(qū)水文溝道位置需人為主觀設(shè)定，雖后續(xù)通過自然穩(wěn)定溝道特征值加以修正，但整體上其地貌求解結(jié)果主觀性較高。NbS(Nature-based Solutions)理論謀求人與生態(tài)系統(tǒng)互惠，與“保護(hù)優(yōu)先，自然恢復(fù)為主”的中國生態(tài)文明建設(shè)理念相契合。在上述研究及理論的指導(dǎo)下，團(tuán)隊內(nèi)諸多學(xué)者均提出師法自然穩(wěn)定地貌的內(nèi)排土場地貌重塑方法。其中，李恒等以內(nèi)蒙古勝利礦為例，通過將周邊自然穩(wěn)定地貌剖面形態(tài)應(yīng)用至排土場邊坡，通過GeoWep模擬論證仿自然形態(tài)邊坡較傳統(tǒng)“臺階式”邊坡更具抗水蝕能力。陳航等以內(nèi)蒙古勝利礦為例，通過提取周邊自然穩(wěn)定地貌特征參數(shù)，在不考慮采復(fù)周期的情景下采用GeoFlue軟件構(gòu)建研究區(qū)內(nèi)排土場近自然重塑地貌，并基于CLiDE模擬驗證其表土抗水蝕能力優(yōu)于傳統(tǒng)內(nèi)排土場地貌重塑結(jié)果。然而，在上述研究中，內(nèi)排土場采排復(fù)一體下復(fù)填區(qū)的可用土方量會隨采復(fù)周期變化。同時，相較于周邊自然穩(wěn)定地貌統(tǒng)計結(jié)果，區(qū)域原始自然穩(wěn)定地貌的特殊性(受區(qū)域氣候及相鄰區(qū)的長期作用影響)也未得到充分考慮。在前置研究中，筆者以新疆黑山露天礦為例，提出基于SAS(Surface Adjustment Surface)的內(nèi)排土場地貌近自然重塑模型，強(qiáng)調(diào)通過局部調(diào)整的方式以盡可能保留原始地貌起伏特征，同時避免模型構(gòu)建中水文溝道的人為主觀設(shè)計，且CLiDE水蝕模擬表明，NNDL_SAS(Near-natural Design Landform with Surface Adjustment Surface)相較于TDL(Traditional Design Landform)在10 a當(dāng)?shù)貧夂驐l件下可減少約39.07%的表土水蝕量。然而，后續(xù)論證表明，上述模型求解形態(tài)局部調(diào)整可能會造成水系斷流、淤積等不可控問題。因此，需要一種既貼合周邊自然地貌，又保留原始地貌及水文特征的內(nèi)排土場地貌重塑方法，以削弱前置研究中地貌設(shè)計結(jié)果水系斷流、淤積等問題。

為此，筆者在前置研究的基礎(chǔ)上，以新疆黑山露天礦為例，提出一種基于HPCS(Hydrological Preserved Curved Surface)的內(nèi)排土場近自然地貌重塑模型，求解并以NNDL_HPCS(Near-natural Design Landform with Hydrological Preserved Curved Surface)為實驗組，NNDL_SAS,TDL和PNL(Pristine Natural Landform)為對照組，結(jié)合GIS空間分析軟件與CLiDE演化模型，模擬并討論NNDL_HPCS在地貌融合、土方調(diào)配和土壤抗水蝕3方面的改良效果。為相同及相似地區(qū)NbS理論下礦山近自然生態(tài)修復(fù)提供方法支持。

1 研究區(qū)概況與數(shù)據(jù)來源

1.1 研究區(qū)概況

研究區(qū)位于新疆黑山露天礦礦權(quán)范圍內(nèi)，地理坐標(biāo)：87.43°E—87.49°E，43.23°N—43.23°N，周邊地形以山間谷地為主，屬大陸半干旱及高寒氣候；區(qū)內(nèi)表土以第四系為主，野生植被類型為典型草原植被。采區(qū)東西長約4 874 m、南北長約1 000 m；采坑預(yù)設(shè)邊界為812 m×1 000 m的矩形，內(nèi)排作業(yè)由東至西(圖1)。

圖1 研究區(qū)采前原始DEM和真彩色衛(wèi)星影像Fig.1 DEM of the natural surface of the study area before acquisition and RGB satellite image of the study area

依據(jù)研究區(qū)煤炭資源儲量核實報告，區(qū)域可采煤層為6,7，8,9，11,12-1,12-2,13-1和13-2號。其中13-2號煤層為全區(qū)可露采、其余為局部可露采(圖2)。內(nèi)排土場臺階高度均設(shè)置為10 m。其中西幫為開采幫，幫坡角為20°；南北兩側(cè)邊幫幫坡角均為13°；東幫為內(nèi)排幫，幫坡角為8°。

1.2 數(shù)據(jù)來源

研究區(qū)采前原始DEM(Digital Elevation Model)、煤層頂板DEM和煤層底板DEM均來源于前置研究，空間分辨率均為30 m×30 m；內(nèi)排土場開采和復(fù)填子區(qū)空間分布數(shù)據(jù)與前置研究相同，均基于內(nèi)排土場采復(fù)子區(qū)位置識別方法獲??；復(fù)填子區(qū)內(nèi)的可用土方量數(shù)據(jù)通過可用土方量計算式求解獲?。粡?fù)填區(qū)HPCS數(shù)據(jù)通過HPCS預(yù)構(gòu)建模型獲??；復(fù)填子區(qū)NNDL_HPCS通過Matlab軟件基于土方平衡下區(qū)域坡度最緩原則定向篩選不同控制點系數(shù)下預(yù)構(gòu)建的HPCS獲得；復(fù)填子區(qū)NNDL_SAS來源于前置研究成果中“近自然設(shè)計DEM”；研究區(qū)PNL數(shù)據(jù)來源于前置研究成果中“自然原始地貌”；研究區(qū)內(nèi)排土TDL與前置研究成果中的“傳統(tǒng)設(shè)計DEM”一致(30°斜坡加平臺式“階梯狀”布局)；地貌融合效果中，各DEM相對于PNL的高程方差值通過ArcGIS空間疊加分析方法，逐柵格累計統(tǒng)計DEM與PNL間柵格高程差的平方值獲??；地貌融合效果中，自然水文保留率指標(biāo)通過ArcGIS坡向模塊和空間疊加分析方法，提取并對比研究區(qū)NNDL_HPCS,NNDL_SAS和TDL三者圖像柵格坡向與PNL空間分布一致性獲??；土方調(diào)配中，NNDL_HPCS,NNDL_SAS和TDL的土方平均運(yùn)距均由Matlab基于每立方土體空間直線位移距離統(tǒng)計獲??；最后運(yùn)用CLiDE演化模型和Matlab軟件，分別模擬并統(tǒng)計研究區(qū)NNDL_HPCS,NNDL_SAS,TDL和PNL四者10 a的土壤水蝕量，以評估NNDL_HPCS的抗水蝕能力。

圖2 研究區(qū)各煤層厚度及最底層煤層底板標(biāo)高空間分布Fig.2 Spatial distribution of the coal seam thickness and bottom plate elevation of the lowest coal seam in the study area

2 研究方法

2.1 基于HPCS的內(nèi)排土場近自然地貌重塑模型

2.1.1 采復(fù)子區(qū)空間位置識別

依據(jù)露天礦采復(fù)周期，獲取開采子區(qū)m和復(fù)填子區(qū)f的水平投影空間位置(圖3)。

圖3 礦區(qū)采復(fù)周期n中采前和復(fù)填后豎直剖面示意((a)，(b))以及采復(fù)子區(qū)網(wǎng)格單元示意(c)Fig.3 Schematic diagram of the vertical section before and after mining in the mining area under the open-pit mining-refill cycle n， schematic diagram of the grid unit of the mining-refill sub-area

如圖3所示，依據(jù)露天礦實際采復(fù)周期及其時間先后順序，沿開采方向?qū)⒙短斓V區(qū)劃分為個子區(qū)。其中，開采子區(qū)和復(fù)填子區(qū)分別對應(yīng)區(qū)域～m和區(qū)域～f。在露天礦內(nèi)排構(gòu)建的采復(fù)過程中，采坑每前進(jìn)一個開采子區(qū)m對應(yīng)修復(fù)后側(cè)一個復(fù)填子區(qū)f(圖3(a)，(b))，直至結(jié)束內(nèi)排過程，礦坑由區(qū)域不斷迭代前進(jìn)至采礦終止邊界。

2.1.2 可用土方量計算

基于輸入數(shù)據(jù)柵格空間分辨率，將開采及復(fù)填子區(qū)均細(xì)分為row行col列的網(wǎng)格單元組(圖3(c))。其中，網(wǎng)格單元大小與柵格分辨率一致。在采復(fù)周期中，開采子區(qū)m可用于復(fù)填子區(qū)f的土方量計算式為

(1)

式中，(m)為開采子區(qū)m可用于復(fù)填子區(qū)f的土方量，其值是開采子區(qū)m總巖土體積(最底層底板以上所有物質(zhì))與其可采礦層總物質(zhì)體積之差；為土體修正系數(shù)，是煤層受剝離損失、夾矸、巖土膨脹等多因素綜合影響下實際土方量與剝離土方量之比；，分別為采復(fù)子區(qū)內(nèi)網(wǎng)格單元的行、列序號，定義域分別為1～row與1～col；(m)為開采子區(qū)網(wǎng)格單元m(,)所在區(qū)域的采前地面平均標(biāo)高；(m)為開采子區(qū)網(wǎng)格單元m(,)所在區(qū)域最底層礦層的底板平均標(biāo)高；(m)為開采子區(qū)網(wǎng)格單元m(,)所在區(qū)域的可采礦層平均厚度；(m)為開采子區(qū)網(wǎng)格單元m(,)的水平投影面積，是行列分辨率的乘積，其余參數(shù)解釋同上。

2.1.3 復(fù)填區(qū)HPCS預(yù)構(gòu)建

自然穩(wěn)定地貌作為區(qū)域長時序演化的結(jié)果，其不僅形態(tài)穩(wěn)定，且地表水文特征與上下流域間高度匹配。為構(gòu)建與周邊自然地貌相銜接的重塑地貌，依據(jù)NbS理論，應(yīng)對其水文特征進(jìn)行保留，以實現(xiàn)重塑地貌與周邊自然地貌相融合。HPCS模型是在上述理論基礎(chǔ)上，通過對水文溝道和坡頂進(jìn)行識別，令其中水文溝道和邊界處高程不變，以保證地貌重塑前后內(nèi)外水系相通，匯水面關(guān)系不變；通過調(diào)整區(qū)域內(nèi)各坡頂高程，以實現(xiàn)土方量平衡；并通過目標(biāo)篩選，以獲取最優(yōu)地貌設(shè)計形態(tài)的地表高程調(diào)整模型。在數(shù)值上，其是重塑地貌DEM較采前自然地貌DEM高程變化的空間分布數(shù)據(jù)。

在HPCS求解過程中，首先需構(gòu)建一批形態(tài)各異的HPCS預(yù)構(gòu)建模型，以為后續(xù)土方控制和坡度緩和優(yōu)化提供可供操作的基底數(shù)據(jù)。故基于其模型定義與NbS理論，為盡可能保留復(fù)填子區(qū)f原始水文特征，通過ArcGIS水文分析模塊與Matlab軟件，創(chuàng)建圖4所示HPCS預(yù)構(gòu)建模型。

圖4 HPCS預(yù)構(gòu)建示意Fig.4 Schematic diagram of pre-construction of HPCS

如圖4所示，對原有自然穩(wěn)定地貌進(jìn)行概化，以減少NNDL_HPCS求解結(jié)果與PNL間的匯水面分布差異。其中，依據(jù)概化要素于原有地貌水文所起作用，劃分為分割匯水面的地形突出點(坡頂點)和匯水面內(nèi)水流集中的溝道線(水系或者溝道)。其中，水文溝道網(wǎng)格提取結(jié)果通過ArcGIS水文分析模塊獲取，其最小匯水面積參數(shù)通過均值變點法(OTSU)確定。原始地貌的地形突出點則基于式(2)提?。?/p>

(2)

其中，(,)為復(fù)填子區(qū)網(wǎng)格單元(，)的取值，凸點取值為1，其余為0；area為以網(wǎng)格單元(，)為中點的*高程搜索窗口，area(2，2)為窗口中心網(wǎng)格單元(，)所對應(yīng)的高程值，其余參數(shù)解釋同上。當(dāng)搜索網(wǎng)格單元中心點單元area(2，2)高程為其窗口內(nèi)最大值，定義該中心網(wǎng)格為凸點。通過Matlab遍歷填挖后復(fù)填子區(qū)PNL，識別凸點網(wǎng)格單元所在位置。在此基礎(chǔ)上，將復(fù)填子區(qū)凸點網(wǎng)格識別結(jié)果與水文溝道網(wǎng)格提取結(jié)果進(jìn)行空間疊加，構(gòu)建式(3)判別式對HPCS網(wǎng)格單元進(jìn)行預(yù)賦值：

(3)

其中，為HPCS網(wǎng)格單元預(yù)賦值結(jié)果;(,)為其網(wǎng)格單元(,)的取值，凸點網(wǎng)格賦值為w，其余為0；w為凸點網(wǎng)格(，)所對應(yīng)的無量綱系數(shù)，定義域為[0,1]，其余參數(shù)解釋同上。在此基礎(chǔ)上，基于式(4)所示模糊算法，預(yù)構(gòu)建HPCS：

(4)

2.1.4 曲面土方量控制及坡度緩和優(yōu)化

同SAS中的土方動態(tài)平衡方法，通過豎直伸縮變化使預(yù)構(gòu)建的HPCS符合土方平衡要求，其計算式為

(5)

′(,)=(,)

(6)

′(f)=′(,)+(f)

(7)

式中，為伸縮變化系數(shù)，其值是復(fù)填子區(qū)f剝離的總巖土(包含礦層)體積減去由開采子區(qū)m獲取的膨脹后的可用土方量之后，與預(yù)構(gòu)建的HPCS總體積比值的相反數(shù)；′(,)為預(yù)構(gòu)建HPCS伸縮變換后網(wǎng)格單元(,)的值，是復(fù)填子區(qū)f最終設(shè)計結(jié)果(NNDL_HPCS)與PNL間高程值的變化；′(f)為NNDL_HPCS中復(fù)填子區(qū)網(wǎng)格單元f(，)所對應(yīng)的地表標(biāo)高，其值等于該處的自然原始標(biāo)高與伸縮變換后的HPCS值之和。

在此基礎(chǔ)上，為篩選出區(qū)域表土抗水蝕能力較強(qiáng)(坡度平緩)且與周邊地貌融合的NNDL_HPCS，構(gòu)建相鄰網(wǎng)格單元間坡角為對象的區(qū)域坡度評分標(biāo)準(zhǔn)。在坡度評分標(biāo)準(zhǔn)設(shè)置中，為避免區(qū)域不同坡度組合間的分值混淆，同時減少軟件運(yùn)算中非必要的內(nèi)存占用，以1°為步長，令相鄰步長間高坡度分值是低坡度分值乘上相鄰網(wǎng)格數(shù)(即8個)的值再加1。其計算式為

(8)

其中，()為坡度夾角(°)所對應(yīng)的坡度評分值，其值域∈；為正整數(shù)，定義域在2°～90°。在區(qū)域坡度評分標(biāo)準(zhǔn)構(gòu)建的基礎(chǔ)上，構(gòu)建區(qū)域坡度評分模型，其計算式為

(9)

其中，(f)為復(fù)填子區(qū)f的區(qū)域坡度總評分值；p()為子區(qū)f內(nèi)網(wǎng)格單元(,)沿某一朝向的坡度評分值;∑p為子單元(,)與周邊8個網(wǎng)格單元間坡角所對應(yīng)評分值的總和。為確保最終篩選對象與周邊已知地貌無縫拼接，除復(fù)填子區(qū)f與礦坑鄰接處網(wǎng)格單元外，其余網(wǎng)格單元都需統(tǒng)計與周邊8個網(wǎng)格單元的坡度評分值。

將設(shè)為w的最小變化步長，通過上述HPCS預(yù)構(gòu)建模型，由0至1遍歷獲取一批不同w組合下的預(yù)構(gòu)建HPCS，進(jìn)而獲取一組形態(tài)各異的復(fù)填子區(qū)′。在此基礎(chǔ)上，以區(qū)域坡度評分最小為目標(biāo)，定向篩選復(fù)填子區(qū)′求解結(jié)果，最終獲取原有自然形態(tài)表面下，土方量動態(tài)平衡且坡度平緩的復(fù)填子區(qū)f內(nèi)NNDL_HPCS。

2.2 水文保留率

坡向是地貌的重要表征參數(shù)，對地表徑流起著直接導(dǎo)向作用，同時受區(qū)域光照、溫度、雨量、風(fēng)速等因子的綜合影響，間接決定區(qū)域植被類型及其生長狀況?；诨疑到y(tǒng)理論，以自然穩(wěn)定地貌坡向空間特征為塑形導(dǎo)向。將水文保留率定義為重塑地貌與采前自然地貌坡向空間分布上一致性，即2者相同坡向區(qū)域的面積占總區(qū)域面積的比例。通過ArcGIS坡向組件，提取PNL,NNDL_HPCS,NNDL_SAS和TDL的坡向空間柵格數(shù)據(jù)。以PNL坡向空間柵格數(shù)據(jù)為基準(zhǔn)，通過Matlab軟件逐柵格進(jìn)行空間疊加分析，計算NNDL_HPCS,NNDL_SAS,TDL三者與PNL的區(qū)域坡向一致性，其計算式為

(10)

其中，為地貌相對于PNL的水文保留率，其定義域在[0,1]，若2者坡向空間一致性越高(即相同位置2者間的坡向夾角越小)，其取值越大；row，col分別為地貌坡向柵格數(shù)據(jù)的總行、總列數(shù)；為地貌的坡向空間柵格數(shù)據(jù)，其中(,)為其第行第列的柵格的坡向值，其取值大小符合羅馬盤規(guī)則；為PNL的坡向空間柵格數(shù)據(jù)，其中(,)是其第行第列柵格的坡向值，其取值大小也符合羅馬盤規(guī)則。在ArcGIS坡向模塊中，平臺坡向賦值為-1，考慮到平臺坡向的萬向性(可視作任意方向)，默認(rèn)將平臺處的坡向認(rèn)定為一致。

2.3 土方直線運(yùn)距

基于露天采排復(fù)一體化理念，令開采工作面土方由上層臺階逐層向下剝離，后側(cè)復(fù)填工作面由煤層剝離后底板處逐層至上復(fù)填。其中，復(fù)填區(qū)當(dāng)且僅當(dāng)相鄰低處臺階的可用土方量(考慮膨脹系數(shù))完全滿足復(fù)填所需時，才進(jìn)行相鄰高處臺階的復(fù)填作業(yè)。由于實際采排復(fù)一體化工作中，土方運(yùn)移受運(yùn)輸通道等諸因素影響，其難以精確計算，故依據(jù)上述思路，按照研究區(qū)邊幫設(shè)計參數(shù)，以10 m為分段間隔，通過Matlab建立采復(fù)前后土方的空間對應(yīng)索引，并依據(jù)其每m土體中心點的位移，計算其土方平均運(yùn)距，以此簡要評價NNDL_HPCS,NNDL_SAS和TDL的土方調(diào)配特征。其計算式為

(11)

其中，為土方平均運(yùn)距；為運(yùn)移土方按時間的排序，取值為以1為步長從0至(m)的遞進(jìn)值；(m,)，′(m,)分別為開采子區(qū)m中第運(yùn)移的立方土體運(yùn)移前后中心點所對應(yīng)的三維坐標(biāo)，其余參數(shù)解釋同上。通過統(tǒng)計各采復(fù)周期中各立方土體運(yùn)移前后位移所對應(yīng)的模長均值，簡要評價重塑地貌土方運(yùn)移特征。

2.4 抗水蝕能力評價

CLiDE模型是基于CAESAR模型開發(fā)的景觀地貌演化模型。相比于傳統(tǒng)的侵蝕預(yù)測模型，在獲取侵蝕和沉積的時序數(shù)據(jù)的同時，還能獲取最終發(fā)育的地貌形態(tài)。并可模擬地表徑流-地貌-地下水3者間的交互作用和補(bǔ)給關(guān)系，極其適用于地貌修復(fù)效果評價。為評價內(nèi)排土場設(shè)計結(jié)果的抗水蝕能力，將研究區(qū)NNDL_HPCS設(shè)為實驗組，NNDL_SAS,TDL和PNL設(shè)為對照組，套用前置研究設(shè)定參數(shù)(氣象數(shù)據(jù)來源于中國氣象網(wǎng))構(gòu)建CLiDE模型，模擬4者10 a的土壤侵蝕過程。其中，依據(jù)研究區(qū)典型草原植被類型，將土地利用類型統(tǒng)一設(shè)置為第5類(適用于半野生草地、粗糙的草地、中性草地等)，并基于區(qū)域植被以1 a生為主，將其植被恢復(fù)時間統(tǒng)一設(shè)置為1 a。并基于Matlab軟件統(tǒng)計各組地貌演變前后土方變化量，以評價其抗水蝕效果。其中為排除礦坑最終位置對土壤水蝕模擬結(jié)果的干擾，在PNL對應(yīng)位置處設(shè)置相同大小的礦坑。其計算式為

第一，專業(yè)人才匱乏。當(dāng)下養(yǎng)老服務(wù)存在的主要困境是養(yǎng)老專業(yè)人才隊伍的匱乏，導(dǎo)致服務(wù)專業(yè)水平不足，影響服務(wù)質(zhì)量和服務(wù)成效。在養(yǎng)老服務(wù)提供方的人才隊伍中，主要包括志愿者、居委會工作人員，以及非專業(yè)出身的其他領(lǐng)域轉(zhuǎn)行而從事養(yǎng)老服務(wù)的人員等，在養(yǎng)老服務(wù)的能力和水平等方面都有待考量，這對于前期智慧養(yǎng)老服務(wù)模式的建設(shè)是一個亟待突破的限制性因素。更值得注意的是，智慧養(yǎng)老的運(yùn)營模式不僅需要的是專業(yè)的養(yǎng)老服務(wù)人才，更是需要懂大數(shù)據(jù)技術(shù)、大數(shù)據(jù)工具的數(shù)據(jù)專業(yè)領(lǐng)域的養(yǎng)老服務(wù)復(fù)合型人才。

(12)

式中，為土壤水蝕總量；(f)，(f)分別為復(fù)填子區(qū)網(wǎng)格單元f(,)中地貌演化前后的地表標(biāo)高均值；(f)為復(fù)填子區(qū)網(wǎng)格單元f(,)的水平投影面積，為行列分辨率的乘積。

3 結(jié)果與分析

3.1 數(shù)據(jù)準(zhǔn)備結(jié)果及分析

研究區(qū)可采煤層共計9層，依據(jù)儲量報告，采用ArcGIS克里金插值模塊獲取研究區(qū)最底層煤層底板標(biāo)高空間分布圖(圖2(j)，對應(yīng)研究區(qū)13-2號煤層))及其各煤層厚度平面分布圖(圖2(a)～(i)，分別對應(yīng)研究區(qū)6，7，8，9，11，12-1，12-2，13-1，13-2號煤層)。為避免數(shù)據(jù)空間分辨率差異造成的干擾，將影像空間分辨率統(tǒng)一設(shè)置為30 m×30 m，數(shù)據(jù)形式均為162列32行的柵格影像。

考慮到研究區(qū)可采煤層非惟一(9層)，且其厚度與底板標(biāo)高存在差異，故依據(jù)三維地質(zhì)資料，以13-2號煤層的底板標(biāo)高為研究區(qū)煤層底板標(biāo)高，并采用空間疊加方法疊加各煤層厚度數(shù)據(jù)，以獲取研究區(qū)煤層頂板標(biāo)高等效數(shù)據(jù)，實現(xiàn)本研究地貌重塑模型于多礦層復(fù)雜情景下的應(yīng)用。

3.2 近自然設(shè)計地貌構(gòu)建結(jié)果及分析

依據(jù)2.1.1節(jié)采復(fù)子區(qū)空間位置識別方法，將研究區(qū)劃分為45個一一對應(yīng)的采復(fù)子區(qū)，與已有研究保持一致。其中，采復(fù)子區(qū)網(wǎng)格行(row)列(col)數(shù)均為32行3列，其水平投影為東西長90 m，南北長960 m的矩形。通過2.1節(jié)內(nèi)排土場近自然地貌重塑模型，將w的最小變化步長取0.01，土體修正系數(shù)取1.1，求解研究區(qū)NNDL_HPCS(圖5(b))。

圖5 研究區(qū)內(nèi)排土場DEMFig.5 DEM of the internal dump in the study area

3.2.1 地貌融合結(jié)果與分析

在視覺效果上，相較于TDL(圖5(d))“階梯狀”地貌空間形態(tài)，NNDL_HPCS(圖5(b))和NNDL_SAS(圖5(c))兩者與PNL(圖5(a))較為相似，存在復(fù)雜的自然地貌紋理特征；在表土高程變化上，受土方平衡約束，NNDL_HPCS,NNDL_SAS,TDL三者較PNL高程變化幅度分別為-112～171,-205.90～50.03,-160.21～73.60 m，變化區(qū)間大小分別為NNDL_HPCS(283 m) > NNDL_SAS(250.93 m)>TDL(233.81 m)；在空間形變特征上，受“開采子區(qū)地表原始標(biāo)高”、“復(fù)填子區(qū)煤層底板標(biāo)高”、“復(fù)填子區(qū)煤層頂板等效標(biāo)高”、“膨脹系數(shù)”4個要素的直接影響，本研究NNDL_HPCS與已有成果NNDL_SAS均較PNL發(fā)生空間上的不規(guī)則形變。例如，在東西方向上，復(fù)填子區(qū)～NNDL均較PNL有所抬升，而在其他區(qū)域均較原始地貌有所下降。與前兩者相似，TDL較PNL也具有明顯的空間形變特征(圖5(d)，6(d))，然而與前兩者不同，受斜坡設(shè)計影響，TDL與周邊原始地貌間“分割”明顯，具“斜坡-平臺-斜坡”3段式的地貌形態(tài)特征(圖5(d)，6(d))。

如圖6所示，以研究區(qū)西邊界為起點，以PNL為基準(zhǔn)，由西向東逐列分別統(tǒng)計NNDL_HPCS,NNDL_SAS,TDL相對于PNL的高程方差。整體上，NNDL_HPCS,NNDL_SAS和TDL三者相對PNL的高程總方差分別為288 942.8,582 248.8和465 452.8 m。其中，NNDL_HPCS相對于PNL的高程方差最小，相較NNDL_SAS與TDL分別低50.37%與37.92%，在絕對高程上更接近于PNL，形態(tài)保留效果最優(yōu)。隨距研究區(qū)西邊界距離變化，在0～810 m段，NNDL_HPCS,NNDL_SAS,TDL,PNL四者均位于礦坑位置，列相對高程方差值均為0；在810～1 800 m段和4 590～4 890 m段，NNDL_HPCS列相對高程方差遠(yuǎn)高于NNDL_SAS和TDL，其均值分別為后2者均值的11.83倍與6.99倍，和17.16倍與6.43倍；在1 800～4 590 m段，NNDL_HPCS列相對高程方差值則遠(yuǎn)低于NNDL_SAS和TDL，其均值分別為后兩者均值的0.13倍與0.18倍。由地質(zhì)資料可知，在810～1 800 m段和4 590～4 890 m段，隨著前側(cè)開采子區(qū)煤層底板走低，可采煤層厚度變化較少，導(dǎo)致后側(cè)復(fù)填可用土方量不斷增加，NNDL_HPCS,NNDL_SAS和TDL對應(yīng)復(fù)填子區(qū)均形成“凸包地形”；在1 800～4 590 m段，隨著前側(cè)開采子區(qū)煤層底板走低，可采煤層厚度增大，導(dǎo)致后側(cè)復(fù)填可用土方量穩(wěn)中有減，NNDL_HPCS,NNDL_SAS和TDL對應(yīng)復(fù)填子區(qū)均形成“凹陷地形”(圖6)。即在“凸包地形”段，NNDL_HPCS與PNL柵格絕對高程重合效果上遠(yuǎn)弱于NNDL_SAS和TDL，在“凹陷地形”段則遠(yuǎn)優(yōu)于2者。

圖6 研究區(qū)PNL,NNDL_HPCS,NNDL_SAS,TDL三維示意(a)～(d)，及其復(fù)填子區(qū)物質(zhì)等效厚度隨距研究區(qū)西邊界距離的 折線統(tǒng)計(e)和3者相對于PNL的列高程方差隨距研究區(qū)西邊界距離的折線統(tǒng)計(f)Fig.6 Schematic diagram of the PNL，NNDL_HPCS，NNDL_SAS,TDL(a)～(d)，the broken line statistics chart of the equivalent thickness of the material in the refill subregion with the distance from the western boundary(e) and they relative to the PNL column elevation variance with the distance from the western boundary in the study area(f)

圖7 研究區(qū)NNDL_HPCS,NNDL_SAS,TDL水文保留 率隨距研究區(qū)西邊界距離的折線統(tǒng)計Fig.7 Hydrological retention rate of NNDL_HPCS，NNDL_SAS， and TDL in the study area is broken-line statistics with the distance from the western boundary of the study area

如圖7所示，以研究區(qū)西邊界為起點，以PNL為基準(zhǔn)，由西向東逐列分別統(tǒng)計NNDL_HPCS,NNDL_SAS,TDL相對于PNL的水文保留率。整體上，NNDL_HPCS,NNDL_SAS和TDL三者水文保留率分別為73.53%,70.68%和65.00%(不包括礦坑部分)。其中，NNDL_HPCS水文保留率最高，相較NNDL_SAS與TDL分別高2.85%與8.53%，在邊坡朝向上更接近于PNL，水文保留效果相較最優(yōu)。隨距研究區(qū)西邊界距離變化上，與高程方差類似，在0～810 m段，NNDL_HPCS,NNDL_SAS和TDL三者均位于礦坑位置，水文保留率均為1；在810～1 800 m段和4 590～4 890 m段，NNDL_HPCS列水文保留率(59.19%和60.60%)低于NNDL_SAS和TDL，其均值分別較后2者均值低23.27%與9.31%，和11.84%與5.23%；在1 800～4 590 m段，NNDL_HPCS列水文保留率(79.44%)高于NNDL_SAS和TDL，其均值分別較后者均值高13.75%和16.52%。與上述高程方差類似，在“凸包地形”段，NNDL_HPCS與PNL柵格水文保留率弱于NNDL_SAS和TDL，在“凹陷地形”段則優(yōu)于兩者(圖6)。

在土方運(yùn)移上，依據(jù)采復(fù)子區(qū)土方對應(yīng)關(guān)系，運(yùn)用Matlab計算獲取NNDL_HPCS,NNDL_SAS及TDL三者的土方平均運(yùn)距(PNL不存在采挖，故不進(jìn)行計算)。計算結(jié)果表明，3者土方平均運(yùn)距NNDL_HPCS(863.23 m/m)>NNDL_SAS(822.20 m/m)>TDL(821.71 m/m)。其中，NNDL_SAS土方平均運(yùn)距較TDL僅增加約0.06%(0.49 m/m)，NNDL_HPCS則較前者約增加5.05%(41.03 m/m)，其施工成本勢必較前2者有所增加。

3.2.3 土壤抗水蝕能力求解結(jié)果與分析

地貌穩(wěn)定性主要受氣候、地質(zhì)、植被等諸多要素的影響。其組成物料強(qiáng)度直接影響邊坡的穩(wěn)定性，是其抗侵蝕能力的主要指標(biāo)。根據(jù)地質(zhì)調(diào)查報告，確定研究區(qū)內(nèi)排土場結(jié)構(gòu)及巖土物理力學(xué)參數(shù)，并結(jié)合中國天氣網(wǎng)(http://www.weather.com.cn/)獲取研究區(qū)2007—2017年逐日氣象參數(shù)。通過CLiDE演化模型，分別模擬NNDL_HPCS,NNDL_SAS,TDL,PNL四者10 a的土壤水蝕過程。其中，為保證模擬結(jié)果的可比性，使NNDL_HPCS,NNDL_SAS和TDL三者的采復(fù)子區(qū)空間位置保持一致，最終輸入范圍如圖6所示。同時，為避免邊緣效應(yīng)對模型模擬結(jié)果的影響，將研究區(qū)周邊1 000 m范圍共同納入運(yùn)算(圖6)。最終通過Matlab軟件進(jìn)行疊加對比分析，統(tǒng)計研究區(qū)內(nèi)排土場NNDL_HPCS,NNDL_SAS,TDL和PNL四者10 a的土壤水蝕總量(表1)。

表1 研究區(qū)NNDL_HPCS,NNDL_SAS,TDL,PNL 土壤10 a總水侵蝕量Table 1 10-year total soil water erosion scale of the NNDL_HPCS，NNDL_SAS，TDLand PNL of the dump site

由表1可知，整體上，本研究內(nèi)排土場10 a土壤水蝕總量在4.098 3×10～23.681 5×10m。其中，PNL土壤抗水蝕能力最強(qiáng)，模擬情景下其土壤10 a水蝕總量僅占4者平均水平(12.175 6×10m)的33.66%。在人工地貌設(shè)計下，本研究NNDL_HPCS和已有研究NNDL_SAS相較TDL分別可減少約39.07%和72.58%的土壤水侵總量，內(nèi)排土場土壤抗水蝕能力提升明顯。并且，相較于已有研究中NNDL_SAS土壤抗水蝕效果，本研究NNDL_HPCS可進(jìn)一步減少約54.99%的土壤水蝕總量。

如圖8，9所示，以研究區(qū)內(nèi)排土場西邊界為起點，以PNL10 a土壤水蝕為基準(zhǔn)(將其列統(tǒng)計值定義為對應(yīng)列的0值)，由西向東逐列分別統(tǒng)計NNDL_HPCS,NNDL_SAS,TDL相對于PNL的10 a列土壤水蝕量。與水文保留率和絕對高程方差的統(tǒng)計結(jié)果不同，在810～1 800 m段，10 a土壤水蝕量中NNDL_SAS(-6.78×10m)

圖8 NNDL_HPCS,NNDL_SAS,TDL相對于PNL的土壤水蝕量隨距研究區(qū)西邊界距離變化Fig.8 Soil water intrusion of NNDL_HPCS，NNDL_SAS，and TDL relative to PNL varies with the distance from the west boundary of the study area

圖9 NNDL_HPCS,NNDL_SAS,TDL和PNL的 土壤水蝕量隨距研究區(qū)西邊界距離變化直方圖Fig.9 Histogram of soil water erosion of NNDL_HPCS， NNDL_SAS，TDL and PNL with the distance from the western boundary of the study area

即整體上，研究區(qū)PNL在全段均發(fā)生較低水平的土壤水蝕；NNDL_SAS土壤水蝕主要發(fā)生在初始段(4 590～4 890 m段)與中段(1 800～4 590 m段)，土壤堆積主要發(fā)生在終段(810～1 800 m段)；NNDL_HPCS土壤水蝕主要發(fā)生在起始段與終段，土壤堆積主要發(fā)生在中段；TDL則在全段均發(fā)生較高水平的土壤水蝕。

4 討 論

煤炭作為我國重要的能源礦產(chǎn)，其露采比例隨社會經(jīng)濟(jì)的發(fā)展不斷增加。內(nèi)排土場作為煤炭露采區(qū)的重要組成部分，其地貌形態(tài)不僅直接影響礦區(qū)人員設(shè)施及生產(chǎn)安全，也直接和間接影響區(qū)域生態(tài)環(huán)境系統(tǒng)穩(wěn)定。受傳統(tǒng)復(fù)填方式影響，內(nèi)排土場往往表現(xiàn)為“梯田式”人工堆墊地貌，從而使原地表要素及水文特征喪失，并在集中降雨時容易引發(fā)嚴(yán)重的表土水蝕等問題。諸多實踐與研究均表明，NbS理念下基于自然穩(wěn)定特征的排土場地貌重塑具有重要的理論與實踐意義，其在有效提高區(qū)域表土穩(wěn)定性的同時所需維護(hù)費(fèi)用較低。與上述結(jié)論相同，團(tuán)隊前置研究均表明，排土場近自然地貌的表土穩(wěn)定性優(yōu)于傳統(tǒng)“梯田式”地貌。且前置研究中，基于SAS的采排復(fù)一體下內(nèi)排土場NNDL較TDL 10 a水蝕量可減少約39.07%，初步驗證了NNDL_SAS在理論上具備可行性。在此基礎(chǔ)上，本研究進(jìn)一步提出了基于HPCS的內(nèi)排土場地貌近自然重塑方法。與前者類似，本文所提方法綜合考慮了采排復(fù)一體化下土方平衡約束、區(qū)域坡度平滑控制與邊界地貌融合3方面約束，在地貌融合效果上不僅優(yōu)于NNDL_SAS，還進(jìn)一步提高了研究區(qū)的抗水蝕能力(10 a土壤水蝕量分別較NNDL_SAS和TDL減少約55%和72.58%)。然而與前者不同的是，本研究基于HPCS的內(nèi)排土場地貌近自然重塑方法不再僅著眼于地貌重塑結(jié)果與周邊自然地貌間的銜接，而是致力于將復(fù)填可用土方再分配于原有的各個水文細(xì)小單元之中，在灰色系統(tǒng)理論支撐下通過保留或微調(diào)整此類要素以盡可能維護(hù)區(qū)域原有生態(tài)系統(tǒng)。

在露天礦采排復(fù)一體化過程中，已有坡面重塑方法能夠較好地應(yīng)用于傳統(tǒng)內(nèi)排土場復(fù)填區(qū)地貌重塑，并在縮短土方運(yùn)距、減小占地面積、提升技術(shù)適用性等方向得到了充分發(fā)展。然而在地貌形態(tài)上，上述內(nèi)排土場多呈現(xiàn)出“棱角分明”的“梯田式”空間布局，存在邊坡表土易發(fā)侵蝕，景觀明顯破碎化等問題。雖然后續(xù)部分學(xué)者提出在原有施工工藝的基礎(chǔ)上通過“弱化”排土場斜坡與平臺間“分界”的方式提高區(qū)域表土穩(wěn)定性及其融合效果，但實際上仍是“去除棱角”的“臺階狀”地貌，未有效解決上述問題。并且，受采復(fù)周期、可采煤層厚度等諸要素的影響，復(fù)填區(qū)可用土方隨采復(fù)周期動態(tài)變化，各周期間內(nèi)排形態(tài)存在差異。本研究和前置研究中，筆者均采用采復(fù)子區(qū)一一對應(yīng)的方式簡化采復(fù)過程中復(fù)填可用土方的動態(tài)變化問題，并相較于TDL在土方運(yùn)距小幅增加(0.06%～5.05%)的情況下大幅減少(39.07%～72.58%)表土水蝕量和獲得更高的地貌融合效果。然而值得注意的是，相較于傳統(tǒng)規(guī)則狀內(nèi)排土場地貌，NNDL形態(tài)較為復(fù)雜。在實際施工作業(yè)中，為避免近自然重塑地貌不規(guī)則起伏所帶來的安全風(fēng)險，應(yīng)在排土場次高層平盤上運(yùn)用卡車運(yùn)輸加單斗挖掘的形式進(jìn)行近自然地貌重塑，其余臺階均應(yīng)嚴(yán)格遵照傳統(tǒng)內(nèi)排土場地貌設(shè)計規(guī)范；同時應(yīng)充分運(yùn)用RTK等定位技術(shù)，確保排土場近自然重塑地貌特征要素(例如2.1.3節(jié)內(nèi)溝道、坡定點)位置符合設(shè)計，通過“點線構(gòu)面”的思路指導(dǎo)地貌近自然重塑。

受地貌發(fā)育影響，自然原始地貌并不總是穩(wěn)定的(例如黃土高原沖刷地貌)。在特定條件下，傳統(tǒng)復(fù)墾后礦區(qū)或?qū)⒏牧荚械孛病Ｓ捎诒狙芯縃PCS模型以PNL為學(xué)習(xí)對象，因此在使用該模型前首先需保證PNL地貌穩(wěn)定(例如本研究PNL 10 a的土壤侵蝕量最小，僅為TDL的28.4%)。并且，本研究各樣本10 a土壤水蝕模擬量隨距西邊界距離變化的統(tǒng)計結(jié)果表明，在局部地區(qū)基于PNL的NNDL構(gòu)建結(jié)果表土抗蝕能力甚至?xí)佑趥鹘y(tǒng)方式(圖9中，1 800～4 590 m段，NNDL_SAS劣于TDL；4 590～4 890 m段，NNDL_HPCS劣于TDL)。因此，需充分論證NNDL,TDL和PNL三者的穩(wěn)定性，從而確定有利于因地制宜的內(nèi)排土場地貌重塑方式。

5 結(jié) 論

(1) 地貌設(shè)計上，基于HPCS的內(nèi)排土場近自然地貌重塑模型可依據(jù)預(yù)設(shè)采復(fù)周期下采復(fù)子區(qū)采前地表、礦層頂板和礦層底板3者的DEM，自動求解復(fù)填可用土方量，獲取研究區(qū)NNDL_HPCS。

(2) 地貌融合上，NNDL_HPCS較NNDL_SAS和TDL視覺上更接近于PNL，其相對PNL的高程總方差較后2者分別低50.37%和37.92%，且其自然水文保留率較后2者分別高出2.85%和8.53%。

(3) 土方調(diào)配上，NNDL_HPCS土方運(yùn)移成本小幅增加，其土方平均運(yùn)距分別較NNDL_SAS和TDL分別高出約4.99%和5.05%。

(4)土壤抗水蝕上，NNDL_HPCS較NNDL_SAS 10 a可減少約55%的表土水蝕總量，較TDL減少約72.6%，其可顯著提高內(nèi)排土場表土抗水蝕能力，與PNL相逼近。