夏嘉南，李根生，卞正富，雷少剛，3，宮傳剛，李 恒

（1.中國礦業(yè)大學 礦山生態(tài)修復教育部工程研究中心，江蘇 徐州 221116；2.中國礦業(yè)大學 公共管理學院，江蘇 徐州 221116；3. 山東省煤田地質(zhì)局第二勘探隊 山東省采煤塌陷地和采空區(qū)治理工程研究中心，山東 濟寧 272000）

0 引 言

內(nèi)排土場作為露天礦區(qū)重要的組成部分，其占地面積隨采復工作面推進不斷增大[1]。在內(nèi)排過程中，傳統(tǒng)地貌重塑方法往往將原始微起伏地貌景觀轉(zhuǎn)變?yōu)槠脚_和邊坡相間的“梯田式”人工堆疊景觀，使原始地表要素及水系特征消失[2]。雖然在后續(xù)復墾設計中，內(nèi)排土場內(nèi)布設輸水管網(wǎng)以滿足復墾植被生長需求，但整體上未能實現(xiàn)其與周邊地貌融合，易引發(fā)土地退化等問題[3-4]。因此，需嘗試在采復過程中保留原始自然地貌特征以破解上述困境。

近年來，不少專家學者就內(nèi)排土場規(guī)劃及內(nèi)排方案開展了研究。例如，王東等[5]以內(nèi)排土場空間利用最大化為目標，基于采場工作線及邊坡角、內(nèi)排土場邊坡角、基底傾角對排土線的影響規(guī)律，構(gòu)建了采場工作幫、內(nèi)排土場、端幫和基底的空間幾何模型；栗嘉彬[6]以技術(shù)最佳適用性為目標，根據(jù)物料用途，構(gòu)建了露天礦表土最優(yōu)采集周期、最佳覆土厚度等求解模型；OGGERI等[7]以最小復填成本為目標，根據(jù)場地、運輸?shù)戎T因素，提出了多情境下露天礦采復處置方法。然而上述研究中，未充分考慮重塑地貌結(jié)果與周邊自然地貌起伏貼合情況，導致內(nèi)排土場往往被重塑為平臺加斜坡的“簡單梯田式”形態(tài)，人工地貌形態(tài)明顯，易引發(fā)嚴重的景觀破碎化問題[8-10]。相較于人工規(guī)則地貌，自然地貌作為長期演變的結(jié)果，具有高度的區(qū)域適宜性與穩(wěn)定性，可作為修復區(qū)地貌重塑的參考對象[3,10-11]。且國外實踐證明[12-13]，基于自然的地貌重塑能更好地適應當?shù)厮畾鈼l件，在有效提高區(qū)域抗水蝕能力的同時僅需要較少的維護費用[14]。因此，迫切需要一種貼合周邊自然地貌的內(nèi)排土場重塑方法，用于解決內(nèi)排土場設計結(jié)果與周邊自然地貌間貼合的問題[15]。

基于此，筆者以新疆黑山露天礦為例，提出一種采排復一體下基于調(diào)整曲面的內(nèi)排土場近自然地貌重塑模型，并結(jié)合地理空間分析技術(shù)與MATLAB軟件，構(gòu)建內(nèi)排土場近自然地貌重塑結(jié)果。并通過CLiDE景觀演化模型，模擬了研究區(qū)10 a間原始地貌、傳統(tǒng)設計地貌及近自然設計地貌三者的土壤水蝕情況，對比評價了近自然設計地貌抗水蝕效果，為相同及相似地區(qū)內(nèi)排土場地貌形態(tài)設計提供借鑒。

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

1.1 研究區(qū)概況

研究區(qū)位于新疆托克遜縣西北約90 km處，地理坐標：87°26′04″～87°29′40″E，43°13′48″～43°14′01″N，周邊地貌以山間谷地為主；年均降水量約160.3 mm，年均氣溫約1.6 ℃，冬寒夏涼，屬大陸干旱及高寒氣候。區(qū)內(nèi)表土以第四系為主，野生植被類型為典型草原植被。

