山前冲积扇上输气管道冲沟侵蚀分形特征

2018-10-31吕擎峰王庆栋王生新姜璐莎

中国地质灾害与防治学报 2018年5期

吕擎峰,王庆栋，郭轩,王生新，马博，姜璐莎

(1.兰州大学西部灾害与环境力学教育部重点实验室，甘肃兰州 730000；2.南京工业大学数理科学学院，江苏南京 211816；3.甘肃省科学院地质自然灾害防治研究所,甘肃兰州 730000)

分形的概念是美籍数学家Mandelbrort首先提出的，1975年，他创立了分形几何学。在此基础上形成了研究分形性质及应用的学科，称为分形理论。分形具有的特征是自相似性、自仿射性、无标度性[1-2]。

分形理论逐渐成为解决复杂问题的有力工具，广泛应用于水系的研究中。早期何隆华等提出了水系分维的含义，对河长和河网结构的分维进行了探讨指出水系的分维反映了地貌侵蚀发育的阶段[3-4]；雷会珠等[5]通过分形理论，统计分析了流域边界周长、长短轴比、汇合角等地貌因子与流域面积的关系；张少文等[6]较系统的论述了分形理论在水系河网结构和流域地形地貌及其演变、河床表面形态、径流过程等水文水资源中的应用；后期国内外学者通过试验，依据分形理论，对流域地貌形态特征与侵蚀关系进行了研究[7-9]；李精忠等[10]通过在汇水量的无标度区间内通过谷底线长度与汇水量阈值的线性拟合，计算谷地的长度分维；随着GIS模块功能的不断完善，将GIS引入水系研究中成为人们关注的热点，范林峰等[11]基于GIS技术，提出了水系三维分形计盒维数的计算方法，并验证了该方法的可靠性。分形理论应用于水系的研究已经较多，而结合已有油气管道研究冲沟侵蚀的相对较少。

西部输气管道常常穿越山前大型冲积扇，冲沟侵蚀十分发育，占所有地质灾害的40.9%，冲沟形态复杂，呈网状、树枝状[12-13]。冲沟侵蚀与坡面地貌有着密切的联系，哪些指标对冲沟侵蚀的影响起主要作用或次要作用，用传统的定性方法很难找到答案。由于冲沟是由主沟及支沟组成的网络系统，具有典型的自相似性，因此可以进行分形分析。本文利用ArcGIS软件，结合对研究区输气管道灾害点的现场调查，确定水系的阈值，获得冲积扇沟网结构，提取坡度及冲积扇的汇水流域。通过累计和变换分形理论，研究了冲沟侵蚀与坡面地貌因子的关系。

1 研究区概况

研究区位于青海省乌兰县茶卡镇西南约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，为典型的高原大陆性气候，以四季不分明、低温多变、干燥多风、蒸发量大、高寒缺氧、日照长、辐射强为其特点。

按照地层和空间展布特点，该扇位于柴达木北缘分区。地层主要有震旦系全吉群、长城系、蓟县系，主要岩性为花岗岩、石英砂砾岩、白云岩、片麻岩、页岩等。地处残山断褶带间，构造作用强烈。冲积扇表面主要为第四纪间歇性河相沉积，由扇顶向扇缘颗粒逐渐变小，扇顶堆积主要为卵、砾，分选性差，随着水流的搬运能力的减弱，扇前缘段主要为细砂、粉砂、黏土。地表被水流冲刷，有些沟道贯穿整个冲积扇，这样的沟道侵蚀程度较大；有些沟道则由坡面局部汇流形成，其侵蚀程度相对较小。

西气东输管道西一线、西二线穿越该冲积扇，穿越长度约27 km，两线相距18 m。冲沟对管道进行不同程度的破坏，本文对西一线47个灾害点进行了调查。

2 地貌因子

对冲积扇输气管道冲沟侵蚀进行调查，结合已有的研究[14-17]，利用ArcGIS软件的分析功能，本文主要从坡度(p)、沟长(L)、汇水域面积(s)这3项指标研究冲沟侵蚀与地貌的分形关系。

2.1 坡度提取

利用GIS软件提取坡度，需要以米为单位的投影坐标系，因此所有的栅格数据均需坐标转换，采用同一坐标投影(WGS-1984-UTM-Zone-47N)[18-19]。利用GIS软件，加载数据，进行地理配准，几何纠正，然后利用坡度提取工具，得到冲积扇的坡度(p)(图1)。

2.2 河网沟长提取

图1 研究区坡度Fig.1 Slope of study area

