岩溶区土壤—表层岩溶带具有复杂的三维空间结构, 在垂直方向上, 大气降水以地表径流、竖井水流、壤中流、表层岩溶带水流、深部溶隙、裂隙水流等形式进入地下河；水平方向上表现为裂隙、管道、洞穴的复杂组合^[1-3]。因此, 相比非岩溶区, 岩溶水系统中“五水”(大气降水、地表水、土壤水、表层岩溶带水、地下水)转化过程更加错综复杂, 难于研究^[4-6]。在岩溶区, 地表有土壤层覆盖和岩石出露, 具有非均匀性；岩溶地下河系统是由岩溶管道、岩溶孔隙、岩溶裂隙和岩溶洞穴等多种岩溶空隙介质体组成的多重复合体系, 具有高度非均质性, 很难查清岩溶裂隙、管道等的分布和组成, 地下河水流常呈现出达西流与非达西流并存的现象, 很难准确描述岩溶水的运动^[7-9]。

岩溶区坡地降水-产流特征及其影响因素对于“五水”转化规律的认识以及水土流失防治的研究具有重要意义, 坡地植被可以通过林冠截留、根系提高入渗来减少地表径流量与流速, 拦截、过滤、沉积地表径流的土壤颗粒, 防止水土流失^[10-12]。坡地产流过程是降水与地表诸多因素(植被覆盖率、土壤前期含水量、基岩裸露率等)共同作用的结果, 形成机理较为复杂^[13-14]。岩溶区的坡地产流主要包括地表径流、界面横向流(壤中流、裂隙流)和深层渗漏^[15-16]。姜光辉^[17]将岩溶山区的产流分为5种产流模式：大气-岩石界面的超渗产流；大气-土壤界面的超渗产流或者饱和产流；土壤-岩石界面的壤中流；表岩溶带-包气带界面的表层岩溶带产流；包气带-饱水带界面的地下径流。前人等的研究表明^[18-19], 径流小区坡地产流、土壤剖面水分对降雨的响应受前期土壤含水量、降水量、降雨强度的影响外, 还与土壤所处的地形地貌有关。但岩溶区坡地产流特征以及产流机制研究的比较少, 且对岩溶槽谷区坡地产流特征及其影响因素的研究报道尚少。

本研究以重庆主城北部中梁山背斜顶部龙凤槽谷和龙车槽谷为研究基地, 根据不同土地利用方式设立4个标准径流小区, 在线监测降水量、土壤含水量和地表径流、壤中流及裂隙流产流量, 对比不同土地利用方式下径流小区产流特征, 运用回归分析研究降水因子和土壤含水率对岩溶槽谷区坡地产流的影响, 同时讨论隧道工程对坡地产流的影响, 为岩溶区确定合理的土地利用方式、地下水资源评价与利用提供理论基础。

1 研究地区和方法 1.1 研究区概况

龙凤、龙车槽谷位于重庆市北碚区境内观音峡背斜中梁山岩溶槽谷区(图 1), 该槽谷区属于川东平行岭谷华蓥山帚状褶皱束南沿部分, 海拔高程480—640m。背斜轴向为NNE走向, 核部的碳酸盐岩被两侧的非可溶岩地层包围, 核部为下三叠统飞仙关组(T_lf)地层, 向两翼地层逐渐过渡为下三叠统嘉陵江组(T₁j)、中三叠统雷口坡组(T₂l)及上三叠统须家河组(T₃xj)。背斜核部发育有东、西两个NNE走向的溶蚀槽谷, 形成“一山三岭两槽”的“笔架”式地形。观音峡背斜北、南部分别被嘉陵江、长江切割, 嘉陵江以南海拔逐渐升高, 至凤凰村附近成为地表分水岭, 形成一个相对独立的水文地质单元即龙凤槽谷区, 流域面积约11.7km²。龙凤槽谷南部, 以凤凰村与青坪村之间的山脊为分水岭, 从龙车寺至中梁镇上形成一个相对独立的水文地质单元, 即龙车槽谷, 面积约26.8km²。3条大型隧道位于龙凤槽谷区下方, 沿西北-东南方向横穿中梁山底部, 分别为G75渝武高速北碚隧道、G5001绕城高速施家梁隧道及轨道交通6号线中梁山隧道, 连通重庆主城与北部区县。隧道贯通后, 两侧隧道口出现大量岩溶地下水排泄, 可能影响上部岩溶槽谷区水文流场。