依據(jù)地質(zhì)勘探資料，選取吉布拉克煤礦探礦權(quán)范圍以南為研究區(qū)（圖1）。區(qū)域東西長約4 874 m、南北長約1 000 m。其中，預設采坑邊界東西長約812 m、南北長約1 000 m，采礦作業(yè)由東至西；區(qū)域可采煤層分別為6、7、8、9、11、12-1、12-2、13-1和13-2號煤層，其中13-2號煤層為全區(qū)可露采、其余為局部可露采。

圖1 研究區(qū)示意Fig.1 Schematic of the study area

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

研究區(qū)采前自然原始數(shù)字高程數(shù)據(jù)（DEM）來源于地理空間數(shù)據(jù)云，拍攝時間2011年12月，條帶號87，行編號43，空間分辨率30 m×30 m；研究區(qū)煤層底板、頂板標高及水平坐標數(shù)據(jù)來源于《新疆托克遜縣黑山礦區(qū)黑山露天煤礦資源儲量核實報告》。此外，基于上述報告通過ArcGIS克里金插值模塊獲取煤層底板及頂板空間分布數(shù)據(jù)；內(nèi)排土場采、復子區(qū)空間分布數(shù)據(jù)基于內(nèi)排土場采復子區(qū)位置識別方法獲??；復填子區(qū)內(nèi)的可用土方量數(shù)據(jù)通過可用土方量計算式求解獲取；復填區(qū)地表調(diào)整曲面數(shù)據(jù)通過地表調(diào)整曲面預構(gòu)建模型獲??；復填子區(qū)地表近自然DEM通過MATLAB軟件基于動態(tài)土方量一致、且區(qū)域斜率較緩的原則定向篩選不同B樣條控制點位置下的地表調(diào)整曲面獲得；研究區(qū)自然原始DEM數(shù)據(jù)通過ArcGIS空間疊加提取獲得；研究區(qū)內(nèi)排土場傳統(tǒng)設計DEM數(shù)據(jù)基于土方量動態(tài)平衡下“30°斜坡加平臺”的設計思路獲取；最后運用CLiDE演化模型與MATLAB軟件，分別模擬并統(tǒng)計研究區(qū)原始DEM、傳統(tǒng)設計DEM及近自然設計DEM三者10 a的土壤水蝕量，以評估近自然整形結(jié)果的抗侵蝕能力。

2 研究方法

2.1 內(nèi)排土場近自然地貌重塑模型構(gòu)建

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

依據(jù)露天礦采復周期，獲取圖2所示各周期下開采子區(qū)amn和復填子區(qū)afn水平投影空間位置。

圖2 露天礦采復周期下開采及復填子區(qū)空間位置提取示意Fig.2 Schematic of extraction of spatial location of the mining sub-area and the refilling sub-area under the open-pit mining-refill cycle

如圖2所示，依據(jù)露天礦實際采復周期，沿開采方向?qū)⒙短斓V區(qū)按采復周期時間先后順序劃分為n個子區(qū)。其中，開采子區(qū)和復填子區(qū)分別對應區(qū)域amn和區(qū)域afn。在露天礦內(nèi)排構(gòu)建的采復過程中，采坑每前進一個開采子區(qū)amn對應修復后側(cè)一個復填子區(qū)afn，直至內(nèi)排過程結(jié)束，礦坑由區(qū)域Ab移至區(qū)域Aa。

2.1.2 可用土方量計算

為便于后續(xù)土方量計算，將開采及復填子區(qū)均細分為rowc行colc列的網(wǎng)格單元組（圖2）。在采復周期n中，開采子區(qū)amn可用于對應復填子區(qū)afn的土方量，計算式如下：

（1）

