SWAT水文模型

SWAT水文模型介绍
1概述
SWAT(Soil and Water Assessment Tool)模型就是美国农业部(USDA)农业研究局(ARS)开发的基于流域尺度的一个长时段的分布式流域水文模型。

它主要基于SWRRB模型,并吸取了CREAMS、GLEAMS、EPIC与ROTO的主要特征。

SWAT具有很强的物理基础,能够利用GIS与RS提供的空间数据信息模拟地表水与地下水的水量与水质,用来协助水资源管理,即预测与评估流域内水、泥沙与农业化学品管理所产生的影响。

该模型主要用于长期预测,对单一洪水事件的演算能力不强,模型主要由8个部分组成:水文、气象、泥沙、土壤温度、作物生长、营养物、农业管理与杀虫剂。

SWAT模型拥有参数自动率定模块,其采用的就是Q、Y、Duan 等在1992年提出的SCE-UA算法。

模型采用模块化编程,由各水文计算模块实现各水文过程模拟功能,其源代码公开,方便用户对模型的改进与维护。

2模型原理
SWAT模型在进行模拟时,首先根据DEM把流域划分为一定数目的子流域,子流域划分的大小可以根据定义形成河流所需要的最小集水区面积来调整,还可以通过增减子流域出口数量进行进一步调整。

然后在每一个子流域内再划分为水文响应单元HRU。

HRU就是同一个子流域内有着相同土地利用类型与土壤类型的区域。

每一个水文响应单元内的水平衡就是基于降水、地表径流、蒸散发、壤中流、渗透、地下水回流与河道运移损失来计算的。

地表径流估算一般采用SCS径流曲线法。

渗透模块采用存储演算方法,并结合裂隙流模型来预测通过每一个土壤层的流量,一旦水渗透到根区底层以下则成为地下水或产生回流。

在土壤剖面中壤中流的计算与渗透同时进行。

每一层土壤中的壤中流采用动力蓄水水库来模拟。

河道中流量演算采用变动存储系数法或马斯金根演算法。

模型中提供了三种估算潜在蒸散发量的计算方法—Hargreaves、Priestley-Taylor与Penman-Monteith。

每一个子流域内侵蚀与泥沙量的估算采用改进的USLE方程,河道内泥沙演算采用
改进的Bagnold泥沙运移方程。

植物吸收的氮采用供需方法计算,植物的氮日需求量就是植物与生物量中氮浓度的函数。

土壤中向植物供给氮,当需求超过供给时,出现营养物压力。

地表径流、壤中流与渗透过程运移的硝态氮量由水量与土壤层中的平均硝态氮浓度来估计。

泥沙中运移的有机氮采用McElroy et al、开发的负荷方程,后经进一步改进。

该负荷方程基于土壤表层的有机氮浓度、泥沙量与富集率来估计径流中的有机氮损失。

植物吸引的磷采用与氮相似的供需方法。

径流中带走的可溶解磷采用土壤表层中的不稳定磷、径流量与磷土分离系数来计算。

泥沙运移的磷采用与有机氮运移相同的方程。

河道中营养物的动态模拟采用QUAL2E模型。

2、1产汇流模型
模型中采用的水量平衡表达式为:
式中:为土壤最终含水量,mm;为土壤前期含水量,mm;t为时间步长,d;为第i天降水量,mm;为第i天的地表径流,mm;为第i天的蒸发量,mm;为第i天存在于土壤剖面底层的渗透量与侧流量,mm;为第i天地下水含量,mm。

SWAT模型水文循环陆地阶段主要有水文、天气、沉积、土壤温度、作物产量、营养物质与农业管理等部分组成。

模型径流量产生计算流程图如下图2-1所示。

图2-1 SWAT模型产流计算流程图
2、1、1地表径流
当落到地表的降水量多余入渗量时产生地表径流。

SWAT模型采用SCS径流曲线法计算。

SCS曲线方程自上世纪50年代逐渐得到广泛使用,属于经验模型,就是对全美小流域降水与径流关系20多年的研究成果。

模型能反应不同土壤类型与土地利用方式及前期土壤含水量对降雨径流的影响,它就是基于流域的实际入渗量(F)与实际径流量(Q)之比等于流域该场降雨前的最大可能入渗量(S)与最大可能径流量()之比的假定基础上建立的
SCS模型的降雨-径流基本关系表达式如下:
式中:假定潜在径流量()为降水量(P)与由径流产生前植物截留、初渗与填洼蓄水构成的流域初损()的差值。