1.2 研究方法 1.2.1 坡面径流小区设置

中梁山岩溶槽谷区, 土地的主要类型分别是林地、耕地及建设地, 三种土地利用类型面积所占中梁山岩溶槽谷区土地面积比例分别是58.1%, 35.8%, 6.1%。根据岩溶区地质地貌特点、生态环境等方面的差异性, 结合当地土地利用方式的实际情况, 在中梁山岩溶槽谷区选择4种典型的土地利用方式(耕地、林地、果园地、竹林地)构成4个标准径流小区(20m×5m), 分别为位于龙凤槽谷(受隧道建设影响区)的耕地、林地、果园地, 以及位于龙车槽谷(非隧道影响区)的竹林地, 其中林地、竹林地保持自然状态(如表 1所示)。径流小区的修建方法是采用四周开挖深沟到达基岩后浇筑钢筋混凝土隔水墙, 墙体高出地面约20cm, 保证了土壤-表层岩溶带原有结构不受干扰, 同时保证了径流小区内部与外界不发生水分交换。根据岩溶区产流特点, 除了对地表径流进行了观测外, 还对壤中流、裂隙流进行了观测, 在径流小区围墙下端设有三个集流池, 每个集流池安装有不锈钢三角堰。在径流小区下方断面上土表的地表径流直接流入对应集流池, 断面土壤-表层岩溶带、表层岩溶带-基岩界面处分别接入导水管, 收集壤中流以及裂隙流, 导水管通过混凝土挡墙过水孔将壤中流和裂隙流引入对应集流池。

1.2.2 数据采集

降水量采用自动气象监测站(DAVIS Vantage Pro2, 天津气象仪器厂有限公司, 中国)实时监测, 精度为±0.2mm, 测量间隔为15min。由于研究径流小区相对集中, 降水量变化不大, 故自动气象站安装于具有代表性的凤凰村村委会楼顶。

土壤含水率采用数据采集器(RR-1016型, 北京雨根科技有限公司, 中国)进行实时监测, 土壤水分的传感器(AV-EC5)分别置于0—20cm和20—40cm土壤深度处, 精度为±0.1%, 量程为0—100%, 测量间隔为15min。

径流小区集流池水位采用HOBO水位计(北京智阳科技有限公司, 中国)进行实时监测, 测量间隔为15min。

1.2.3 数据分析方法

采用Excel、SPSS、Origin统计分析软件, 利用2017年6月1日至2017年10月31日径流小区的降水、土壤含水率、径流等监测数据, 对重庆市中梁山典型岩溶槽谷区不同土地利用方式的产流特征进行了分析研究。

采用线性方程、多项式、指数函数、对数函数、幂函数等五种函数对各径流小区地表径流与降水参数(降水量、降水强度、降水历时)进行了回归分析, 根据拟合参数的显著性来确定最优拟合方程。采用多元线性回归方法分析了各径流小区的地表径流/径流系数与降水参数及土壤含水率之间的相关关系。

2 研究结果 2.1 降水特征

本次监测从2017年6月1日开始到2017年10月31日截止, 如图 2所示, 该时段总降水量为752.4mm, 占年总降水量的70%以上。观测期间共有34次降水过程(根据实测产流特点, 将降水过程中时间间隔超过6h定义为两次降水过程), 观测到22次有效降水过程, 即每个径流小区均有产流发生。4个径流小区有产流过程发生时次最小降水量为14.6mm, 按照气象部门对降水等级的划分：小雨日降雨量为0.1—9.9mm, 中雨日降水量为10—24.9mm, 大雨日降水量为25—49.9mm, 暴雨日降水量≥50mm, 整个研究时间段有产流过程发生时, 共有中雨13次, 大雨7次, 暴雨2次。

2.2 前期土壤含水率

观测期间前期土壤含水率(降水前3h的土壤含水率)如图 3所示, 耕地在22.0%—47.6%之间, 林地在22.8%—38.6%之间, 果园地在18.9%—41.5%之间, 竹林地在10.9%—37.8%之间。平均前期土壤含水率, 耕地为37.6%, 林地为30.1%, 果园地为34.9%, 竹林地为27.1%。土壤含水率和降水特征密切相关, 降水-产流过程发生之前一段时间内无降水或降水少, 则前期土壤含水率低, 降水过程中土壤达到饱和的时间长；降水-产流过程发生之前一段时间内有降水则前期土壤含水率相对较高, 土壤更容易达到饱和。土壤理化性质的决定了不同土地利用方式径流小区产生地表径流时土壤含水率的差异, 地表径流产生时平均土壤含水率耕地＞果园地＞竹林地≈林地(图 3), 这和前期土壤含水率的大小顺序(耕地＞果园地＞林地＞竹林地)基本一致。