式中:V（amn）為開采子區(qū)amn可用于對應復填子區(qū)afn的土方量，其值是開采子區(qū)afn最下層可采礦層底板以上的總物質(zhì)體積減去區(qū)域內(nèi)所有可采礦層總物質(zhì)體積；i,j分別為開采及復填子區(qū)內(nèi)網(wǎng)格單元的行、列序號，取取值分別在1～rowc和1～colc之間；H（amnij）為開采子區(qū)amn網(wǎng)格單元（i,j）所在區(qū)域的采前地面標高；hd（amnij）為開采子區(qū)amn網(wǎng)格單元（i,j）所在區(qū)域的最底層礦層底板標高；ht（amnij）為開采子區(qū)amn網(wǎng)格單元（i,j）所在區(qū)域的可采礦層總厚度；S（amnij）為開采子區(qū)amn子單元（i,j）的水平投影面積，其余參數(shù)解釋同上。

2.1.3 復填區(qū)地表調(diào)整曲面預構(gòu)建

依據(jù)區(qū)域原有自然穩(wěn)定地貌，構(gòu)建圖3所示調(diào)整曲面模型，以盡可能保留復填子區(qū)afn原有地貌起伏特征。

圖3 地表調(diào)整曲面預構(gòu)建示意Fig.3 Schematic of pre-built adjustment surface

如圖3所示，將反映復填子區(qū)最終地貌設計形態(tài)較其原始自然地貌空間高程變化的曲面定義為地表調(diào)整曲面。為保證地表調(diào)整曲面能使復填子區(qū)表面與周邊自然地貌“無縫”融合，選取邊界處高程變化小，中間變化大，且與區(qū)域自然地貌相似的連續(xù)放緩曲面。參照“船體”構(gòu)建思路，采用連續(xù)性較優(yōu)的B樣條曲線[16]為“龍骨”，三角函數(shù)為側(cè)方“支架”，構(gòu)建復填子區(qū)afn所對應的地表調(diào)整曲面。圖中虛線箭頭為開采方向；x,y,z三軸正方向分別對應開采方向、水平垂直開采方向及豎直方向；fb和ft分別對應B樣條曲線函數(shù)和三角函數(shù)。其中，B樣條曲線fb位于復填子區(qū)afn水平投影中心線位置，采用控制點構(gòu)建，其計算式如下：

fb（j）=Bezier（v,tb,ts,te）

（2）

（3）

其中:fb（j）為復填子區(qū)afn第j列的B樣條曲線值，定義域在0～1；Bezier為豎直平面上的2維B樣條曲線函數(shù)，形態(tài)由起始點0,y1、終止點1,yNC和中間NC-2個控制點共同控制；NC為Bezier函數(shù)中起點和終點在內(nèi)所有控制點的個數(shù)，取值為>2的整數(shù)；v為NC行2列的矩陣，由上至下行數(shù)據(jù)分別對應起點、控制點1至控制點NC-2、終點的坐標；tb為起點所對應網(wǎng)格單元的列序號；ts為Bezier函數(shù)取值間隔，是終點與起點所對應網(wǎng)格單元的列序號之差，取值為大于0的正數(shù)；te是終點所對應網(wǎng)格單元的列序號；x1～xNC-2分別為控制點1至控制點NC-2在中心線上以起點為原點的水平投影距離，其中x1≤x2≤…≤xNC-2；y1～yNC分別是起點、控制點1至控制點NC-2、終點的預設高程值，定義域皆在0～1之間，其余參數(shù)解釋同上。在確定復填子區(qū)B樣條曲線函數(shù)形態(tài)的基礎上，通過三角函數(shù)確定地表調(diào)整曲面，其計算式如下：

（4）

式中:ft（i,j）為復填子區(qū)afn地表調(diào)整曲面網(wǎng)格單元（i,j）所對應的預設高程值，取值在0～1；π為圓周率，其余參數(shù)解釋同上。

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

為保證采復周期中可用土方量動態(tài)平衡，通過豎直伸縮變化使地表調(diào)整曲面符合土方平衡要求，其計算式如下：

（5）

（6）

（7）

在此基礎上，為篩選出區(qū)域整體坡度較緩的復填子區(qū)地表高程模型，以提高其表土抗水蝕能力，構(gòu)建以網(wǎng)格單元與周邊網(wǎng)格單元間坡角為對象的區(qū)域坡度評分標準，其計算式如下：

（8）

