傅舜宇，张青，陶宇，鲍依群

（绍兴市水环境科学研究院有限公司，浙江 绍兴 312000）

香榧Torreya grandis‘Merrillii’是我国特有的珍稀干果，主产于浙江会稽山区的绍兴、诸暨、嵊州、东阳、磐安五个县市[1]。根据笔者的实地走访、调查，约从2014 年开始，在绍兴南部、嵊州北部、诸暨西部的丘陵山区，大面积的香榧造林采用先在平地种植香榧苗，待成年后再移栽至山坡的种植方式。由当地乡镇政府提供的统计数据表明，2017 年，绍兴市柯桥区的稽东镇，嵊州市的谷来镇、竹溪乡，香榧的种植面积占这些乡镇土地总面积的20%以上。为了方便讨论，本文将树龄≤3 a，在平地（多为水田或旱地改造）以起垄的方式栽种的香榧林称为香榧苗圃；将树龄3～30 a，被移栽至山坡种植的香榧林称为香榧幼林；将树龄＞30 a 的次生香榧林和香榧自然林称为香榧成林。根据2017 年的遥感影像解析，结合数月的现场实地走访和资料查阅，得知香榧林（包括香榧成林、香榧幼林、香榧苗圃）约占汤浦水库流域面积的10.2%，茶Camellia s inensis园占4.3%，其他农林地类占流域面积的比例均低于3%。

据笔者2016—2017 年对汤浦水库上游流域96 个行政村的实地走访和问卷调查显示，香榧幼林和香榧成林的施肥方式主要有两种：（1）一年两次，一次在2—3 月，一次在9—10 月；（2）一年三次，一次在2 月，一次在4—5 月，一次在9—10 月。每次施肥的氮施入量如下：香榧幼林为2.3～18.0 kg·hm-2、香榧成林为6.8～27.0 kg·hm-2、香榧苗圃为0.8～4.5 kg·hm-2。另有少数农户在1、2、3、4、5、9、10、11 月分别对香榧幼林进行施肥，每次施入氮为2.3～5.3 kg·hm-2，（均折合为纯氮用量计量）。此外，部分农户每年会施1～2 次有机肥，每次施用量为100～500 kg·hm-2。

香榧作为松杉纲Coniferopsida 植物，耐胁迫能力强，能够在陡峭而贫瘠的山坡上生长。在浙北红壤丘陵地区，一般认为不适合农林业开发的坡度25°以上的坡地，甚至坡度大于40°的坡地，被大量开垦为香榧幼林。在谷来镇、竹溪乡、稽东镇，多处坡度60°～70°的极陡坡也被开垦为香榧幼林。有研究表明，在红壤丘陵地区，覆盖草本植物的坡地的地表径流、壤中流、基流造成的氮流失均显著低于清除地表草本植被的裸土坡地[2]。另一方面，已有很多研究表明，过量施氮或不合理施氮将加速土壤酸化、破坏土壤结构[3-4]，在一定条件下可能因淋溶导致大量氮营养盐进入水体[5-9]。

过去，国内外对香榧的研究集中在生长习性[5]、生态习性[6-7]、香榧对土壤生态的影响[1,8-10]这几个领域。2020年，有学者研究了香榧人工种植林在面对极端降雨时，土壤总氮、土壤总磷、土壤总有机碳的损失[11]。目前，国内外还没有运用SWAT 模型或径流小区试验研究香榧对流域的非点源氮污染的影响。本研究将着眼于香榧林，特别是人工种植的香榧幼林对于汤浦水库氮输入的影响，旨在为水库水质保护和非点源污染防治提供依据。

1 研究区概况

