王文章，程 艳，敖天其,2，黎小东

(1.四川大学 水利水电学院，成都 610065；2.四川大学 水力学与山区河流开发保护国家重点实验室，成都 610065)

0 引 言

随着点源污染得到了有效治理与控制，面源污染已成为国内水环境污染治理研究的热点[1]。SWAT模型以扩展模块的形式嵌入与Arcgis软件中，能够在不同土壤类型和土地利用的大尺度流域内模拟流域产流、产沙以及营养物负荷等的变化情况,已成功应用于国内多个流域的非点源污染模拟研究[2]。

子流域划分是分布式水文模型模拟的最基本步骤，也是保证模拟精度的前提。流域的建模过程和模拟结果都会受到子流域大小、规模和数量的影响，不同的子流域划分的参数敏感性差异也较大。张雪松等[3](2004)认为模型输入参数的空间集总程度会随子流域划分数量大小而改变，从而对流域模拟结果产生一定影响; 喻晓等[4](2017)的研究表面，较低的亚流域划分水平会造成模型结果不稳定；王国波等(2016)[5]探讨了不同子流域数量与流域关键污染源区(CSAs)的相关关系，表明相对于SWAT模型自动确定的集水阈值，运用合理的方法确定研究流域的最佳集水阈值更能准确的识别CSAs。目前子流域划分的最佳集水阈值确定方法有水系河网密度法、流域宽度分布推求、水系分形维数法等，其中河网密度法、流域宽度分布推求法适用于较大流域，而水系分维值反映流域地形地貌特征，更加适用于小流域的最佳集水阈值推求[6]。考虑到研究区流域面积，因此本文运用分形理论确定流域最佳集水阈值，基于SWAT模型进行合理的子流域划分与流域面源污染模拟，以期能够精确的识别流域重点治理区域。

1 研究区概况

古蔺县坐落于四川盆地的南缘，地处东经105°30′-106°20′，北纬27°41′-28°20′。古蔺河是古蔺县境内最大的河流，是赤水河右岸一级支流，流域占古蔺县各乡镇总面积的大部分，集水面积为965 km2，干流河长70.75 km，河道平均比降1.84%。根据古蔺县境内国家水文站1959-2007年实测的径流资料计算，流域多年的平均年径流深为422 mm，多年的平均流量5.07 m3/s，年径流量1.6 亿m3。径流在年内的分配不均匀，丰水期(5-9月)的径流占年径总流的65.8%，枯水期11月-翌年3月约占年径总流的21.7%。最枯月为12月，仅占年径流的3.7%。

研究区包含的8个乡镇，包括太平镇、石屏乡、鱼化乡、龙山镇、永乐镇、护家乡、古蔺镇、德跃镇；总人口为344 328，占全县人口的4.01% ，其中农业人口283 408，城镇人口60 920，各乡镇的平均人口密度为360人/km2，远远超过了中国平均人口密度(130 人/km2)。研究区包含的8个乡镇2013年种植小麦7 027 hm2，水稻5 087 hm2，玉米6 840 hm2，饲养牛36 795头，生猪151 802头，羊17 221只，鸡406 714只，其中包含集中式畜禽养殖与散养畜禽养殖。

图1 古蔺河流域Fig.1 Gulin River Basin

2 研究方法

2.1 模型数据库建立

(1)数字高程DEM：古蔺河流域的地形数据来源于中国科学院计算机网络信息中心国际科学数据镜像网站(http:∥datamirror.csdb.cn)提供的30m分辨率的DEM(Digital Elevation Model)数据。在DEM输入SWAT之前在Arcgis中进行处理，通过填洼处理避免水流逆流或者汇流[7-8]，SWAT中输入的研究区DEM见图2。

图2 数字高程Fig.2 DEM

(2)气象数据：古蔺河流域有古蔺水文站一个国家水文站点，位于古蔺镇，本次研究收集到2000-2013年逐时降雨数据，由于分布式水文模型中要使用日降雨数据，将其整理为日降雨，并收集到古蔺水文站2000-2010年的逐月径流数据；研究区的太阳辐射、日相对湿度、风速等气象数据使用SWAT内置的天气发生器模拟获得。

