呂擎峰,王慶棟，郭 軒,王生新，馬 博，姜璐莎

(1.蘭州大學(xué)西部災(zāi)害與環(huán)境力學(xué)教育部重點(diǎn)實(shí)驗(yàn)室，甘肅 蘭州 730000；2.南京工業(yè)大學(xué)數(shù)理科學(xué)學(xué)院，江蘇 南京 211816；3.甘肅省科學(xué)院地質(zhì)自然災(zāi)害防治研究所,甘肅 蘭州 730000)

分形的概念是美籍?dāng)?shù)學(xué)家Mandelbrort首先提出的，1975年，他創(chuàng)立了分形幾何學(xué)。在此基礎(chǔ)上形成了研究分形性質(zhì)及應(yīng)用的學(xué)科，稱為分形理論。分形具有的特征是自相似性、自仿射性、無(wú)標(biāo)度性[1-2]。

分形理論逐漸成為解決復(fù)雜問(wèn)題的有力工具，廣泛應(yīng)用于水系的研究中。早期何隆華等提出了水系分維的含義，對(duì)河長(zhǎng)和河網(wǎng)結(jié)構(gòu)的分維進(jìn)行了探討指出水系的分維反映了地貌侵蝕發(fā)育的階段[3-4]；雷會(huì)珠等[5]通過(guò)分形理論，統(tǒng)計(jì)分析了流域邊界周長(zhǎng)、長(zhǎng)短軸比、匯合角等地貌因子與流域面積的關(guān)系；張少文等[6]較系統(tǒng)的論述了分形理論在水系河網(wǎng)結(jié)構(gòu)和流域地形地貌及其演變、河床表面形態(tài)、徑流過(guò)程等水文水資源中的應(yīng)用；后期國(guó)內(nèi)外學(xué)者通過(guò)試驗(yàn)，依據(jù)分形理論，對(duì)流域地貌形態(tài)特征與侵蝕關(guān)系進(jìn)行了研究[7-9]；李精忠等[10]通過(guò)在匯水量的無(wú)標(biāo)度區(qū)間內(nèi)通過(guò)谷底線長(zhǎng)度與匯水量閾值的線性擬合，計(jì)算谷地的長(zhǎng)度分維；隨著GIS模塊功能的不斷完善，將GIS引入水系研究中成為人們關(guān)注的熱點(diǎn)，范林峰等[11]基于GIS技術(shù)，提出了水系三維分形計(jì)盒維數(shù)的計(jì)算方法，并驗(yàn)證了該方法的可靠性。分形理論應(yīng)用于水系的研究已經(jīng)較多，而結(jié)合已有油氣管道研究沖溝侵蝕的相對(duì)較少。

西部輸氣管道常常穿越山前大型沖積扇，沖溝侵蝕十分發(fā)育，占所有地質(zhì)災(zāi)害的40.9%，沖溝形態(tài)復(fù)雜，呈網(wǎng)狀、樹(shù)枝狀[12-13]。沖溝侵蝕與坡面地貌有著密切的聯(lián)系，哪些指標(biāo)對(duì)沖溝侵蝕的影響起主要作用或次要作用，用傳統(tǒng)的定性方法很難找到答案。由于沖溝是由主溝及支溝組成的網(wǎng)絡(luò)系統(tǒng)，具有典型的自相似性，因此可以進(jìn)行分形分析。本文利用ArcGIS軟件，結(jié)合對(duì)研究區(qū)輸氣管道災(zāi)害點(diǎn)的現(xiàn)場(chǎng)調(diào)查，確定水系的閾值，獲得沖積扇溝網(wǎng)結(jié)構(gòu)，提取坡度及沖積扇的匯水流域。通過(guò)累計(jì)和變換分形理論，研究了沖溝侵蝕與坡面地貌因子的關(guān)系。

1 研究區(qū)概況