研究区为汤浦水库坝址流域（120°28'～120°49' E，29°37'～29°56' N），地处亚热带季风区，在行政区划上分属柯桥区稽东镇、王坛镇、平水镇和嵊州市谷来镇、竹溪乡、王院乡。汤浦水库位于曹娥江支流小舜江的上游，水库水域面积为13.64 km2，流域集水面积约为460 km2；流域年均降水量为1 564.4 mm；年均径流深为795.2 mm。主要入库河道为南溪（长约为43.8 km，流域面积约为215 km2）、北溪（长约为25.8 km，流域面积约为132 km2）、王化溪（流域面积约为61 km2）、万乔溪（流域面积约为14 km2）。汤浦水库自1997 年建库以来，水质总体良好，大部分指标符合GB3838—2002 的I 类水标准。不过，库区水的总氮浓度在2013—2015 年时有所升高，虽然在2016—2017 年时略有下降，但仍高于2011—2012 年时的水平。

研究区地貌以低山丘陵为主，地势从西南向东北倾斜。红壤是绍兴市柯桥区南部和嵊州市北部的丘陵地区的主要土壤类型[12-13]。按照土壤发生分类法，研究区内的土壤以黄泥土、黄泥砂土、黄砾泥等红壤为主[13]。研究区内的香榧幼林、成林的土壤总氮平均含量为1 569.0±481.8 mg·kg-1，有效氮平均含量为154.0±41.0 mg·kg-1。

2016—2018 年，本研究在对汤浦水库流域的调查中发现，南溪的榆树支流上游、竹溪支流上游是香榧幼林的集中种植区。这些地区坡度40°以上的坡地被大面积开垦为香榧幼林。这些香榧幼林的郁闭度普遍小于0.1。当地农民在将香榧从香榧苗圃移栽至山坡时，普遍使用草甘膦等除草剂清除所有林下植被。这些因素导致了多处山体的土壤近似于裸露。

2 研究方法

2.1 SWAT 模型的构建

2.1.1 数据库的构建 前人在汤浦水库流域的研究已证明SWAT 模型适用于该流域[14-17]。以往的研究均采用20世纪90 年代第一次土地普查所制的土地利用类型分布图和土壤类型分布图，均未对SWAT 模型的数据库进行本土化，也均未涉及香榧这一本地特色作物。本研究根据高分一号卫星遥感影像，结合实地走访，手动绘制了全流域土地利用类型图（图1），并对土壤属性、土壤类型等数据库进行了本土化处理。本研究构建SWAT 模型的输入数据来源见表1。

表1 SWAT模型的空间与属性数据Table 1 The spatial and attribute data of SWAT model

图1 S WAT模型的土壤类型分布图（a）和土地利用方式图（b）Table 1 Soil type (a) land-use (b) of SWAT model

研究区内的土地利用类型被分为16 类。其中，稻田、旱地、茶园、自然林等大多数土地利用类型在SWAT模型自带的数据库中有对应的土地利用类型代码。但是，本研究重点关注的香榧成林、香榧幼林、香榧苗圃在SWAT 模型自带数据库中，没有相应的土地利用类型数据。因此，本研究在SWAT模型的土地利用类型数据库中，自建了香榧幼林、香榧成林、香榧苗圃3 种土地利用类型代码。对稻田、旱地、茶园、自然林等其他土地利用类型，本研究采用SWAT 模型自带数据库中对应的土地利用类型数据。

2.1.2 子流域和水文响应单元（HRU）的划分 在数据、信息的收集工作完成后，进行空间离散化操作。汤浦水库上游流域被划分成31个子流域（图2）和1 150个HRU。SWAT模型只允许存在一个总出水口，而汤浦水库库区湖面面积大，不能将湖面视为一个子流域。此外，汤浦水库4条主要入库河道中，万乔溪、王化溪、除占岙支流外的北溪和南溪下游的土地利用类型的分布相对传统。它们和有香榧幼林集中分布的北溪占岙支流、南溪上游榆树、竹溪等支流存在较大差异。因此，本研究先选取香榧集中种植区所在的几个小流域进行模拟，再将模型推广到其他河道流域。

图2 汤浦水库流域的子流域划分Figure 2 The sub streams of Tangpu Reservoir basin

2.1.3 模型的评价 SWAT模型的模拟结果采用决定系数（R2）和Nash-Sutcliffe模拟效率系数（Ens）来评价。R2和Ens的计算方法如下：