由此推到上式有:
初损受土地利用、耕作方式、灌溉条件、冠层截留、下渗、填洼等因素的影响,它与土壤最大可能入渗量S呈一定的正比关系,美国农业部土壤保持局在分析了大量长期的实验结果基础上,提出了二者最合适的比例系数为0、2,即:
由此可得SCS方程为:
流域当时最大可能滞留量在空间上与土地利用方式、土壤类型与坡度等下垫面因素密切相关,模型引入的值可较好地确定,公式如下:
CN就是一个无量纲参数,CN值就是反映降雨前期流域特征的一个综合参数,它就是前期土壤湿度、坡度、土地利用方式与土壤类型状况等因素的综合。

2、1、2蒸散发
模型考虑的蒸散发就是指所有地表水转化为水蒸气的过程,包括树冠截留的水分蒸发、蒸腾与升华及土壤水的蒸发。

蒸散发就是水分转移出流域的主要途径,在许多江河流域,蒸发量都大于径流量。

准确地评价蒸散发量就是估算水资源量的关键,也就是研究气候与土地覆盖变化对河川径流影响的关键。

(1)潜在蒸散发
模型提供了Penman-Monteith、Priestley-Taylor与Hargreaves三种计算潜在蒸散发的方法,另外还可以使用实测资料或已经计算好的逐日潜在蒸散发资料。

一般采用Penman-Monteith方法来计算流域的潜在蒸散发。

(2)实际蒸散发
实际蒸散发以潜在蒸散发为计算基础。

在计算流域实际蒸散发量的时候,模型首先计算植物冠层截留水分的蒸发,然后计算最大蒸腾量、最大升华量与最大
土壤蒸发量,最后计算实际的升华量与土壤水分蒸发量。

(3)冠层截留蒸发量
模型在计算实际蒸发时假定尽可能蒸发冠层截留的水分,如果潜在蒸发量小于冠层截留的自由水量,则:
式中:为某日流域的实际蒸发量,mm;为某日冠层自由水蒸发量,mm;为某日的潜在蒸发量,mm;为某日植被冠层自由水初始含量,mm;为某日植被冠层自由水终止含量,mm。

如果潜在蒸发量大于冠层截留的自由水含量则:
当植被冠层截留的自由水被全部蒸发掉,继续蒸发所需的水分就会从植被与土壤中得到。

(4)植物蒸腾
假设植物生长在一个理想的条件下,植物蒸腾可用以下表达式计算:
当0≤LAI≤3、0时
当LAI＞3、0时
式中:为某日最大蒸腾量,mm;为植被冠层自由水蒸发调整后的潜在蒸发,mm;LAI为叶面积指数。

因为没有考虑到植物下面图层的含水量问题,由此公式计算处的蒸腾量可能比实际蒸腾量要大一些。

(5)土壤水分蒸发
在计算土壤水分蒸发时,首先区分出不同深度土壤层所需要的蒸发量,土壤深度层次的划分决定土壤允许的最大蒸发量,可由下式计算:
式中:为z深度处蒸发需要的水量,mm;z为地表以下土壤的深度,mm。

表达式中的系数就是为了满足50%的蒸发所需水分来自土壤表层10mm,以及95%的蒸发所需水分来自0～100mm土壤深度范围内。

土壤水分蒸发所需要的水量就是有土壤上层蒸发需水量与土壤下层蒸发需水量决定的:
式中:为ly层的蒸发需水量,mm;为土壤下层的蒸发需水量,mm;
为土壤上层的蒸发需水量,mm。

土壤深度的划分假设50%的蒸发需水量由0～10mm内土壤上层的含水量提供,因此100mm的蒸发需水量中50mm都要由10mm的上层土壤提供,显然上层无法满足需要,这就需要建立一个系数来调整土壤层深度的划分,以满足蒸发需水量,调整后的公式可以表示为:
式中:esco为土壤蒸发调节系数,该系数就是SWAT为调整土壤因毛细作用与土壤裂隙等因素对不同土层蒸发量二提出的,对于不同的esco值对应着相应的土壤层划分深度。