研究區(qū)位于青海省烏蘭縣茶卡鎮(zhèn)西南約20 km處的沖積扇，地處98°49′～99°09′E，36°28′～36°42′N，東西寬約31 km，面積約495 km2。該扇地形起伏變化較大，海拔3 070～4 840 m，年均降雨量17.6～764 mm，為典型的高原大陸性氣候，以四季不分明、低溫多變、干燥多風(fēng)、蒸發(fā)量大、高寒缺氧、日照長(zhǎng)、輻射強(qiáng)為其特點(diǎn)。

按照地層和空間展布特點(diǎn)，該扇位于柴達(dá)木北緣分區(qū)。地層主要有震旦系全吉群、長(zhǎng)城系、薊縣系，主要巖性為花崗巖、石英砂礫巖、白云巖、片麻巖、頁(yè)巖等。地處殘山斷褶帶間，構(gòu)造作用強(qiáng)烈。沖積扇表面主要為第四紀(jì)間歇性河相沉積，由扇頂向扇緣顆粒逐漸變小，扇頂堆積主要為卵、礫，分選性差，隨著水流的搬運(yùn)能力的減弱，扇前緣段主要為細(xì)砂、粉砂、黏土。地表被水流沖刷，有些溝道貫穿整個(gè)沖積扇，這樣的溝道侵蝕程度較大；有些溝道則由坡面局部匯流形成，其侵蝕程度相對(duì)較小。

西氣東輸管道西一線、西二線穿越該沖積扇，穿越長(zhǎng)度約27 km，兩線相距18 m。沖溝對(duì)管道進(jìn)行不同程度的破壞，本文對(duì)西一線47個(gè)災(zāi)害點(diǎn)進(jìn)行了調(diào)查。

2 地貌因子

對(duì)沖積扇輸氣管道沖溝侵蝕進(jìn)行調(diào)查，結(jié)合已有的研究[14-17]，利用ArcGIS軟件的分析功能，本文主要從坡度(p)、溝長(zhǎng)(L)、匯水域面積(s)這3項(xiàng)指標(biāo)研究沖溝侵蝕與地貌的分形關(guān)系。

2.1 坡度提取

利用GIS軟件提取坡度，需要以米為單位的投影坐標(biāo)系，因此所有的柵格數(shù)據(jù)均需坐標(biāo)轉(zhuǎn)換，采用同一坐標(biāo)投影(WGS-1984-UTM-Zone-47N)[18-19]。利用GIS軟件，加載數(shù)據(jù)，進(jìn)行地理配準(zhǔn)，幾何糾正，然后利用坡度提取工具，得到?jīng)_積扇的坡度(p)(圖1)。

2.2 河網(wǎng)溝長(zhǎng)提取

圖1 研究區(qū)坡度Fig.1 Slope of study area

在ArcGIS軟件的水文擴(kuò)展分析模塊下，首先對(duì)DEM數(shù)據(jù)進(jìn)行洼地填平，得到無(wú)洼地的DEM數(shù)據(jù)。然后利用水文分析工具，生成流向柵格。其原理是通過(guò)計(jì)算DEM數(shù)據(jù)與其相鄰的8個(gè)柵格單元間的落差值，然后按照坡度最陡原則確定水流方向。最后計(jì)算匯流累積量[20-22]，當(dāng)匯流累積量達(dá)到某一閾值時(shí)，該柵格就會(huì)出現(xiàn)在河道上，出現(xiàn)在河道的柵格賦值為1，小于特定閾值的賦值為0，這樣就生成了連續(xù)的水流線，從而得到河流網(wǎng)絡(luò)柵格。匯流量的閾值需要結(jié)合實(shí)際管道災(zāi)害調(diào)查資料，進(jìn)行多次嘗試，確定適當(dāng)?shù)拈撝??？紤]到坡面局部匯流，也能對(duì)管道運(yùn)營(yíng)造成影響，因此閾值確定為80。對(duì)確定的河網(wǎng)矢量化，并進(jìn)行平滑處理。最后利用Upstream選項(xiàng)，計(jì)算水流長(zhǎng)度。在本文中，水流沖刷形成溝道，因此水流的長(zhǎng)度代表溝長(zhǎng)(L)(圖2)。