式中，n为观测次数，Qm为观测值，QS为模拟值，Qavgm为观测值平均值，Qavgs为模拟值平均值。当R2>0.6、Ens>0.5 时，通常认为模型的模拟效果较好，具有较高的可信度。

2.2 香榧径流小区试验

香榧幼林的产流率、氮输出负荷并无文献可查。另一方面，香榧幼林是本研究在SWAT 模型数据库中自建的土地利用类型，SWAT 模型对这一地类的模拟需要验证。因此，本研究于2017 年10 月—2018 年11 月开设了径流小区试验。香榧幼林集中分布在南溪的榆树支流上游，因此，本研究在榆树支流上游的长蛇水库北岸山坡（29°42' N，120°33' E）的香榧幼林设置径流小区试验场地。小区从坡顶延至坡底，长约为23 m，宽约为5 m，坡度约为40°，投影面积约为100 m2。设两组平行试验。两个小区种植的香榧树树龄为8～10 a，株高为1.2～1.5 m，郁闭度约为0.05，株行距为4 m×4 m，均无林下植被。每组小区的边界用水泥板围成。水泥板露出地面20 cm，埋入地下30 cm，水泥板表面和板之间用防水涂料进行处理。径流小区的2 个集水槽及其连同的集水桶，分别用来收集地表径流出水和0～30 cm 土层土壤中的流出水。另外，在径流小区旁设取样器收集雨水水样。径流小区于2018 年3 月下旬、5 月下旬分别施入氮3.75 kg·hm-2，于9 月上旬施入氮2.25 kg·hm-2（施肥量全部已折算成纯氮用量）。

2.3 土壤调查

本研究在香榧幼林径流小区及研究区内若干代表性地段设置样线。每条样线采用蛇形布点法设16 个采样点。每个采样点之间相距5 m 左右。在香榧幼林的样线，在每个采样点内随机选3 棵香榧，在每棵香榧的树干半径50 cm 距离内取样作为冠内组；将相邻的香榧树沿斜对角相连，选3 处交叉处取样作为冠外组。在其他样线，每个采样点内按照十字交叉线取5 个土样。土样用土壤采样器分0～10 cm 的表层土壤和＞10～20 cm 的深层土壤采集。将每条样线中相同土层的土样混合，用四分法将多余的土壤舍弃；可得2 个混合样，测定所有样线混合土样的土壤全氮含量。

3 结果与分析

3.1 参数率定的结果

以2011 年为缓冲年，2012 年1 月—2015 年12 月作为校准期，2016 年1 月—2018 年12 月作为验证期，以汤浦水库入库水量对模型模拟的流域径流量进行校准。再以2015 年为校准期，以2016—2018 年为验证期，以南溪入库口的断面20 次水质调查的总氮数据作为观测值，将观测值和模型对总氮的模拟值分别进行营养盐率定。再用校准参数后的模型对北溪、王化溪、万乔溪的入库口断面的总氮实测值进行验证。本研究选用的率定参数见表2。

表2 径流率定、营养盐率定的参数及其范围Table 2 Parameter calibration of runoff and eutrophication and their range and optimum value

模型模拟的校准期径流量和同期实测值的R2=0.96、Ens=0.96；模拟的验证期径流值和同期实测值的R2=0.90、Ens=0.89，均符合评判标准R2>0.6、Ens>0.5 的要求（图3）。模型对南溪入库口校准期、验证期的总氮模拟，也均符合R2>0.6、Ens>0.5 的要求。模型对榆树支流等多个香榧种植区、骆村支流等多个混合种植区、北溪入库口、王化溪入库口、万乔溪入库口的总氮输出量的模拟结果，和对应河道断面同期水质调查的实测值的评价指标，亦均符合R2>0.6、Ens>0.5 的要求（图4）。这表明本模型对汤浦水库上游流域总氮输出的模拟效果较好，模型模拟的可信度高。

图3 汤浦水库流域的入库径流Figure 3 Inflow runoff into Tangpu Reservoir