2、1、3土壤水
渗入到土壤中的水有多种不同运动方式。

土壤水可以被植物吸收或蒸腾而损耗,可以渗透到土壤底层最终补给地下水,也可以在地表形成径流,即壤中流。

由于主要考虑径流量的多少,因此对壤中流的计算简要概括。

模型采用动力储水方法计算壤中流。

相对饱与区厚度计算公式为:
式中:为土壤饱与区内可流出的水量,mm;为山坡坡长,m;为土壤可出流的孔隙率;表示土壤层总孔隙度,即与土壤层水分含量达到田间持水量的孔隙度之差。

山坡出口断面的净水量为:
式中:为出口断面处的流速,mm/h。

其表达式为:
式中:为土壤饱与导水率,mm/h;为坡度。

总结上面表达式,模型中壤中流最终计算公式为:
2、1、4地下水
模型采用以下表达式来计算流域地下水:
式中:为第i天进入河道的地下水补给量,mm;为第(i-1)天进入河道的地下水补给量,mm;为时间步长,d;为第i天蓄水层的补给流量,mm;为基流的退水系数。

其中补给流量由下式计算:
式中:为第i天蓄水层补给量,mm;为补给滞后时间,d;为第i天通过土壤剖面底部进入地下含水层的水分通量,mm/d;为第(i-1)天蓄水层补给量,mm。

2、2土壤侵蚀模型
泥沙生成量用MUSLE方程来预测,计算渠道泥沙输移量的公式为:
式中:为输移能力,t/m3;为流速,m/s;a与b就是常数。

根据天气条件,泥沙输移量可以高于或者低于输移能力,导致沉积过量的泥沙通过渠道侵蚀再悬浮输移泥沙。

流速方程为:
式中:为流量,m3/s;为渠道宽度,m;为径流深,m。

对于低于齐岸深度的径流,径流深使用Manning方程来计算,假定渠道宽度远大于深度:
式中:为渠道曼宁系数;为渠道坡度,m/m。

由于降水与径流产生的土壤侵蚀就是用MUSLE方程来计算的,MUSLE就是修正的通用土壤流失方程(USLE)。

USLE方程就是通过降水动能函数预测年均侵蚀量,而在MUSLE中,用径流因子代替降水动能,改善了泥沙产量的预测,这样就不需要泥沙输移系数,并且可以将方程用于单次暴雨事件,因为径流因子就是先行湿度与降水动能的函数。

USLE中需要输移系数就是因为降水动能因子表示的能量只在作用流域内起作用。

修正的通用土壤流失方程为;
式中:为土壤侵蚀量,t;为地表径流,mm/h;为洪峰流量,m3/s;为水文响应单元的面积,hm2;为土壤侵蚀因子;为植被覆盖与管理因子;
为水土保持措施因子;为地形因子;为粗碎屑因子。

2、2、1土壤侵蚀因子
当其她影响侵蚀的因子不变时,因子反映不同类型土壤抵抗侵蚀力的高低。

它与土壤物理性质的影响,如机械组成、有机质含量、土壤结构、土壤渗透性等有关。

当土壤颗粒粗、渗透性大时,值就低,反之则高;一般情况下值得变幅在0、02 ~0、75之间。

值得直接测定方法就是:在标准小区(坡长为22、1m,宽为1、83m,坡度为9%)上没有任何植被,完全休闲,无水土保持措施,降水后收集由于坡面径流而冲蚀到集流槽内的土壤,烘干、称重,由公式计算得到值。

试验测算值既费时又费力,1971年一些学者发展了一个通用方程来计算土壤侵蚀因子值,该方程在土壤黏土与壤土组成少于70%时适用。

式中:为颗粒尺度参数;为有机物含量百分比,%;为土壤分类中的结构代码;为土壤剖面可渗透性类别。

1995年另一位学者提出了一个替换方程:
式中:为粗糙砂土质地土壤侵蚀因子;为粘壤土土壤侵蚀因子;为土壤有机质因子;为高砂质土壤侵蚀因子。

各因子的计算公式如下:
式中:为粒径在0、05～2、00mm沙粒的百分含量;为粒径在0、002～0、05mm的淤泥、细砂百分含量;为粒径小于0、002mm的粘土百分含量;为各土壤层中有机碳含量,%。