圖2 坡面溝長(zhǎng)Fig.2 Length of surface gully

2.3 匯水域的生成

在GIS中對(duì)于匯水域的劃分，都是基于DEM數(shù)據(jù)計(jì)算的(非指定外部?jī)A瀉點(diǎn))，一般用作小流域的提取。匯水域的意義在于水文的相關(guān)性，流域面作為空間范圍，可以圈定該空間位置的一切地理要素、時(shí)空要素與水系之間的關(guān)系。匯水區(qū)域是由分水嶺分割而成的集水區(qū)，可利用水流方向確定出相互連接并處于同一流域盆地的柵格區(qū)域，通過(guò)Basin工具可以計(jì)算匯水區(qū)域。再通過(guò)柵格轉(zhuǎn)面工具，就可以將匯水域矢量化。圖3為矢量河網(wǎng)數(shù)據(jù)與匯水區(qū)域在同一窗口顯現(xiàn)。

圖3 匯水區(qū)域分布圖Fig.3 Distribution area of watershed

3 沖溝侵蝕分形特征

3.1 累計(jì)和變換分形理論

在自然界中，不規(guī)則現(xiàn)象廣泛存在，但幾乎沒(méi)有完全滿足常維分形的，因此，我國(guó)學(xué)者付昱華[23]提出了“變換形成的分形”和“累計(jì)和系列變換”的方法。任意函數(shù)都可以通過(guò)一系列1階、2階、…、n階的變換，分別建立各階累計(jì)和的變維分形模型，然后篩選出最優(yōu)的效果，作為其分形維數(shù)。累計(jì)和變換分形主要通過(guò)以下四個(gè)步驟：

步驟1：確定原始數(shù)據(jù)(Ni，ri)(i=1，2，…，n)，建立lnNi和lnri的雙對(duì)數(shù)坐標(biāo)系，利用式(1)計(jì)算相鄰兩點(diǎn)間的斜率的相反數(shù)，即為研究對(duì)象的變段分維值Di,i+1。

Di,i+1=ln(Ni/Ni+1)/ln(ri+1/ri)

(1)

步驟2：

{Ni}={N1，N2，N3，…}

(2)

在(2)式的基礎(chǔ)上，按式(3)、(4)構(gòu)造各階累計(jì)和序列：

{S1i}={N1，N1+N2，N1+N2+N3，…}(i=1，2，3，…，n)

(3)

{S2i}={S11，S11+S12，S11+S12+S13，…}(i=1，2，3，…，n)

(4)

式中，S1，S2，S3，…是指1，2，3，…階累計(jì)和序列。

步驟3：建立各階累計(jì)和的變維分形模型，并將各階累計(jì)和數(shù)據(jù)點(diǎn)繪制在雙對(duì)數(shù)坐標(biāo)系下。

表1 沖溝侵蝕與坡度的分段變維維數(shù)計(jì)算結(jié)果Table 1 Computed segmentation variable dimension value between gully erosion and slope

利用式(1)計(jì)算相鄰點(diǎn)間的斜率，取其相反數(shù)即為該階累積和的變段分維數(shù)。

步驟4：確定最優(yōu)模型，即直線擬合最好的一次累計(jì)，并計(jì)算斜率，即為分維值D。

3.2 主溝各地貌指標(biāo)與沖溝侵蝕分形關(guān)系

在變換分形過(guò)程中，同一地貌指標(biāo)進(jìn)行分類或分級(jí)，統(tǒng)計(jì)該狀態(tài)下侵蝕溝道的數(shù)目，對(duì)各狀態(tài)下的數(shù)量按照由大到小進(jìn)行排列，分別編號(hào)r=1，2，3，…，7。表1為沖溝侵蝕與坡度的分段變維維數(shù)計(jì)算結(jié)果，沖溝侵蝕主要發(fā)育在0°～14°坡度區(qū)間，占總數(shù)的61.7%。規(guī)定用N表示溝道數(shù)目(個(gè))，一階和二階累計(jì)和分別用S1、S2表示。圖4～圖6分別為沖溝侵蝕與坡度(p)、溝長(zhǎng)(L)、匯水域面積(s)的變維分形關(guān)系。