2.3 坡地产流特征

从图 4、图 5可以看出, 4个径流小区的坡地径流由地表径流、壤中流和裂隙流组成, 地表径流是主要产流形式。在降水过程中, 随着降水过程的持续, 先后产生了地表径流、壤中流和裂隙流(图 4)。22次产流过程4个径流小区的地表径流、壤中流、裂隙流所占比例如图 5所示, 地表径流占总径流量的比例耕地为90.2%；林地为80.6%；果园地为81.2%；竹林地为73.2%。4个径流小区壤中流和裂隙流在总径流中占比均较小, 耕地壤中流占6.7%, 裂隙流占3.1%；林地壤中流占9.9%, 裂隙流占9.5%；果园地壤中流占7.2%, 裂隙流占11.5%；竹林地壤中流占17.8%, 裂隙流占9.0%。

不同土地利用方式由于受不同程度的人类活动影响, 导致植被、土壤状况等的不同, 其对地表径流、壤中流、裂隙流的径流组成的影响也有差异。

总径流量与径流系数：在相同的降水条件下, 4个不同土地利用方式的径流小区总径流量耕地最大, 为3696.9L, 其次为果园地, 产流总量为3657.2L, 再次为竹林地, 产流总量为2922.9L, 产流总量最少为林地, 总产流量为2211.1L。总径流量耕地＞果园地＞竹林地＞林地(图 4)。

4个径流小区中, 总径流系数在3.1%—5.2%之间, 和非岩溶区相比, 降水量转化成径流量的比例较低, 非岩溶区径流系数大于20%^[20-21]。径流系数耕地为5.2%；果园地为5.2%；竹林地为4.1%；林地为3.1%。在岩溶槽谷区, 到达地面的降水除转化为地表径流、壤中流、裂隙流外, 绝大部分降水由土壤层入渗后经过岩溶裂隙、岩溶管道等最终汇入地下河网；在非岩溶区, 到达地面的降水几乎全部转化为坡面径流和壤中流, 因此, 非岩溶区的总径流系数大于岩溶区。

地表径流与地表径流系数：由图 4可知, 在相同的降水条件下, 4个不同土地利用方式的径流小区地表径流表现为耕地＞果园地＞竹林地＞林地, 和总径流量一致。从6月到10月, 耕地的地表径流最大, 为3334.2L, 其次为果园地, 地表径流为2970.5L, 再次为竹林地, 地表径流为2138.8L, 地表径流最少为林地1781.8L, 径流系数耕地为4.7%；果园地为4.2%；竹林地为3.0%；林地为2.5%。地表径流系数最大为4.7%, 远低于非岩溶区^[21-22]。

壤中流：由界面产流理论, 包气带在垂直方向上的水力传导度的差异为侧向壤中流的产生提供了条件^[23], 在一定的降水条件下, 岩溶槽谷坡地除了产生地表径流外, 在土壤中水力传导度出现差异的地方, 将会有壤中流的产生。该研究表明, 相同的降水条件下, 4个径流小区壤中流也有所差异：耕地为248.9L；林地为219.7L；果园地为264.6L；竹林地为521.1L。壤中流表现为竹林地＞果园地＞耕地＞林地, 径流系数竹林地为0.7%, 果园地为0.4%, 耕地为0.4%, 林地为0.3%。在22场降水过程中, 4个径流小区观测到的壤中流产流次数如下：耕地9次；林地6次；果园地8次；竹林地11次。

不同土地利用方式壤中流的大小差异主要与土壤层的厚度、植物覆盖率差异有关^[24]。竹林地根系较浅, 根系的阻挡作用可有效增加土壤入渗, 从而使得竹林地中壤中流的产流量大于果园地、耕地及林地^[25]。壤中流的产生较地表径流的出现有一定的滞后时间, 壤中流的出现比地表径流的时间晚60min以上, 每场降水过程中壤中流出现的滞后时间和降水量、降水强度以及土层厚度等都有密切关系。