2、2、2植被覆盖因子
植被覆盖与管理因子表示植物覆盖与作物栽培措施对防止土壤侵蚀的综合效益,其含义就是在地形、土壤、降水条件相同的情况下,种植作物或林草地的土地与连续休养地土壤流失量的比值,最大取值为1、0、由于植被覆盖受植物生长期的影响,SWAT模型通过下面的方程调整植被覆盖与管理因子:
式中:就是最小植被覆盖与管理因子值;就是地表植物残留量,kg/hm2。

最小C因子可以由已知年平均C值,通过以下方程计算。

式中:表示不同植被覆盖的年均C值。

2、2、3水土保持措施因子
水土保持因子就是指有保持措施的地表土壤流失与不采取任何措施的地表土壤流失的比值,这里的水土保持措施包括等高耕作、带状种植与梯田。

等高耕作对于中低强度的降水侵蚀具有保护水土流失的作用,但对于高强度的降水其保护作用则很小,等高耕作对坡度为3%～8%之间的土地非常有效。

2、2、4地形因子
地形因子的计算公式如下:
式中:为坡长;为坡长指数;为坡度。

坡长指数的计算公式如下:
式中:为水文相应单元的坡度,。

2、2、5因子
因子就是通过下面公式计算的:
式中:为第一层土壤中砾石的百分比,%。

2、3污染负荷模型
SWAT模型可以模拟不同形态氮的迁移转化过程,地表径流流失、入渗淋失、化肥输入等物理过程,有机氮矿化、反硝化等化学过程以及作物吸收等生物过程,氮可以分为有机氮、作物氮与硝酸盐氮三种化学状态,氮的生物固定、有机氮向无机氮的转化以及溶解性氮随侧向壤中流的迁移等过程,有机氮又被划分为活泼有机氮与惰性有机氮两种状态,以及铵态氮挥发过程的模拟。

2、3、1硝态氮污染负荷模型
硝态氮主要随地表径流、侧向流或渗流在水体中迁移,要计算随水体迁移的硝态氮量必须先计算自由水中硝态氮浓度,用这个浓度乘以各个水路流动水的总量,即可得到从土壤中流失的硝态氮总量。

自由水部分的硝态氮浓度可用下面公式计算:
式中:为自由水中硝态氮浓度,kg/mm;为土壤中硝态氮的量,kg/hm2;为土壤中自由水的量,mm;为孔隙度;为土壤饱与含水量。

(1)通过地表径流流失的溶解态氮计算公式:
式中:为通过地表径流流失的硝态氮,kg/hm2;为硝态氮渗流系数;
为地表径流,mm。

(2)通过侧向流流失的溶解态氮的量计算公式:
对于地表10mm图层:
对于10mm以下的土层:
式中:为通过侧向流流失的硝态氮,kg/hm2;为硝态氮渗流系数;为自由水的硝态氮浓度,kg/mm;为侧向流,mm。

(3)通过渗流流失的溶解态氮量计算公式:
式中:为通过渗流流失的硝态氮,kg/hm2;为自由水的硝态氮浓度,kg/mm;为渗流,mm。

2、3、2有机氮污染负荷模型
有机氮通常就是吸附在土壤颗粒上随径流迁移的,这种形式的氮负荷与土壤流失量密切相关,土壤流失量直接反映了有机氮负荷。

有机氮随土壤流失的输移
量计算公式为:
式中:为有机氮流失量,kg/hm2;为有机氮在表层土壤中的浓度,kg/t;为土壤流失量,t;为水文相应单元的面积,hm2;为氮富集系数,氮富集系数就是随土壤流失的有机氮浓度与土壤表层有机氮浓度的比值。

计算富集系数的公式如下:
式中:为地表径流中泥沙含量。

的计算公式如下:
式中:为土壤流失量,t;为水文相应单元面积,hm2;为地表径流,mm。

2、3、3溶解态磷污染负荷模型
溶解态磷在土壤中的迁移主要就是通过扩散作用实现的,扩散就是指离子在微小尺度下由于浓度梯度而引起的溶质迁移,由于溶解态磷不很活跃,所以由地表径流以溶解态形式带走的土壤表层的磷很少,地表径流输移的溶解态磷可由下面公式计算:
式中:为通过地表径流流失的溶解态磷,kg/hm2;为土壤中溶解态磷,kg/hm2;为土壤溶质密度,mg/m3;为表层土壤深度,mm;为土壤磷分配系数,表层土壤中溶解态磷的浓度与地表径流中溶解态磷浓度的比值。