其中：fp（a）為坡度夾角a所對應的坡度評分值，其值為大于0的整數(shù)；au為正整數(shù)，其取值在2°～90°，其余參數(shù)解釋同上。在區(qū)域坡度評分標準構(gòu)建的基礎上，構(gòu)建區(qū)域坡度評分模型，其計算式如下：

（9）

式中：p（afn）為復填子區(qū)afn的區(qū)域坡度總評分；fpij（a）為子區(qū)afn內(nèi)網(wǎng)格單元（i,j）沿某一朝向的坡度評分值；∑fpij（a）為子單元（i,j）與周邊8個網(wǎng)格單元間坡角所對應評分值的總和。

依據(jù)上述評分模型，通過MATLAB軟件，調(diào)整2.1.3節(jié)中B樣條曲線起點、控制點1至控制點NC-2、終點的坐標位置，遍歷構(gòu)建一組形態(tài)各異的復填子區(qū)地表高程模型。并以區(qū)域坡度總評分最小為目標，定向篩選復填子區(qū)Ht′求解結(jié)果，最終獲取原有自然形態(tài)表面下，土方量動態(tài)平衡且坡度平緩的復填子區(qū)afn地表高程近自然設計形態(tài)。

2.2 設計結(jié)果抗水蝕能力評價

CLiDE模型是基于CAESAR模型開發(fā)的景觀地貌演化模型[17]。與傳統(tǒng)的侵蝕預測模型相比，其不僅能獲取侵蝕和沉積的時序數(shù)據(jù)，還能獲取最終發(fā)育的地貌形態(tài)。同時嵌入了Lisflood-FP流體力學相關模塊與地下水分布式水文模塊，可以模擬地表徑流-地貌-地下水三者間的交互作用和補給關系[18]，故極其適用于地貌修復效果評價[10]。為評價內(nèi)排土場設計結(jié)果的抗水蝕能力，將研究區(qū)內(nèi)排土場近自然設計結(jié)果為試驗組，原始自然地貌及傳統(tǒng)復填方法下內(nèi)排土場設計結(jié)果為對照組，結(jié)合實地調(diào)查情況構(gòu)建CLiDE模型，模擬其10 a的土壤侵蝕過程。并基于MATLAB軟件統(tǒng)計各組地貌演變前后土方變化量，以評價其抗水蝕效果。

3 結(jié)果與分析

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

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

圖4 13-2號煤層底板標高平面分布與研究區(qū)各煤層厚度平面分布Fig.4 Plane distribution map of floor elevation of No. 13-2 coal seam and plane distribution map of the thickness of each coal seam

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

依據(jù)2.1.1節(jié)采復子區(qū)空間位置識別方法，在30 m×30 m柵格影像數(shù)據(jù)下，按采復周期將研究區(qū)劃分為45個一一對應的采復子區(qū)。其中，采復子區(qū)網(wǎng)格行數(shù)為32，列數(shù)為3，即復填子區(qū)和開采子區(qū)均設計為水平投影東西長90 m，南北長960 m的水平矩形區(qū)域，以使采復工作前后采坑的水平投影大小不變（圖5）。

圖5 研究區(qū)開采子區(qū)及復填子區(qū)空間位置示意Fig.5 Schematic of the spatial location of the mining sub-area and the refill sub-area

在此基礎上，依據(jù)2.1節(jié)所示內(nèi)排土場地表近自然設計方法，通過MATLAB軟件以0.01為控制點數(shù)值最小變化間距，膨脹系數(shù)k取均值1.1，求解研究區(qū)內(nèi)排土場地貌近自然設計DEM結(jié)果（圖6）。其中，考慮到礦區(qū)可采煤層非惟一（9層），且厚度及底板標高存在差異，依據(jù)研究區(qū)三維地質(zhì)資料，以13-2號煤層的底板標高為研究區(qū)煤層底板標高，并采用空間疊加方法疊加各煤層厚度數(shù)據(jù)，獲取研究區(qū)煤層頂板標高等效數(shù)據(jù)，以實現(xiàn)本研究內(nèi)排土場近自然地貌重塑模型于多礦層復雜情景下的應用。

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