裂隙流：4个径流小区中, 相同的降水条件下, 裂隙流果园地(422.1L)＞竹林地(262.9L)＞林地(209.6L)＞耕地(113.7L), 裂隙流径流系数果园地为0.6%；竹林地为0.4%；林地为0.3%；耕地为0.2%。

2.4 降水特征和前期土壤含水率与产流关系

采用线性方程、多项式、指数函数、对数函数、幂函数等五种函数对各径流小区地表径流与降水参数包括降水历时(T)、降水量(P)、15min最大雨强(I₁₅)、平均降水强度(I)之间的关系进行了单因子回归分析, 根据拟合参数的显著性原则, 多项式函数能更好的反映降水参数与产流量的回归关系。由三项式回归方程拟合结果(图 6), R²均大于0.567, P＜0.01, 拟合效果好。耕地、果园地的地表径流与15min最大雨强呈显著的相关关系, 林地、竹林地的地表径流与降水量存在显著相关关系。

坡面径流除了受到降水量、降水强度、降水历时等降水因子的影响, 还会受到土壤含水率的影响, 径流系数综合反映了一个地区地质土壤等地表状况对径流的影响^[26], 为了消除降水这一对径流量影响最大的因子, 采用地表径流系数数据与降水前期土壤含水率建立相关关系, 如图 7所示。通过线性方程、多项式、指数函数、对数函数、幂函数等五种函数对各径流小区坡面径流系数与前期土壤含水率关系进行了单因子回归分析, 根据拟合参数的显著性原则, 多项式函数能更好的反映坡面径流系数数据与前期土壤含水率的回归关系。由三项式回归方程拟合结果, 耕地、林地、果园地R²均大于0.470, 拟合效果比较好, 耕地、林地、果园地的地表径流系数与前期土壤含水率呈现比较显著相关关系, 竹林地与前期土壤含水率相关关系不明显。

为了研究降水因子、前期土壤含水率共同对地表径流的影响, 采用多元线性回归, 建立各径流小区的地表径流与降水历时(T)、降水量(P)、15min最大雨强(I₁₅)、平均降水强度(I)、前期土壤含水率(Sm)之间的回归方程, 发现降水因子对坡面产流影响大于前期土壤含水率。

将不同土地利用方式下径流小区的回归方程模型系数标准化：

耕地：R₁=-235.225+0.445P+0.671I₁₅     (R²=0.850)

林地：R₂=-112.439+0.720P+0.414I₁₅     (R²=0.836)

果园地：R₃=-177.466+0.484P+0.533I₁₅    (R²=0.656)

竹林地：R₄=-175.332+0.852P        (R²=0.725)

耕地、林地、果园地、竹林地回归方程中R²大于0.656, 说明拟合结果好。通过对上述模型进行了F检验, t检验(表 2)和残差的正态性检验(图略), F检验以及t检验的显著性值均小于0.05, 残差直方图表明模型学生化残差基本成标准正太分布。说明模型具有可靠性, 各自变量与因变量之间存在显著的线性回归关系。

通过以上模型发现, 不同土地利用方式径流小区的地表径流主要受降水因子的影响^[27], 与降水因子都呈显著相关关系。系数标准化的模型中, 各指标的系数大小反应该指标的贡献率, 由模型可知, 耕地、果园地15min最大雨强的贡献率最大, 对地表径流的影响也最显著；林地、竹林地降水量的贡献率最大, 是地表径流的主要影响因素。在耕地、林地、果园地及竹林地中, 前期土壤含水率未通过显著性检验, 说明前期土壤含水率对地表径流的影响小于降水因子。

从以上回归模型可以得出, 降水因子对地表径流的影响远大于前期土壤含水率。

3 讨论