(3)土壤类型数据：研究区土壤数据是从FAO网上下载的HWSD 数据，并自行剪切出中国区的土壤栅格图，其数据值(Value)有效范围是11 000～11 935，古蔺河流域的土壤数据就是从中国区的土壤栅格图裁剪而来。在下载的HWSD数据库中打开 HWSD_DATA 表，可以查找土壤的部分属性值，其余属性值通过《四川土种志》查找与用SPAW软件计算得到[9]，并新增到SWAT模型的土壤数据库中。

表1 水文气象站点概况Tab.1 Hydrometeorology station

(4)土地利用数据：古蔺流域的土地利用/植被覆盖类型数据来源于四川省古蔺县林业局的相关土地利用信息，利用ARCGIS软件制成了古蔺县土地利用图。原始数据中包含了水田、旱地、有林地、草地、村庄、河流水面等各种土地利用分类。在分布式水文模型中，主要将土地利用类型分为9类，其中林地占54.08%(FRST),耕地占35.14%(AGRL),草地占7%(PAST)，因此古蔺县主要以林业和农业种植为主。由于美国的土地利用类型与我国的属性相似，因此在本次研究中土地利用属性直接采用SWAT模型系统默认的土地利用类型属性数据[10]。

图3 土地利用图Fig.3 Landuse

图4 土壤类型图Fig.4 Soil

(4)污染源数据：古蔺河流域的污染源中按照不同类型又可划分为若干种类；其中点源主要包括工业污染、城镇居民所产生的垃圾废水，非点源包括散养畜禽污染、化肥污染、农村居民所产生的垃圾废水，此类统计数据来源于古蔺县环保局、农业局、统计局等。本次研究假设各统计因子在行政区内等密度分布，在实际计算时，采用面积比例法，获得古蔺河流域的统计数据信息。

2.2 水系分维数计算

计盒维数法是求解分形维数的计算方法之一[11]，其基本思路是用不同尺寸的盒子覆盖分形体，并统计对应尺寸的盒子数量。当水系作为分形体时，可以利用Arcgis的空间分析模块选定某一集水阈值Si，提取研究区的栅格河网，然后不断地改变栅格边长r的大小，统计不同边长的栅格的数量N，当r→0 时，可以得到研究区水系分形维数Di：

(1)

3 结果与分析

3.1 最佳集水阈值确定

本文选取300 hm2到4 500 hm2之间共15个集水阈值，计算对应于10种栅格边长下的分维数。图5(a)为利用计盒维数法求得的分维数与集水阈值关系曲线，可以看出，随着集水阈值的不断增大，水系分维值开始减小，变化速度也逐渐降低。

图5 水系分维值与集水阈值关系Fig.5 The relationship of fractal dimension and threshold of drainage area

由于集水阈值越大，河网密度越小，所能反映的流域信息就越少，而集水阈值过小，又会导致流域出现过多的伪河道，使模拟结果与实际情况相差过大，因此需要保证流域不出现过多伪河道的情况下，河网密度尽量的大。在此基础上，建立分维数与集水阈值的二阶导函数关系(即分维值的变化加速度)如图5(b)，当集水阈值小于1 500 hm2时，分维值变化加速度极大，1 500 hm2后分维值变化加速度趋于稳定，说明该值即是相关研究中提到的阈值“拐点”[12,13]，因此本文选取集水阈值为1 500 hm2进行合理的子流域划分。

3.2 模型参数率定与验证

通过模型自带的敏感性分析模块对研究区参数进行敏感性分析，确定了流域中影响径流的6个敏感性参数，以及影响氮磷负荷敏感的6个参数，具体见表2。本次研究收集到古蔺河流域古蔺站2000-2010年实测月径流和古蔺河水质监测断面2012-2013年的水质数据。水质模拟需要控制断面的径流数据，因此首先率定和验证研究区的径流过程，调试模型使其精度满足模拟精度要求，一般要求校准期以及验证期径流模拟值和实测值的Nash效率系数Nash≥0.5，决定系数R2≥0.6[6]。