2、3、4有机磷与矿物质磷污染负荷模型
有机磷与矿物质磷通常就是吸附在土壤颗粒上通过径流迁移的,这种形式的
磷负荷与土壤流失量密切相关,土壤流失量直接反映了有机磷与矿物质磷负荷,有机磷与矿物质磷随土壤流失输移量计算公式为:
式中:为有机磷流失量,kg/hm2;为有机磷在表层土壤中的浓度,kg/t;为土壤流失量,t;为水文相应单元的面积,hm2;为磷富集系数。

2、4河道演算模型
模型中定义的河道均就是明渠流,SWAT用曼宁公式来定义河道糙率与水流流速。

水流在河道中演进过程使用变储量演算法或马斯金根法,两种方法都就是动力波方程。

在模拟中一般使用马斯金根法来计算。

(1)马斯金根法
马斯金根法假设河道内水体形状就是由一个楔形蓄水体与一个棱形蓄水体组成,如下图2-2所示。

图2-2 河道水体示意图
当洪水波行进到某个河段槽,入流量大于出流量便形成了楔形蓄水体。

当洪水波退去,在河段槽便出现了出流量大于入流量的负楔蓄体。

另外对于楔蓄水体,河段槽内始终包含一个体积为流域长度上横截面不变的棱柱状水体。

总的蓄水容量为:
式中:为蓄水容量,m3;为入流量,m3/s;为出流量,m3/s;为稳定流情况下的河段传播时间;为流量比重因素。

该公式可以重新整理为如下形式:
流量比重因素的下限为0、0,上限为0、5。

这个因子就是楔蓄量的函数。

对于水库式蓄水,没有楔蓄,;而对于一个完整的楔蓄,;对于河流,落在0、0与0、3之间,其平均值接近0、2。

对于蓄水容量的定义可以加入连续公式并简化为:
式中:为该时间段开始时的入流量,m3/s;为该时间段结束时的入流量,m3/s;为该时间段结束时的出流量,m3/s;为该时间段结束时的出流量,m3/s。

其中,。

为了用体积单位表示所有值,在蓄水容量公式两端乘以该时间段得到:
为了保持数值稳定与避免出现负出流量的计算,必须满足以下条件:
流量比重因素X的值由使用者输入,蓄水时间常数的值估算如下:
式中:为稳定流情况下的河段传播时间,s;与为权重系数,由使用者输入;为稳定刘情况下渠道蓄满水的河段传播时间,s;为渠道蓄满1/10水量时河段传播时间,s。

要计算与,Cunge于1969年提出一个公式:
式中:为渠道长度,km;为指定深度处的波速,m/s。

波速的计算公式为:
其中流速由曼宁公式求解,将曼宁公式代入上式得:
式中:为湿周,m;为河段坡度,%;为曼宁系数;为流速,m/s。

(2)变储量演算法
对于一个给定的河段,储量演算基于连续方程,可写为:
式中:为稳定流情况下的河段传播时间,s;为该时间段开始时的入流量,m3/s;为该时间段结束时的入流量,m3/s;为该时间段开始时的出流量,m3/s;为该时间段结束时的出流量,m3/s;为该时间段开始时的蓄水容量,m3;为该时间段结束时的蓄水容量,m3。

演进时间就是由渠道中的水容量除以水流流量:
式中:为蓄水容量,m3;为出流量,m3/s。

联合以上两式,简化为:
其中,为蓄水系数,;为入流平均流量,
,m3/s。

整理上式得:
为用体积单位来表达所有变量,公式两边都乘以时间段,则有:
3软件操作
3、1输入数据
要创建SWAT数据集,ArcSWAT需要访问提供流域信息的准确类型的ArcGIS 的栅格、矢量与数据文件。

这些必需的地图主题与数据文件必须在运行SWAT之前准备好。

3、1、1空间数据
空间数据包括必需的ARCSWAT空间数据与可选的ARCSWAT空间数据。

SWAT 空间数据集可以以任何投影类型来创建,但所有的地图必须用同一投影。