如圖6所示，視覺效果上，研究區(qū)內(nèi)排土場近自然設計DEM與原始自然地貌及紋理特征較為相似，且高程值均在+2 000.26～+2 959.00 m；表土高程上，受土方平衡約束，近自然設計DEM較自然原始DEM存在差異，變化幅度在-205.90～50.03 m。一方面，受研究區(qū)“開采子區(qū)地表原始標高”、“復填子區(qū)煤層底板標高”和“復填子區(qū)煤層頂板等效標高”的直接影響，近自然設計結(jié)果較原始地貌發(fā)生不規(guī)則形變。例如，在東西方向上，復填子區(qū)af1和af14近自然設計DEM較原始地貌有所抬升，而在其他區(qū)域較原始地貌有所下降。另一方面，受近自然地貌重塑模型中地表調(diào)整曲面模型影響（2.1.3節(jié)，地表調(diào)整曲面突出部分集中在區(qū)域中部），使得近自然設計DEM較原始地貌間形變也集中在區(qū)域中部。

與前者相似，相較于區(qū)域自然原始DEM，傳統(tǒng)設計結(jié)果也具有中部形變明顯的空間立體特征（圖6c、圖7c），然而與前者不同，受斜坡設計影響，內(nèi)排土場傳統(tǒng)設計DEM與周邊原始地貌間地貌“分割”明顯。例如，在同一位置由北至南分別做原始DEM、近自然設計DEM及傳統(tǒng)設計DEM三者的豎直剖面，在相同可用土方量控制下，相較于前兩者地貌“波動式”形態(tài)（圖7a、圖7b），傳統(tǒng)設計DEM具有明顯的“斜坡-平臺-斜坡”三段式的地貌形態(tài)特征（圖7c）。

在土方運移上，依據(jù)采復子區(qū)土方對應關系，運用MATLAB計算獲取近自然設計及傳統(tǒng)設計下土方平均運距。其中，近自然設計DEM與傳統(tǒng)設計DEM兩者土方平均運距分別為822.20 m/m3與821.71 m/m3，前者較后者僅高出約0.06%，兩者土方平均運距極為接近（相差0.26 m/m3）。

3.3 設計結(jié)果抗水蝕能力評價

地貌穩(wěn)定性主要受氣候、地質(zhì)、植被等諸要素影響。其組成物料強度直接影響邊坡的穩(wěn)定性，是其抗侵蝕能力的主要指標。依據(jù)地質(zhì)調(diào)查報告，確定研究區(qū)內(nèi)排土場地層結(jié)構(gòu)及巖土物理力學參數(shù)，并結(jié)合中國天氣網(wǎng)獲取研究區(qū)2007—2017年氣候參數(shù)。通過CLiDE演化模型，分別模擬3.2節(jié)所示10 a間內(nèi)排土場自然原始DEM、傳統(tǒng)設計DEM及近自然設計DEM的演變過程。其中，為保證模擬結(jié)果具有可比性，將研究區(qū)內(nèi)排土場傳統(tǒng)設計DEM的采復子區(qū)空間位置與3.2節(jié)近自然設計DEM設定為一致，其最終設計結(jié)果如圖7c所示。同時，為避免模型邊緣效應影響，將研究區(qū)周邊1 000 m范圍共同納入運算（圖7）。最后通過MATLAB疊加及對比分析，統(tǒng)計研究區(qū)內(nèi)排土場自然原始DEM、傳統(tǒng)設計DEM及近自然設計DEM下10 a的土壤水蝕總量（表1）。

如表1所示，整體上，本研究內(nèi)排土場10 a土壤水蝕量在4.098 3×105～23.681 5×105m3。其中，自然原始DEM土壤抗水蝕能力最強，模擬情境下其土壤10 a水蝕總量僅占樣本平均水平（14.069 5×105m3）的29.13%。在人工地貌設計下，本研究近自然設計DEM相較于傳統(tǒng)設計DEM可減少約39.07%的土壤水蝕量，其土壤水蝕總量逼近樣本平均水平。

