常一倫,謝婉麗,2,姬文龍,楊 惠,單 帥
(1.西北大學(xué) 地質(zhì)學(xué)系,陜西 西安 710069;2.西北大學(xué) 大陸動力學(xué)國家重點實驗室,陜西 西安 710069;3.陜西陜煤陜北礦業(yè)有限公司 生態(tài)環(huán)保部,陜西 榆林 719000)
隨著生態(tài)文明被寫入憲法,綠色礦山理念已上升為國家戰(zhàn)略[1-2]。大規(guī)模礦山活動在給我國經(jīng)濟社會注入新鮮活力的同時[3],不合理的礦山開發(fā)也破壞了礦區(qū)的生態(tài)與資源環(huán)境[4],誘發(fā)崩塌、滑坡、泥石流、地面塌陷、地裂縫等系列地質(zhì)災(zāi)害[5]。作為礦山災(zāi)害多發(fā)國之一,我國礦山地質(zhì)災(zāi)害有范圍廣、種類多、損失重、影響大等特點[6]。崩塌是指較陡斜坡上的巖土體在重力作用下以傾倒、滑動、墜落等突然的方式與母巖分離,并以滑動、滾動、跳躍等方式沿斜坡向下運動,最終在坡腳堆積的地質(zhì)現(xiàn)象[7-8]。其發(fā)生具有災(zāi)害頻發(fā)性、時間不確定性、發(fā)生快速、造成危害大等特點。目前,關(guān)于崩塌的研究主要集中在形成與破壞機理、勘察設(shè)計治理等方面[9-12]。在陜北開展煤礦開挖邊坡崩塌特征分析,預(yù)測崩塌的運動特征,對已建工程和擬建工程及區(qū)內(nèi)車輛、設(shè)備、人員的安全防護具有重要意義。
崩塌的研究基礎(chǔ)是現(xiàn)場調(diào)查,傳統(tǒng)現(xiàn)場調(diào)查方式多受地形、技術(shù)[13-14]、人員、成本[15-16]等條件制約。近幾年迅速發(fā)展的UAV技術(shù)具有高精度、低成本、易操作、數(shù)據(jù)易處理等優(yōu)點,很好解決了這一痛點問題[17]。葉震等[18]利用多旋翼無人機航測獲取了邊坡三維模型,并對模型進行球形k均值聚類分析,成功提取邊坡巖體結(jié)構(gòu)面的傾向、傾角參數(shù)。王頌等[19]現(xiàn)場調(diào)查結(jié)合無人機航拍,成功提取了西藏某崩塌三維數(shù)字模型和崩塌危巖區(qū)的巖體結(jié)構(gòu)面信息,進行了崩塌模擬,為沿線工程的崩塌防護提供了試驗依據(jù)。
對崩塌運動特征的研究主要利用數(shù)值模擬,常用的二維崩塌數(shù)值模擬軟件有ROCFALL和Rockfall Alalyst。ROCFALL可模擬崩塌運動路徑,提取高度、能量、速度和到達概率等,其應(yīng)用較為廣泛。SAGHIR[20]對Chalk cliffs進行了基于ROCFALL和現(xiàn)場調(diào)查的崩塌幾何尺寸、崩落點和崩落范圍分析。KISTER等[21]選取崩塌發(fā)育地區(qū)的路堤,進行了ROCFALL和Rockyfor3D的崩塌模擬,并進行了比較。Rockfall Alalyst是一款基于ArcObjects和C#的GIS擴展程序[22],可在ArcGIS下通過參數(shù)設(shè)置,利用無人機航測的高精度DEM影像底圖進行崩塌的2D模擬,并在ArcScene里3D可視化。仉義星等[23]利用Rockfall Alalyst模擬了崩塌區(qū)落石的路徑和危險性分析,所獲得的危險性物理模型結(jié)合易發(fā)性統(tǒng)計模型,完成了區(qū)域崩塌滑坡的易發(fā)性和危險性等級綜合分析。然而,ROCFALL是一款二維軟件,其剖面選取的人為性使得結(jié)果具有很大的偶然性。Rockfall Alalyst雖可以三維可視化,但落石的距離、角度間隔等,也都是人工設(shè)置后固定不變的。以上2種基于確定性的軟件都無法實現(xiàn)更接近真實的落石三維模擬運動。RocPro3D是一款專業(yè)三維崩塌模擬軟件,使用密集點云生成的三角網(wǎng)格作為地形底圖,可模擬落石在三維空間的運動,其功能涵蓋:地面模型生成、地面介質(zhì)屬性定義、落石參數(shù)定義與運動模擬、數(shù)值分析與成圖等,廣泛應(yīng)用于崩塌、泥石流等模擬與危險性評估中[24]。
以張家峁煤礦工業(yè)廣場開挖邊坡崩塌為研究區(qū),調(diào)查歷史崩塌的分布和危巖區(qū)危巖體特征,利用無人機進行航攝測量,獲取DEM、和點云數(shù)據(jù),進而獲取巖體結(jié)構(gòu)面信息。結(jié)合野外統(tǒng)計數(shù)據(jù)和無人機航拍,采用RocPro3D對危巖區(qū)進行三維數(shù)值模擬。旨在反演研究區(qū)崩塌模擬的各項參數(shù),預(yù)測落石的運動軌跡、影響范圍、彈跳高度、能量分布等,以期為崩塌治理、防護設(shè)計提供參考依據(jù)。
張家峁煤礦位于神木市北30 km,坐標N 39°1′15.17″,E 110°24′24.13″。場地建設(shè)過程中對植被的破壞和原地形開挖,加上差異風(fēng)化、卸荷作用等,崩塌災(zāi)害尤為嚴重。崩塌所在邊坡坡底標高+978 m,坡頂標高+1 070 m,走向355°,寬度約500 m,長度約140 m,坡度約35°,巖體破碎度高,節(jié)理裂隙發(fā)育,植被蓋度低?,F(xiàn)有防護措施包括局部主動防護網(wǎng)和被動防護網(wǎng)、擋墻、排水渠處理(圖1),為防止再次發(fā)生崩塌,進行了現(xiàn)場調(diào)查和無人機航拍,將東邊坡分為A~F 6處典型危巖區(qū)(圖2)。分別距坡腳80、80、67、55、27、31 m。坡體和坡腳均有不同數(shù)量的歷史崩塌堆積。
圖1 工業(yè)廣場東邊坡歷史工程治理措施Fig.1 Treatment measures of historical engineering for east slope of industrial site
圖2 工業(yè)廣場東邊坡崩塌概況Fig.2 Overview of collapse of east slope
東邊坡歷史上曾多次發(fā)生崩塌,坡腳和坡體上堆積了大量落石。據(jù)現(xiàn)場調(diào)查將坡腳和坡體上落石劃分為3個主要區(qū)域:北側(cè)、中南側(cè)和南側(cè)。
坡腳采用人工現(xiàn)場量測(圖3),坡頂?shù)入y以到達的地方借助高精度DOM(分辨率4.5 cm/pix)進行基于顏色、形狀、大小、紋理等的人機交互目視解譯(圖4),多次量測取均值方式,統(tǒng)計體積不小于0.1 m3的落石數(shù)量分別為:73、201、75塊(圖5)。按落石與邊坡的位置關(guān)系和粒徑區(qū)間來統(tǒng)計落石數(shù)量,結(jié)果如圖6所示。
圖3 落石粒徑現(xiàn)場調(diào)查Fig.3 Field survey of rockfall size
圖4 落石粒徑目視解譯Fig.4 Visual interpretation of rockfall size
圖5 歷史崩塌堆積分布Fig.5 Distribution of historical collapse
圖6 歷史崩塌落石尺寸-數(shù)量分布特征Fig.6 Size-quantity distribution characteristic of historical rockfalls
由圖6a知:落石崩落粒徑主要集中在0.1~5.0 m3,占落石總量的97.4%;中南側(cè)區(qū)域落石分布最多,共201塊,占比57.6%,北側(cè)、南側(cè)相差不大,占比20.9%和21.5%。由圖6b知:坡體上落石,崩落塊度在0.1~0.5 m3最為集中,占坡體上落石比64.0%,這是因為初始速度一致時,體積更小的落石崩落后動能小,易受坡體摩擦與障礙物阻擋而最先滯留坡體上部,體積中等的具有更大動能到達坡體下部,另外受崩塌物源的節(jié)理面間距、沉積巖層厚度的限制,體積大于5~10 m3的落石本身即較少,因此無論坡體上還是坡腳都數(shù)量最少。由圖6c知:坡腳落石,崩落塊度在0.1~1.0 m3最為集中,占比79.7%。崩塌的運動方式主要有滑動、滾動、躍動和自由崩落等,受限于節(jié)理面間距、沉積巖層厚度,體積大于5 m3的崩塌大多磨圓差,不易滾動,大概率只能以較小的初始速度滑動,難以到達坡腳。總體而言,落石崩落塊度分布區(qū)間,主要受控于結(jié)構(gòu)面數(shù)量、展布特征及間距;坡體上的具體概率分布,則受其粒徑大小、形態(tài)、磨圓、脆性、坡表性質(zhì)(障礙物特征、介質(zhì)特征)、坡度和微地形等共同作用。
坡腳和坡體上落石相對位置與粒徑大小的分布規(guī)律,可以為崩塌模擬中的地面介質(zhì)參數(shù)反演提供一定數(shù)據(jù)參考。
由于研究區(qū)所在的邊坡高度高,坡度大,難以人工進行結(jié)構(gòu)面產(chǎn)狀測量。利用無人機航攝和內(nèi)業(yè)處理得到三維點云數(shù)據(jù),通過Coltop軟件可獲取危巖體結(jié)構(gòu)面產(chǎn)狀(圖7)。
圖7 結(jié)構(gòu)面信息提取統(tǒng)計Fig.7 Structure information extraction statistics
調(diào)查發(fā)現(xiàn)危巖區(qū)主要受3組結(jié)構(gòu)面控制,從6處危巖區(qū)分別均勻選取若干產(chǎn)狀平直光滑的結(jié)構(gòu)面,每組結(jié)構(gòu)面提取200個點的XYZ坐標值,并通過計算得到每組結(jié)構(gòu)面的200個產(chǎn)狀,將這些產(chǎn)狀信息作圖成極射赤平投影極點圖(圖8)。
J1:結(jié)構(gòu)面1;J2:結(jié)構(gòu)面2;J3:結(jié)構(gòu)面3圖8 巖體結(jié)構(gòu)面極射赤平投影圖Fig 8 Polar stereographic projection of rock structure surface
如圖8所示,每組結(jié)構(gòu)面產(chǎn)狀數(shù)據(jù)相對集中,離散程度小。計算每組結(jié)構(gòu)面產(chǎn)狀的平均值,可獲得6處危巖區(qū)對應(yīng)的3組主要結(jié)構(gòu)面產(chǎn)狀平均值,見表1,并在圖8中作出3組結(jié)構(gòu)面產(chǎn)狀平均值對應(yīng)的極射赤平投影。
表1 危巖區(qū)巖體結(jié)構(gòu)面產(chǎn)狀統(tǒng)計
受巖體結(jié)構(gòu)面的切割,6處危巖區(qū)均為塊狀結(jié)構(gòu),結(jié)構(gòu)面間距直接控制了危巖體的塊度。通過無人機航攝處理獲取的正射影像圖,統(tǒng)計每組結(jié)構(gòu)面間距,并計算對應(yīng)各危巖區(qū)的危巖體積(表2)。
由表2知,E、F處危巖體塊度最大,A、B和D處較小,C處居中。由前面小節(jié)統(tǒng)計數(shù)據(jù)可知,北側(cè)、中南側(cè)、南側(cè)3個區(qū)域的坡體上、坡腳的巖體優(yōu)勢粒徑分別在0.1~0.5、0.1~0.5、0.1~5和0.5~5、0.1~1、0.1~0.5 m,而A~D 4處危巖物源區(qū)——除E、F 2處點狀危巖區(qū)(尚未發(fā)生崩塌)外,優(yōu)勢體積分別為2、2、7.5、3.75 m3。
表2 危巖區(qū)主要結(jié)構(gòu)面間距與體積
由圖9可見歷史崩塌體積在邊坡的整體分布主要集中在0.1~5.0 m3,與物源區(qū)危巖體結(jié)構(gòu)面控制下的平均體積具有很強關(guān)聯(lián)性,說明歷史崩塌最主要源自A~D危巖區(qū)。
研究歷史崩塌的體積和空間分布特征,可以在崩塌模擬的調(diào)參過程中為參數(shù)范圍的確定提供參考,而物源區(qū)危巖體的體積分布與展布特征,可以為數(shù)值模擬中物源區(qū)幾何參數(shù)的選取提供更可靠的演示數(shù)據(jù)。
圖9 歷史崩塌優(yōu)勢體積在坡體的分布概況Fig.9 Distribution of dominant volume of historical collapse on slope
危巖體發(fā)生破壞的原因主要是由于主控結(jié)構(gòu)面的端部巖石破裂所致[25]。目標區(qū)危巖體主要受控3組結(jié)構(gòu)面而被切割,且每組結(jié)構(gòu)面方向與巖體臨空面方向一致,因此可根據(jù)巖體結(jié)構(gòu)面之間的切割關(guān)系,初步進行危巖體的穩(wěn)定性判斷。如圖8所示,如選取J1結(jié)構(gòu)面作為臨空面,則J2和J3交于A點,交線OA傾向與J1結(jié)構(gòu)面傾向一致,交線傾角大于J1結(jié)構(gòu)面傾角,說明危巖區(qū)在J1臨空面是相對穩(wěn)定的;如選取J2組結(jié)構(gòu)面作為臨空面,則J1和J3交于B點,交線OB傾向與J2結(jié)構(gòu)面傾向一致,交線傾角小于J2組結(jié)構(gòu)面傾角,說明危巖區(qū)在J2臨空面是不穩(wěn)定的;如選取J3結(jié)構(gòu)面作為臨空面,則J1和J2交于C點,交線OC傾向與J3結(jié)構(gòu)面傾向一致,交線傾角小于J3組結(jié)構(gòu)面的傾角,說明危巖區(qū)在J3臨空面是不穩(wěn)定的。但由于沉積層差異風(fēng)化作用和其他次級結(jié)構(gòu)面,危巖體也存在以J1結(jié)構(gòu)面為臨空面失穩(wěn),發(fā)生崩塌的可能存在。
傾向J1臨空面的危巖體受巖體結(jié)構(gòu)面切割后呈塊狀、板狀,危巖體近乎水平,在危巖體下部抗風(fēng)化能力較弱巖層發(fā)生掏蝕后,容易發(fā)生錯斷式崩塌。傾向J2和J3臨空面的危巖體傾角多在65°~75°,受結(jié)構(gòu)面切割呈高角度傾斜塊狀,在重力或其他外力作用下,有沿著臨空面發(fā)生向下滑移、或沿著危巖體底部某一支點傾倒的趨勢。據(jù)現(xiàn)場調(diào)查歷史崩塌分布和無人機航拍,A~D 4處帶狀危巖區(qū)現(xiàn)狀欠穩(wěn)定,推測主要以J2、J3結(jié)構(gòu)面失穩(wěn),E、F 2處點狀危巖現(xiàn)狀欠穩(wěn)定,其下部巖層差異風(fēng)化嚴重,推測主要以J1結(jié)構(gòu)面失穩(wěn)產(chǎn)生錯斷式崩塌。
危巖體失穩(wěn)后,按所發(fā)生崩塌體的運動過程,可分為滑動、滾動、躍動、自由落體或復(fù)合等多種形式[26],主要與危巖體的密度、磨圓、塊徑、高差及坡形、坡度、坡面粗糙度等有關(guān)。同一物源區(qū)發(fā)生2次崩塌,其運動路徑、能量等也絕非一致,故需要借助數(shù)值手段進行模擬。
采用RocPro3D進行崩塌數(shù)值模擬時主要得確定2類參數(shù):危巖體特征參數(shù)和巖土體介質(zhì)表面參數(shù)。其中危巖體參數(shù)主要有shape、d、d/h、density、v和H,shape為危巖體形狀,d為柱狀危巖體截面直徑,d/h為柱狀危巖體直徑與柱狀危巖體母線之比,density為危巖體密度,v為危巖體初始速度,H為危巖體距離坡體表面的初始高度。巖土體介質(zhì)表面參數(shù)主要有Rn、Rt、k和β,Rn為法向恢復(fù)系數(shù),Rt為切向恢復(fù)系數(shù),k為摩擦系數(shù),β為滾動與躍動間的轉(zhuǎn)化角。
據(jù)不同性質(zhì)對巖土體進行類型劃分(圖10),A~F表示危巖區(qū),①~⑥表示原始基巖面、開挖斜坡風(fēng)化基巖面、生活區(qū)、G338國道、黃土路肩、河流水面,不同巖土體區(qū)域其參數(shù)取值不同。據(jù)表2和DOM測量,給出各危巖區(qū)的危巖體特征參數(shù)取值見表3。
表3 危巖體特征參數(shù)選取
首先根據(jù)Rockfor3D軟件指導(dǎo)[27]及前人經(jīng)驗[28-29]選取巖土體介質(zhì)參數(shù)的取值范圍,然后在取值范圍內(nèi)試算,獲取崩塌模擬的空間分布特征,最后對比模擬結(jié)果和歷史崩塌分布特征,確定危巖體特征參數(shù)和巖土體介質(zhì)表面參數(shù)的合理取值。最終確定參數(shù)見表4,模擬得到危巖區(qū)崩塌的空間分布如圖11所示。
表4 巖土體介質(zhì)表面參數(shù)設(shè)置
圖10 巖土體介質(zhì)及危巖區(qū)劃分Fig.10 Geotechnical media and dangerous rock area division
根據(jù)反演結(jié)果,確定模擬次數(shù)、輸入的危巖體體征參數(shù)及巖土體介質(zhì)表面參數(shù),對A~F危巖區(qū)分別進行計算,分析每處危巖區(qū)的崩塌運動特征。
獲得危巖區(qū)落石到達概率如圖12所示,據(jù)各危巖區(qū)崩塌落石到達概率大小,作出各危巖區(qū)崩塌落石最優(yōu)運動路徑(圖13),并給出最優(yōu)路徑下,各危巖區(qū)落石與坡面間的關(guān)系(圖14),從而量化危巖區(qū)發(fā)生崩塌落石對其下方各區(qū)域的影響,統(tǒng)計結(jié)果見表5。
從落石空間分布來說,A區(qū)崩塌主要從坡體正前方通過,大概率對生活區(qū)形成威脅,最遠也僅到達生活區(qū)內(nèi),不會對停車場、G338國道、黃土路肩和河流水面構(gòu)成威脅;B區(qū)崩塌主要從坡體正前方通過,C、D區(qū)崩塌主要從坡體正前方和側(cè)溝谷通過,B、C和D危巖區(qū)大概率對停車場和少部分生活區(qū)形成威脅,最遠可到達靠近G338國道的邊緣處,不會對G338國道、黃土路肩和河流水面構(gòu)成威脅;E、F區(qū)大概率對停車場形成威脅,最遠也僅到達停車場內(nèi),不會對生活區(qū)、G338國道、黃土路肩和河流水面構(gòu)成威脅。
圖11 危巖區(qū)崩塌空間分布Fig.11 Spatial distribution of collapse in dangerous rock area
從運動規(guī)律來說,落石從危巖區(qū)向下運動過程中,隨著落石與坡面的碰撞,損失的勢能轉(zhuǎn)化成內(nèi)能、熱、聲能等。在坡降大的位置,落石距離坡面的豎向距離也大,碰撞后彈跳的高度往往也較大。落石距離坡面的豎向最大距離多位于坡體的中下部,其勢能逐漸減小,動能先增后減。從坡降較小到坡降較大的位置,落石與坡面間的豎向距離持續(xù)增加,其動能通常也隨著下降運動,在發(fā)生下一次碰撞前的最后一刻,達到峰值。以B區(qū)落石最優(yōu)運動路徑為例,落石平均尺寸2 m×1 m×1 m,通過坡面和停車場時,其落石的質(zhì)心距離坡面的距離6.39 m,距離停車場地面距離6.5 m,對停車場的沖擊動能可達1 600~1 800 kJ,因此會對其停車場構(gòu)成威脅。
據(jù)數(shù)值模擬結(jié)合現(xiàn)場調(diào)查的歷史崩塌分布,6處危巖區(qū)物源豐富明確,節(jié)理發(fā)育,巖體風(fēng)化破碎程度較高,屬于高易發(fā)區(qū),將持續(xù)發(fā)育崩塌。從模擬結(jié)果來看,A危巖區(qū)落石最優(yōu)路徑下,其下方的被動防護網(wǎng)可以起到一定的防護作用,B危巖區(qū)落石最優(yōu)路徑下,其下方被動防護網(wǎng)與落石在防護網(wǎng)處的高度基本持平,難以絕對起到防護作用,C、D、E、F等4個危巖區(qū)下部除有少量樹木,無其他防護,有必要進行防治??衫脭?shù)值模擬的結(jié)果為其防治提供依據(jù)。據(jù)圖8—圖10模擬所得崩塌落石影響范圍、堆積位置、最優(yōu)運動路徑,可明確防護設(shè)施的布設(shè)范圍。據(jù)圖14模擬所得最優(yōu)運動路徑下各危巖區(qū)的崩塌落石距離目標區(qū)的高度和動能,可明確防護的高度、角度和強度。
圖13 落石最優(yōu)運動路徑預(yù)測Fig.13 Prediction of optimal movement path of rockfalls
圖14 最優(yōu)路徑下各危巖區(qū)落石與坡面關(guān)系Fig.14 Relationship of rockfall and slope in each dangerous rock area under optimal path
由圖14知,受6個危巖區(qū)影響最大的目標區(qū)主要是生活區(qū)和停車場,G338國道、黃土路肩、河流水面均不受影響,從圖14的A危巖區(qū)可以看出,在落石最優(yōu)運動路徑下,被動防護網(wǎng)能較為有效地攔截崩塌落石;由圖14的B危巖區(qū)可以看出,在落石最優(yōu)運動路徑下,被動防護網(wǎng)攔截崩塌落石的幾率大幅降低。由于被動防護網(wǎng)僅在生活區(qū)東側(cè)坡體中下部有所設(shè)置,因此具體的防護措施布設(shè)范圍應(yīng)選取從生活區(qū)與停車場交界開始的整個停車場區(qū)域,并根據(jù)模擬結(jié)果結(jié)合工程投資等,給出本研究區(qū)更為合理的崩塌防治方案。
1)野外調(diào)查結(jié)合無人機航測技術(shù)、點云結(jié)構(gòu)面參數(shù)獲取技術(shù),可以實現(xiàn)對陜煤張家峁煤礦工業(yè)廣場東邊坡崩塌落石的歷史崩塌分布、巖體結(jié)構(gòu)面產(chǎn)狀和物源區(qū)危巖體特征的快速提取,有效克服了傳統(tǒng)野外崩塌調(diào)查耗時長、難度大、精度低等弊端,為廣泛工況下的崩塌地質(zhì)災(zāi)害的識別、預(yù)測、防治提供了思路和技術(shù)方法。
2)基于無人機自動航攝的外業(yè)和內(nèi)業(yè)流程化處理,提取高精度點云和DEM數(shù)據(jù),并以現(xiàn)場調(diào)查所統(tǒng)計的數(shù)據(jù)為依據(jù),利用參數(shù)反演確定崩塌模擬的各項參數(shù),進行了崩塌的三維數(shù)值模擬。結(jié)果表明模擬結(jié)果與現(xiàn)場較為吻合,為運動學(xué)模擬提供了可靠依據(jù)。
3)利用三維數(shù)值模擬所獲得的崩塌落石在單位柵格上的到達路徑的概率分布,可預(yù)測落石在坡體上的最優(yōu)運動路徑,有效避免了2D崩塌模擬軟件選取剖面時的主觀干擾因素,顯著減小誤差,有效的獲取了崩塌的影響區(qū)域。結(jié)果表明,6處危巖區(qū)主要威脅工業(yè)廣場東北側(cè)的生活區(qū)和東南側(cè)的停車場,不會對工業(yè)廣場中部的G338國道、工業(yè)廣場西側(cè)的黃土路肩和河流水面構(gòu)成威脅;落石距離坡體的最大高度多集中在坡體中下部,其運動的動能先增大,后減?。晃挥谏顓^(qū)東側(cè),坡體中下部的被動防護網(wǎng)能較為有效地防護A危巖區(qū)的落石,但對B危巖區(qū)的落石防護作用不顯著,C、D、E、F等4個危巖區(qū)主要威脅停車場無防護措施。以上,可為崩塌的防治設(shè)計提供一定的指導(dǎo)依據(jù)。