圖4 沖溝侵蝕與坡度的變維分形關(guān)系Fig.4 The relationships of variable dimension fractal between gully erosion and slope

圖5 沖溝侵蝕與溝長(zhǎng)的變維分形關(guān)系Fig.5 Variable dimension fractal curves of gully erosion vs. gully length

3.3 分維值的意義探討

分形維數(shù)是描述分形結(jié)構(gòu)的定量參數(shù)[24]，反映了沖溝侵蝕分布自相似特征的程度，本質(zhì)上能夠體現(xiàn)沖溝侵蝕與坡度(p)、溝長(zhǎng)(L)、匯水域面積(s)各項(xiàng)地貌指標(biāo)在結(jié)構(gòu)上的復(fù)雜程度。具體表現(xiàn)在：隨著變換階數(shù)越高，各指標(biāo)對(duì)沖溝侵蝕的影響機(jī)制越復(fù)雜；分維值越大，影響程度越強(qiáng)。

圖6 沖溝侵蝕與匯水域面積的變維分形關(guān)系Fig.6 Curves of variable dimension fractal of gully erosion varied with area of watershed

沖溝侵蝕程度與坡度(p)、坡長(zhǎng)(L)呈二階累計(jì)和分形分布，變維分形特征復(fù)雜，分維值也較大，表明這兩項(xiàng)指標(biāo)對(duì)沖溝侵蝕影響較大，是控制沖溝侵蝕的主要因素，坡度(p)對(duì)沖溝侵蝕的控制作用最為復(fù)雜，影響程度高于坡長(zhǎng)(L)；沖溝侵蝕程度與匯水面積(s)呈低階累計(jì)和分形分布，分維值相對(duì)較小，變維分形特征相對(duì)簡(jiǎn)單，說(shuō)明該指標(biāo)對(duì)沖溝侵蝕的控制程度較弱(表2)。

表2 沖溝侵蝕與地貌指標(biāo)分形匯總表Table 2 The fractal summary table of gully erosion and geomorphological index

4 結(jié)論

本文通過(guò)利用GIS的水文分析功能，結(jié)合管道現(xiàn)場(chǎng)調(diào)查，獲取研究區(qū)的坡度(p)、溝長(zhǎng)(L)、匯水域面積(s)3項(xiàng)地貌指標(biāo)，運(yùn)用累計(jì)和變換分形理論，研究了沖積扇輸氣管道沖溝侵蝕與地貌因子的變維分形特征，得出以下結(jié)論：

(1)沖溝侵蝕與坡度(p)、溝長(zhǎng)(L)、匯水域面積(s)3項(xiàng)地貌指標(biāo)均有變維分形的特征，這些指標(biāo)在一定程度上控制著沖溝侵蝕的發(fā)育程度。

(2)沖溝侵蝕與坡度(p)、溝長(zhǎng)(L)呈現(xiàn)二階累計(jì)和分形分布，說(shuō)明這兩項(xiàng)指標(biāo)對(duì)沖溝侵蝕的控制作用最為復(fù)雜，且坡度分維值(Dp=1.35)>溝長(zhǎng)分維值(DL=1.23)，表明坡度對(duì)沖溝侵蝕的影響高于溝長(zhǎng)。

(3)沖溝侵蝕與匯水域面積(s)呈一階累計(jì)和分形分布，分維值(Ds=0.39)較小，表明匯水域面積對(duì)沖溝侵蝕的影響較弱。

(4)本文為研究山前沖積扇輸氣管道沖溝侵蝕與地貌因子的復(fù)雜關(guān)系提供一種新思路?？梢愿鶕?jù)對(duì)沖積扇地貌的初判，對(duì)輸氣管道運(yùn)營(yíng)期的主要治理區(qū)域提供一定的依據(jù)。