图4 南溪、北溪、王化溪、万乔溪入库口断面的总氮通量Figure 4 Total nitrogen flux at the entrance sections of Nanxi River,Beixi River,Wanghua River and Wanqiao River

3.2 径流小区试验结果

径流小区于2017 年6 月底完成施工，于8 月完成调试。径流小区试验自2017 年10 月开始，至2018 年9月结束。表3 为试验期间15 次取样的日期和降雨情况。由产流结果显示，香榧幼林的产流情况受到雨量等级和降雨历时两重影响。试验期间的两次大雨、暴雨和一次暴雨的后续降雨中，香榧幼林的地表径流产流量可达到13%～15%。同浙北红壤丘陵地区的研究结果相比，香榧幼林在大雨、暴雨期间的产流率与临安青山湖没有林下植被的竹林（14%～21%）[18]在暴雨期间的产流率相近；高于钱江源地区硬阔林（无径流－8%）[19]和临安青山湖具有林下植被的竹林（4%～7%）[18]在大雨、暴雨期间的产流率，高于新安江、西苕溪的茶园（4%～9%）[20-21]在大雨、暴雨期间的产流率。

图5 历次降雨事件中，香榧幼林径流小区的总氮流失量Figure 5 The TN loss of runoff plot in young T.grandis ‘Merrillii’ stand on rainy days

表3 径流小区历次降雨情况Table 3 Rainfall in runoff plot

3.3 径流小区试验对SWAT 模型的验证

基于香榧幼林径流小区的15 次径流出水、同期雨水和当地其他自然林的径流水样等实测数据，本研究采用输出系数模型对南溪榆树支流上游的总氮输出量进行估算。由输出系数模型计算得到的该子流域总氮输出量，被用来验证SWAT 模型对这一香榧集中种植区总氮输出量的估算结果。本文所用的输出系数模型的简化公式源自尹微琴等的研究[22]。

式中，L表示总氮输出量，Ei表示第I 种土地利用输出系数，Ai表示第i种土地利用类型的面积，Ii表示单位面积总氮输出量，P表示降雨输入的总氮。

在榆树支流上游的香榧集中种植区，输出系数模型计算得到的氮输出量和SWAT 模型的同期模拟值的R2=0.75、Ens=0.61，符合评判标准R2>0.6、Ens>0.5（图6）。这一结果进一步说明，本研究的SWAT 模型对香榧集中种植区的总氮输出的模拟效果较好，模型模拟的可信度高。

图6 基于径流小区试验的输出系数模型和SWAT 模型估算的南溪榆树支流上游的氮输出量Figure 6 Nitrogen f rom u pstream o f Y ushu e stimated by S WAT model and output coefficient model based on runoff plot

3.4 SWAT 模型对流域内总氮的模拟

3.4.1 非点源氮污染的空间分布 SWAT 模型的模拟结果表明，2017 年10 月—2018 年9 月，南溪流域能贡献约75%的入库总氮（表4）。王化溪、万乔溪、库区的总氮输出量的贡献比率则相对较低。其中，南溪的榆树支流上游、竹溪支流上游、北溪的占岙支流上游的总氮年输出量可达3.0 t·km-2（图7）。这些区域都是香榧集中种植区，香榧林的面积可占相应子流域面积的50%～88%。

图7 汤浦水库流域2017 年10 月—2018 年9 月总氮输出的空间分布Figure 7 Spatial distribution of total nitrogen output in Tangpu Reservoir basin from Oct.2017 to Sep.2018

表4 SWAT 模型模拟的汤浦水库入库各河道流域非点源氮通量Table 4 Non point source nitrogen flux in Tangpu Reservoir inflows simulated by SWAT model

3.4.2 不同土地利用方式对非点源氮污染的贡献 如表5显示，约有45%的总氮入库量来自香榧幼林，加上香榧成林和香榧苗圃贡献的总氮入库量，三者的总和比例超过入库总氮的一半。香榧幼林占汤浦水库上游约6%的面积。从单位面积氮负荷来看，香榧幼林远远超过其他几种土地利用方式。