為進一步分析排土場各種地貌設計方式下對區(qū)域整體土壤水蝕的影響，將各DEM演化模擬結(jié)果與其初始DEM求差，獲取研究區(qū)內(nèi)排土場表土高度的空間變化影像（圖8），其中正數(shù)為該區(qū)域表土高度增高高度，負數(shù)為降低高度。

圖7 研究區(qū)內(nèi)排土場自然原始、傳統(tǒng)設計及近自然設計下DEM影像及橫向剖面示意Fig.7 DEM image and cross-sectional schematic of the natural primitive, traditional design and near-natural design of the internal dump

表1 研究區(qū)內(nèi)排土場自然原始、傳統(tǒng)設計及近自然設計下區(qū)域土壤10 a總水侵蝕量

圖8 研究區(qū)內(nèi)排土場演化前后自然原始、傳統(tǒng)設計及近自然設計下DEM表層高度變化Fig.8 The surface height changes of the natural original DEM, the traditional design DEM, and the near-natural design DEM before and after the evolution of the dump

如圖8所示，整體上，各樣本土壤在10 a的演變模擬中均發(fā)生了不同程度的運移，使區(qū)域原有DEM發(fā)生不同程度的改變。在表土高度變化上，各DEM樣本高程變化區(qū)間均在-0.62～2.79 m之間，其變化幅度由大到小依次為：傳統(tǒng)設計DEM（3.39 m）、近自然設計DEM（3.35 m）和自然原始DEM（3.06 m）。排序結(jié)果與各樣本土壤10 a水蝕總量排序一致。在空間分布上，各樣本的土壤堆積及流失區(qū)域均有所不同。其中，傳統(tǒng)設計DEM表土流失及堆積區(qū)與研究區(qū)邊界高度相關，且主要沿邊界（圖8研究區(qū)邊界）呈向內(nèi)規(guī)則分布，視覺上與周邊自然地貌間具有明顯割裂，景觀及水文融合性較低。與前者不同，近自然設計DEM表土堆積及流失區(qū)空間分布與自然原始DEM較為相似，兩者均呈現(xiàn)出高度復雜的“離散化”空間分布。在周邊自然地貌融合上，兩者與周邊自然地貌間無明顯“痕跡”間隔，表土流失及堆積區(qū)域空間分布與研究區(qū)邊界無明顯相關，景觀及水文銜接性較好。并且，以自然原始DEM表土流失區(qū)（演化后高程下降區(qū)域）和堆積區(qū)（演化后高程上升區(qū)域）空間位置為參照的空間疊加分析表明，相較傳統(tǒng)設計DEM（空間重疊率：54.31%），近自然設計DEM土壤流失及堆積區(qū)具有更高的空間重疊率（57.27%）。

4 討 論

內(nèi)排土場作為露天礦山的重要組成部分，其占地面積隨工作面的推進不斷增大。然而，受傳統(tǒng)復填方式影響，內(nèi)排土場區(qū)域原始微起伏地貌往往被替換為斜坡加平臺的“梯田式”人工堆墊地貌，使原有地表要素及水文特征喪失，并在集中降雨條件下，易引發(fā)表土水蝕等問題[2-4]。為此，專家學者提出依據(jù)臨近自然穩(wěn)定地貌特征以實現(xiàn)排土場地貌重塑[3,8,10]。例如，有學者以勝利一號礦排土場邊坡為例，提出以仿自然“反S形”邊坡替代傳統(tǒng)斜坡，減少區(qū)域表土侵蝕[10]；同樣，有學者以勝利一號礦為例，提出依據(jù)內(nèi)排土場周邊自然地貌參數(shù)，基于區(qū)域整體土方平衡，構(gòu)建內(nèi)排土場近自然地貌結(jié)果，實現(xiàn)景觀內(nèi)外融合，且相較傳統(tǒng)設計地貌具備更好的表土穩(wěn)定性[3]。與上述研究相似，本研究同樣以自然地貌為學習對象，且內(nèi)排土場近自然設計地貌較傳統(tǒng)設計地貌具有更強的抗水蝕能力（10年土壤水蝕減少約39.07%）。然而不同的是，上述研究強調(diào)以區(qū)域地貌統(tǒng)計特征為依據(jù)，本研究則強調(diào)以原始自然穩(wěn)定地貌為依據(jù)，即在排土場塑形過程中盡可能保留原始自然穩(wěn)定地貌特征，以實現(xiàn)礦區(qū)內(nèi)外地貌相連。