在ArcGIS软件的水文扩展分析模块下，首先对DEM数据进行洼地填平，得到无洼地的DEM数据。然后利用水文分析工具，生成流向栅格。其原理是通过计算DEM数据与其相邻的8个栅格单元间的落差值，然后按照坡度最陡原则确定水流方向。最后计算汇流累积量[20-22]，当汇流累积量达到某一阈值时，该栅格就会出现在河道上，出现在河道的栅格赋值为1，小于特定阈值的赋值为0，这样就生成了连续的水流线，从而得到河流网络栅格。汇流量的阈值需要结合实际管道灾害调查资料，进行多次尝试，确定适当的阈值。考虑到坡面局部汇流，也能对管道运营造成影响，因此阈值确定为80。对确定的河网矢量化，并进行平滑处理。最后利用Upstream选项，计算水流长度。在本文中，水流冲刷形成沟道，因此水流的长度代表沟长(L)(图2)。

图2 坡面沟长Fig.2 Length of surface gully

2.3 汇水域的生成

在GIS中对于汇水域的划分，都是基于DEM数据计算的(非指定外部倾泻点)，一般用作小流域的提取。汇水域的意义在于水文的相关性，流域面作为空间范围，可以圈定该空间位置的一切地理要素、时空要素与水系之间的关系。汇水区域是由分水岭分割而成的集水区，可利用水流方向确定出相互连接并处于同一流域盆地的栅格区域，通过Basin工具可以计算汇水区域。再通过栅格转面工具，就可以将汇水域矢量化。图3为矢量河网数据与汇水区域在同一窗口显现。

图3 汇水区域分布图Fig.3 Distribution area of watershed

3 冲沟侵蚀分形特征

3.1 累计和变换分形理论

在自然界中，不规则现象广泛存在，但几乎没有完全满足常维分形的，因此，我国学者付昱华[23]提出了“变换形成的分形”和“累计和系列变换”的方法。任意函数都可以通过一系列1阶、2阶、…、n阶的变换，分别建立各阶累计和的变维分形模型，然后筛选出最优的效果，作为其分形维数。累计和变换分形主要通过以下四个步骤：

步骤1：确定原始数据(Ni，ri)(i=1，2，…，n)，建立lnNi和lnri的双对数坐标系，利用式(1)计算相邻两点间的斜率的相反数，即为研究对象的变段分维值Di,i+1。

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

(1)

步骤2：

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

(2)

在(2)式的基础上，按式(3)、(4)构造各阶累计和序列：

{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，…阶累计和序列。

步骤3：建立各阶累计和的变维分形模型，并将各阶累计和数据点绘制在双对数坐标系下。

表1 冲沟侵蚀与坡度的分段变维维数计算结果Table 1 Computed segmentation variable dimension value between gully erosion and slope

利用式(1)计算相邻点间的斜率，取其相反数即为该阶累积和的变段分维数。

步骤4：确定最优模型，即直线拟合最好的一次累计，并计算斜率，即为分维值D。

3.2 主沟各地貌指标与冲沟侵蚀分形关系

在变换分形过程中，同一地貌指标进行分类或分级，统计该状态下侵蚀沟道的数目，对各状态下的数量按照由大到小进行排列，分别编号r=1，2，3，…，7。表1为冲沟侵蚀与坡度的分段变维维数计算结果，冲沟侵蚀主要发育在0°～14°坡度区间，占总数的61.7%。规定用N表示沟道数目(个)，一阶和二阶累计和分别用S1、S2表示。图4～图6分别为冲沟侵蚀与坡度(p)、沟长(L)、汇水域面积(s)的变维分形关系。

图4 冲沟侵蚀与坡度的变维分形关系Fig.4 The relationships of variable dimension fractal between gully erosion and slope

图5 冲沟侵蚀与沟长的变维分形关系Fig.5 Variable dimension fractal curves of gully erosion vs. gully length

3.3 分维值的意义探讨

分形维数是描述分形结构的定量参数[24]，反映了冲沟侵蚀分布自相似特征的程度，本质上能够体现冲沟侵蚀与坡度(p)、沟长(L)、汇水域面积(s)各项地貌指标在结构上的复杂程度。具体表现在：随着变换阶数越高，各指标对冲沟侵蚀的影响机制越复杂；分维值越大，影响程度越强。

图6 冲沟侵蚀与汇水域面积的变维分形关系Fig.6 Curves of variable dimension fractal of gully erosion varied with area of watershed

冲沟侵蚀程度与坡度(p)、坡长(L)呈二阶累计和分形分布，变维分形特征复杂，分维值也较大，表明这两项指标对冲沟侵蚀影响较大，是控制冲沟侵蚀的主要因素，坡度(p)对冲沟侵蚀的控制作用最为复杂，影响程度高于坡长(L)；冲沟侵蚀程度与汇水面积(s)呈低阶累计和分形分布，分维值相对较小，变维分形特征相对简单，说明该指标对冲沟侵蚀的控制程度较弱(表2)。