表5 SWAT 模型对不同土地利用方式的总氮入库量的模拟Table 5 Simulation of total nitrogen inflows into Tangpu Reservoir from different landuse by SWAT model

3.5 香榧林土壤调查结果

香榧种植区普遍存在过量施化肥的现象[1,10]。在研究区内的实地调查表明，香榧集中种植区的施肥量明显高于附近地区其他土地利用方式的施肥量。本研究于2018 年在研究区的土壤调查表明，香榧种植区0～15 cm 土层的土壤总氮含量可达同一子流域的其它自然林和茶园同层土壤的1.6 倍和1.5 倍；采用“一年两次，每次大量”施肥方式下的香榧林15～40 cm 土层的土壤总氮含量可达同一子流域的其它自然林和茶园同层土壤的2.0 倍和1.5 倍（表6）。本研究的土壤调查和戴文圣等2005 年在绍兴、嵊州等五县市香榧林开展的土壤调查[1]情况相似，都反映了浙江红壤丘陵地区的香榧林因为滥施肥导致土壤氮含量过量。

表6 研究区的土壤氮含量Table 6 Nitrogen content in research area

4 结论与讨论

4.1 结论

本研究利用汤浦水库2011—2018 年的入库水量和2015 年多个河道断面的水质检测结果对SWAT 模型的参数进行率定、校准：在上游的香榧幼林集中种植区——榆树支流上游，采用香榧幼林径流小区试验数据对SWAT模型进行了验证；在南溪、北溪的入库断面，采用2016—2018 年的河道断面水质检测结果进行验证。模拟结果表明，本研究所用的SWAT 模型的模拟精度较高，并且该模型能够反映汤浦水库上游流域包括香榧幼林、香榧成林、香榧苗圃在内的多种土地利用方式的非点源氮输出。根据SWAT 模型的模拟结果，占汤浦水库上游流域面积约46%的南溪小流域贡献了75%的入库总氮。南溪的榆树支流上游、竹溪支流、北溪的占岙支流这些香榧幼林集中种植区上游的非点源氮污染较严重。占流域面积约6%的香榧幼林提供了约45%的入库总氮，是贡献最大的土地利用类型。

4.2 讨论

香榧幼林的氮流失很可能与以下因素有重要关联：（1）过低的郁闭度、近乎裸露的土壤。林下植被能够对降雨进行再分配，影响坡地径流的产生，进而影响总氮的流失[23]。林下植被的存在可以减少大雨、暴雨期间大雨滴对土壤颗粒的溅蚀作用，使得坡面流没有足够的潜力冲走更多的土壤颗粒[24]。榆树支流上游、竹溪支流上游等地通常在陡坡开荒种植香榧幼林，并使用草甘膦等除草剂清除全部林下植被。这种种植方式导致山坡土壤裸露。本研究的径流小区试验表明，这种方式形成的香榧幼林在大雨、暴雨期间，坡面会产生比自然林、茶园更多的径流。（2）过陡的种植坡度。在研究区内，大量香榧幼林种植在40°以上的陡坡上，甚至60°～70°的极陡坡上。陡坡的开垦加剧了水土流失。（3）过量的、不合理的施用氮肥。这使香榧集中种植区的土壤氮本底高于周围其他作物的土壤。在浙江临安山区，类似的砍伐山坡的自然林，清除全部林下植被种植经济林木的生产方式，被证明会产生严重的土壤侵蚀和氮素流失[25]。本研究中，大量香榧幼林的种植坡度（40°左右，部分地区可达60°～70°）高于临安山核桃Carya cathayensis的种植坡度（25°左右），因此，可能将面临更严峻的氮流失风险。

目前，汤浦水库上游流域各乡镇政府已不允许在陡坡开垦香榧幼林。未来，有关部门可向当地农户普及水土保持学知识，使他们了解保留林下植被的益处。此外，应提倡在香榧幼林套种林下经济作物，推广等高种植等减少水土流失的种植方式，提倡测土配方施肥等科学施肥方式。