在露天礦采復過程中，受土體膨脹系數(shù)、可采煤層標高、采復周期及采前地表標高等諸要素影響，使得各采復周期間后方復填可用土方量存在差異；且由于復填區(qū)形態(tài)各異，造成內(nèi)排土場復填地貌存在走高、走低等諸多情況。在近自然設計中，已有坡面近自然重塑方法雖然能夠較好地應用于傳統(tǒng)內(nèi)排土場復填區(qū)地貌設計，然而其實際形態(tài)仍以去除邊角的“梯田式”地貌為主[8,10]。后者提出的內(nèi)排土場近自然塑形方法雖具有良好的景觀融合效果，并考慮到開采前后土方平衡，然而設計中未考慮采復周期影響[3]，易造成采復周期中復填區(qū)實際可用土方量無法或超出原有設計表面等問題。受采復周期影響，本研究近自然設計地貌也呈現(xiàn)出部分區(qū)域較自然原始表面抬高、沉降等高程形變特征。然而與已有方法不同的是，本研究內(nèi)排土場近自然地貌重塑方法通過采復子區(qū)一一對應的方式解決采復過程中土方動態(tài)變化的問題，并相較于傳統(tǒng)設計地貌，近自然設計地貌與周邊自然地貌無明顯分界，模擬構(gòu)建上可實現(xiàn)采后地形與周邊自然地形的融合。

值得注意的是，在自然地貌中，區(qū)域原始地貌并不總是穩(wěn)定的[19]。在一定條件下，傳統(tǒng)復墾后礦區(qū)或?qū)⒏欣趨^(qū)域生態(tài)恢復[20]。由于本研究近自然地貌重塑方法以采前原始地貌為學習對象，故在構(gòu)建內(nèi)排土場地貌形態(tài)前，需充分論證原有自然地貌是否穩(wěn)定（例如，本研究自然原始地貌10年水蝕量最小，僅為傳統(tǒng)設計地貌的28.40%），以選擇有利于區(qū)域生態(tài)恢復的內(nèi)排土場設計方法。同時，受地下煤層剝離影響，本研究設計結(jié)果無法在形態(tài)上與原始自然地貌保持一致，并受采后松散土體空間不規(guī)則壓實沉降[21]，設計結(jié)果可能遠離預期。故在后續(xù)研究中，可進一步從力學視角，優(yōu)化本研究基于調(diào)整曲面的內(nèi)排土場近自然地貌重塑模型。

5 結(jié) 論

1）研究區(qū)內(nèi)排土場近自然設計結(jié)果表明，提出的近自然地貌重塑方法可在土方動態(tài)平衡的前提下，依據(jù)采復子區(qū)采前地表、礦層頂板及底板數(shù)字高程數(shù)據(jù)，獲取可用地表標高設計結(jié)果；在土方調(diào)配上，近自然設計地貌土方平均運距為822.20 m/m3，較傳統(tǒng)設計地貌高0.06%，兩者土方平均運距極為接近；且在視覺效果上，與周邊自然景觀“無痕”拼接，對區(qū)域開采前后地貌改變較少。

2）研究區(qū)10 a土壤水蝕模擬結(jié)果表明，本研究內(nèi)排土場近自然整形結(jié)果較傳統(tǒng)“梯田式”減少約39.07%的土壤流失量，表土抗水蝕能力較優(yōu)；在土壤演化運移中，本研究近自然設計地貌與區(qū)域原始地貌在土壤堆積及流失上空間重疊率達57.27%，較傳統(tǒng)設計結(jié)果高2.96%，且視覺上呈“離散化”分布，與周邊自然地貌銜接較好。