表2 敏感性参数表Tab.2 Sensitive parameters

(1)径流模拟：将2000年作为预热期，2001-2007年作为率定期，2008-2010年作为验证期，古蔺站径流模拟精度以及月径流实测值与模拟值对比见表3。可见SWAT模型在古蔺河流域有较好的适用性，率定期Nash效率系数0.78，决定系数R2为0.92；验证期Nash效率系数为0.71，决定系数R2为0.88；模拟结果令人满意。

表3 研究区径流模拟率定及验证结果Tab.3 Results of runoff simulation and calibration

(2)总氮、总磷模拟：将古蔺河流域出口太平渡断面2012年实测总氮、总磷数据(古蔺县环保局提供)用于模型率定，2013年的数据用于模型验证。氮、磷污染负荷模拟的精度要求对径流模拟的相对较低，只要保证氮、磷负荷的模拟值与实测值得Nash≥0.5，决定系数R2≥0.6即可。通过模型的调试，古蔺河流域水质监测断面的氮、磷污染负荷率定与验证结果满足精度要求，曲线拟合较好。TN与TP模拟值与实测值得月平均负荷对比见图6。

图6 总氮、总磷模拟结果Fig.6 results of SWAT in simulating TN and TP

3.3 总氮、总磷输出系数的空间分布

对古蔺河各子流域单位面积内的非点源污染负荷(即各子流域TN、TP的输出系数)进行空间分析(见图7)，由图中可以看出子流域TN、TP的输出系数空间分布特征一致，TN、TP输出系数最大的为6、11、14、18、30子流域，值分别为46.73～67.67、5.72～8.29 kg/hm2。根据研究区的概况及流域划分，11、14、30子流域属于古蔺镇的行政区内，子流域18属于石屏乡的行政区内，子流域6属于德跃镇的行政区内；而古蔺镇是古蔺河流域内面积最大，农业人口、畜禽养殖最多的乡镇，石屏乡的农业人口、畜禽养殖较多，但面积最小，因此造成行政区内的子流域污染物输出系数较大。

图7 各子流域单位面积TN、TP输出Fig.7 Unit area output of TN and TP in Subbasin

3.4 污染源贡献率分析

本次研究中的非点源污染包括农村居民生活所产生的污水、散养畜禽养殖以及化肥，分别模拟有相应污染源输入和无输入时研究区出口断面负荷的变化，根据变化的差值可以估算出各种污染源对研究区污染负荷的贡献。从图8中可以看出，大多数子流域中各类污染源的排列顺序为化肥产生的污染>畜禽养殖产生的污染>居民生活产生的污染；尤其是TP图中可以看出，化肥产生的污染TP产量远远高于畜禽养殖污染和化肥污染，可见化肥污染与畜禽养殖所造成的污染是古蔺河流域非点源污染的主要原因。

图8 各类污染源输出Fig.8 Output of each pollution source

4 结 论

(1)基于分形理论方法，计算出古蔺河流域最佳集水阈值为1 500 hm2，该阈值下子流域划分个数为37个。

(2)在古蔺河流域建立SWAT模型对氮磷负荷进行研究，模型率定期与验证期的径流模拟精度Nash效率系数值大于0.70，决定系数R2大于0.88；总氮、总磷模拟率定期与验证期Nash效率系数值大于0.5，决定系数R2大于0.8，模拟评价参数均满足精度要求，表明SWAT模型在古蔺河流域农村面源污染模拟研究中具有一定的适用性。

(3)古蔺河流域2013年总氮输出量为2 761.72 t，总磷输出量为338.27 t，古蔺镇、德跃镇、石屏乡污染负荷输出强度占前三，是需要重点治理的乡镇。

(4)农业化肥是流域主要污染源，农田中化肥施用的比例失调会大大降低农作物对化肥的利用率，同时增大化肥的流失，所以严格控制农业用肥量，进行测土配方施肥是研究区面源污染治理的需要采取重要措施。

□