不同土地利用方式的径流小区, 土壤状况、植被覆盖等的差异会导致产流的差异^[28-29]。在中梁山岩溶槽谷区, 总径流量为耕地＞果园地＞竹林地＞林地, 地表径流的大小也表现为耕地＞果园地＞竹林地＞林地。4个径流小区的土壤均为嘉陵江组灰岩发育而来的石灰土, 0—20cm和20—40cm土层深度的pH(7.44—7.88)、容重(0—20cm土层深度为1.05—1.13, 20—40cm土层深度为1.11—1.35)及平均土层厚度(40—45cm)等性质比较一致。但4个径流小区土地利用方式不同, 植被覆盖率有所差异, 植被覆盖率低的地区易产生地表径流^[30], 有研究表明：林地植被可对地表径流造成高达89.3%的消减率^[31]。在非隧道影响区不同的土地利用方式的径流小区中, 耕地的植被覆盖率最低, 且土壤扰动频繁, 径流量最高；果园地植被覆盖率比较低, 果园地中种有柑橘树, 缺乏灌草层次结构, 凋落物少地表覆盖层很薄, 地表土壤层紧密且裸露多, 地表径流相对较大；林地植被覆盖率最大, 由于乔木和灌木交错镶嵌分布, 更有效地利用了不同空间的光照、水分和养分, 构成更稳定的复合缓冲层和地表覆盖层, 对降水截留、径流吸持及延缓径流汇聚起决定作用^[32], 因此, 林地的地表径流最小。综上, 在相同降水条件下, 不同的土地利用方式, 由于植被覆盖率不同, 导致地表径流耕地＞果园地＞林地。

在中梁山岩溶槽谷区, 由于隧道建设, 破坏了地下水循环系统, 加大了深层溶隙水的排泄量, 造成了地下水水源枯竭^[33-35], 如G75渝武高速北碚隧道的贯通使隧道年平均渗漏水量为283.9万m³, 造成水源枯竭面积7.2km², 在影响区内的23个泉(井)中, 17个泉(井)完全枯竭, 4个泉(井)水量减少明显仅在夏季有少量泉水涌出, 仅两个涌水量很少的浅层泉基本没受到影响^[36]。在前期的走访调研中也发现, 三条隧道开挖以来, 槽谷区内原有满足生活、农业用水等的泉(井)水几乎完全枯竭, 目前当地村民用水靠区外抽水到蓄水池再定时集中供水来解决；隧道开挖造成原有溶蚀洼地水田失去水源补充, 部分弃耕或被迫改为旱地。对于隧道建设对坡地径流的影响, 对比林地(隧道影响区)和竹林地(非隧道影响区)两个径流小区, 二者地质背景一致, 土壤类型相同, 植被状况相近, 在研究期间, 竹林地总产流量＞林地总产流量, 地表径流也表现出同样的规律。竹林地总径流系数为4.1%, 林地总径流系数为3.1%, 隧道影响区的坡地总径流比非隧道影响区的坡地总径流减少约24.4%, 说明由于隧道的建设, 造成了坡地岩溶水的漏失, 改变了原有岩溶水的流场特征, 进一步加剧了岩溶槽谷区缺水状况。

学者在研究岩溶区坡地产流规律上大多采用模拟降水方式进行^[37-39], 且对坡地产流的研究多集中在地表径流上^[40-41], 该研究选取了位于重庆市中梁山龙凤(隧道影响区)和龙车(非隧道影响区)槽谷不同土地利用方式(耕地、林地、果园地、竹林地)的4个径流小区, 在自然降水过程中, 除了对降水量、地表径流进行监测外, 同时也对壤中流、裂隙流和土壤含水率进行了同步监测, 数据采集过程中, 每15min输出一次数据, 精度更高。该研究利用较高精度的实时监测数据, 对4种不同土地利用方式径流小区产流特征和影响因素进行对比分析, 更贴近岩溶槽谷区不同土地利用方式下坡地自然产流规律, 能更准确的反映出隧道工程对岩溶槽谷区水文流场的影响作用, 这对岩溶区土地管理、地下水资源的合理利用具有理论意义。该研究仅分析了岩溶槽谷区坡面产流特征及影响因素, 该研究区的具体产流过程和产流机制有待进一步研究。

4 结论

(1) 4个不同土地利用方式的径流小区坡地总产流量耕地＞果园地＞竹林地＞林地, 径流系数远低于非岩溶区。各径流小区地表径流大小和总径流量一致, 壤中流和裂隙流产生滞后于地表径流。

(2) 降水因子、前期土壤含水率共同影响地表径流, 但降水因子对地表径流的影响远大于前期土壤含水率。15min最大雨强是影响耕地、果园地的地表径流的主要因素, 降水量是影响林地、竹林地的地表径流的主要因素。

(3) 前期土壤含水率对耕地、林地、果园地地表径流影响较大, 对竹林地地表径流影响较小。