(1)必需的ARCSWAT空间数据
DEM,ESRI GRID Format。

DEM的高程值可以用整型或者实数型。

决定地图分辨率的单位不必与高程的单位保持一致。

例如地图分辨率可以就是米
(m),而高程可以就是英尺(feet)。

地图分辨率单位必须定义为以下几种:
米(m),公里(km),英尺(feet),码(yards),英里(miles)。

高程单位必须定义为以
下几种:米(m),厘米(cm),码(yards),英尺(feet),英寸(inches)。

●Land Cover/Land Use,ESRI GRID,Shapefile,or Feature Class Format。

土地利
用图中的类别需要重新分类为SWAT需要的土地利用类型。

用户可以用
三种方法重新分类土地利用,第一种就是创建地图时用美国地质勘探局
的分类代码;第二种就是加载土地利用图时为每一分类选定SWAT土地利
用类型;第三种就是为土地利用图的不同分类创建一个能够识别4位
SWAT代码的look up表。

●Soil,ESRI GRID,Shapefile,or Feature Class Format。

用户需要输入SWAT土壤
文件,在创建工程之前把图中每一类的土壤数据输入User Soil database
中。

要对土壤图的土壤类型进行重分类,信息必须在界面中手动输入,列
出土壤信息look up表
(2)可选的ARCSWAT空间数据集
●DEM Mask,ESRI GRID,Shapefile or Feature Class Format。

模型可以加载一个
mask在DEM之上。

●STREAMS,Shapefile or Feature Class Format。

在那些地势很低的地方,DEM
地图网格不能精确地判断河流的位置,就需要有河网描绘的shapefile线
性文件。

●User-Defined Watersheds,Shapefile or Feature Class Format。

用户可以自定
义子流域划分加载到模型中去,如果选择自定义子流域也必须自定义河
网。

●User-Defined Streams,Shapefile or Feature Class Format。

用户自定义子流
域后,模型将不会进行河网提取,需要用户自行定义。

3、1、2属性数据
属性数据包括:土壤属性数据、气象资料数据与水文数据。

SWAT模型通过三个个数据库来存储属性数据,其中气象资料数据通过模型自带的数据库进行存储,另外自定义的土壤属性数据与气象站参数数据通过两个附加的数据库进行存储,水文数据用于模拟结果率定。

(1)土壤属性数据
模型需要的土壤属性数据包括两大类:土壤物理属性与土壤化学属性。

土壤的物理属性对土壤剖面中水与气的运动状况起决定作用,并显著影响水文响应单元(HRU)中水文过程的模拟。

土壤物理属性参数主要包括土层厚度、土壤层数、密度、土壤颗粒大小分布与土壤饱与水力传导率等。

由于获得的土壤属性数据难以完全满足构建SWAT模型土壤库的要求,因此,部分无法直接获取的参数如土壤可利用有效水、饱与水力传到率等采用美国农业部开发的土壤水特性程序SPAW进行估算;部分参数通过对前人的研究成果与经验公式获得。

估算结果仅简略反映土壤参数特征,在模型参数校准的过程中还需做进一步调整。

(2)气象、水文资料数据
SWAT模型所需要的气象数据包括气象站点的位置以及各气象站点的实测数据。

气象站点包括雨量站、湿度站、温度站、风速站、太阳辐射站,其实测数据分别为日降水量、相对湿度、日最高/最低气温、风速与太阳辐射。

本文采用北京站的气象资料来构建模型“天气发生器”所需的数据。

SWAT模型自带的天气发生器WXGEN,其功能有二,一就是用于生成气候数据,通过大量基础气象数据构建完成之后,可以用来生成任意年份的气象数据,二就是对缺失的数据进行填补。

3、2输出数据
SWAT模型每次模拟都会生成五个输出文件:输入汇总文件(input、std)、输出汇总文件(output、std)、HRU输出文件(output、hru)、子流域输出文件(output、sub)、河道输出文件(output、rch)。

3、2、1 output、hru
HRU输出文件包含了流域中每一个水文响应单元的信息,文件中输出数据介绍如下:
MON:时间步长,根据用户选择可以为日、月、年;
AREA:HRU的面积,km2;
PRECIP:时间步长内HRU的降雨量,mm;
SNOFALL:时间步长内HRU的降雪量,mm;。