一种大流域尺度水氮迁移耦合模拟方法与流程

本发明涉及寒区氮循环技术领域，尤其涉及一种大流域尺度水氮迁移耦合模拟方法。

背景技术：

由于人口数量的增长以及工业化的加快发展，氮肥施用量、畜禽养殖量、城镇农村生活污水量等也在逐渐增长并且通过人类活动进入到生态系统，破坏了原有生态系统的平衡，造成土壤污染、水体富营养化等问题。在寒区土壤冻结过程中会有大量的面源氮污染物持留，在融冻过程中伴随着冻土和冰雪融化过程以及春季降水进入地表和地下水体，导致水质迅速下降，造成寒冷地区面源污染物流出具有累积性和突发性。

此外，国内主流的流域水文模型之一的大尺度流域水与能量转化过程模拟(waterandenergytranferprocessesinlargeriverbasins，简称wep-l)可以很好的模拟冻融过程中地表、土壤、地下水、坡面沟以及河道中的水分迁移；该模型为水文模型，未考虑土壤中污染物的迁移过程。

技术实现要素：

本发明的目的在于提供一种大流域尺度水氮迁移耦合模拟方法，从而解决现有技术中存在的前述问题。

为了实现上述目的，本发明采用的技术方案如下：

一种大流域尺度水氮迁移耦合模拟方法，包括如下步骤，

s1、根据地表水文过程，结合地表氮素总量，分别计算由于氮沉降和施肥沉降到地表的随地表径流流失、下渗到土壤以及储存在地表洼地的溶解态氮分量；

s2、以土壤水文过程为基础，计算源汇项，并根据源汇项利用对流弥散方程计算土壤各层氮素浓度；

s3、根据地下水水文过程，利用质量守恒方程计算地下水氮素浓度；

s4、根据坡面沟水文过程，利用对流弥散方程计算坡面沟中氮素浓度；

s5、根据河道水文过程，利用对流弥散方程计算河道断面氮素浓度。

优选的，步骤s1具体为，根据地表水文过程，计算由于氮沉降和施肥沉积到地表的氮素随地表径流流失、下渗到土壤以及储存在地表洼地的比例，并利用比例乘以地表氮素总量分别计算随地表径流流失、下渗到土壤以及储存在地表洼地的溶解态氮量；计算公式分别如下，

min＝mfer+mwet

mrun＝raterun·min

minf＝rateinf·min

mstore＝ratestore·min

其中，rsuf_c为地表径流量；inf_c为地表下渗到土壤的水量；rstore_c为洼地储留量；raterun为地表径流流失的水量占地表总水量的比例；rateinf为下渗到土壤的水量占地表总水量的比例；ratestore为储存在地表的水量占地表总水量的比例；min为地表氮的输入总量；mfer为地表施肥量；mwet为地表氮沉降量；mrun为随地表径流流失的溶解态氮量；minf为下渗到土壤的溶解态氮量；mstore为储存在地表的溶解态氮量。

优选的，步骤s2具体为，根据土壤水文过程，计算源汇相，源汇项包括由于土壤吸力引起的上下迁移量、氮素转化量以及通过壤中流流失量；进而利用源汇项使用对流弥散方程计算土壤各层氮素的浓度；

利用下式计算在dt时间内的源汇项

其中，v为土壤体积；mfer为施氮量；mwet为每日氮沉降量；mupdown为由于土壤吸引力引起的氮素在相邻两层的迁移量；mtran为氮素不同成分之间的转化量，mtran＝-μici，μi为各种反应系数，ci表示土壤中污染物的浓度；muptake为植物根系吸收量；mrunoff为壤中流损失量；

在冻融期，施氮量和植物根系吸收量为零；土壤分为11层，只在第一层土壤考虑氮沉降与施肥，层与层之间考虑由于重力势引起的下渗，由于吸引力引起的土壤水分在上下层的迁移量，进而引起氮素在上下层的迁移；

利用对流弥散方程计算土壤各层氮素浓度的计算公式如下，

其中，cdis为土壤水中溶解态氮的浓度；csorb为吸附态氮的浓度，csorb＝kd·cdis，kd表示固相、液相分配系数；·为土壤容重；v为土壤水流通量；z为表示各层土壤相对于基准面的深度；θ为t时刻位于深度z处的土壤体积含水率；d为水动力弥散系数。

优选的，步骤s3具体为，根据地下水水文过程，统计土壤下渗到地下水的氮量、该等高带地下水与周围等高带地下水氮素交换量、地下水与河道的交换量以及地下水氮素转化量，并利用质量守恒方程计算地下水氮素浓度；部分溶解态有机氮、氨氮和硝态氮会进一步淋溶到地下水中，其在地下水中的迁移转化过程采用质量守恒方程描述，

其中，p为饱和含水率；h为饱和带的厚度；qr为非饱和土壤带向饱和土壤带的渗漏速率；cri为非饱和带渗漏水中溶解态氮的浓度；qriver为地下水与河道交换水量；qdgd为地下水与周围等高带地下水的交换水量。

优选的，步骤s4具体为，根据坡面沟水文过程，统计由地表径流以及壤中流进入坡面沟的氮量以及坡面沟中氮素转化量，利用对流弥散方程计算坡面沟中氮素浓度；除对流扩散过程之外还包括旁侧入流、一个等高带的坡面沟到下一个等高带的坡面沟以及不同组分之间的相互转化作用，具体计算公式为，

其中，c′dis为河道中溶解态氮的浓度；v为水流的流速；ex为扩散系数；s为有机氮、氨氮或硝态氮的源汇项，用于表示旁侧入流以及不同组分之间的相互转化作用，具体指代哪一种氮素的源汇项，需要根据计算的氮素种类进行确定。

优选的，对于有机氮，矿化作用会使有机氮转化为氨氮从而减少，有机氮在坡面沟中或河道中的源汇项可表示为，

其中，m1为有机氮的旁侧入流；dx为河流的长度；a为河流的横截面积；k1为有机氮的矿化速率；c1为河流中有机氮的浓度；

对于氨氮，矿化作用会使有机氮转化为氨氮从而使河流中氨氮含量增加，硝化作用会使氨氮转化为硝态氮从而使河流中氨氮含量降低；则氨氮在坡面沟中或河道中的源汇项可表示为：

其中，m2表示氨氮的旁侧入流；k2表示氨氮的硝化速率；c2表示河流中氨氮的浓度；

对于硝态氮，硝化作用会使氨氮转化为硝态氮从而使河流中硝态氮含量增大，反硝化作用可使硝态氮转化为n2o，使河流中硝态氮含量降低；则硝态氮在坡面沟中或河道中的源汇项可表示为：

式中，m3表示硝态氮的旁侧入流；k3表示硝态氮的反硝化速率；c3表示河流中硝态氮的浓度。

优选的，步骤s5具体为，根据河道水文过程，统计坡面沟流入河道的量以及地下水与河道交换量，计算河道中氮素转化量，利用对流弥散方程计算河道断面氮素浓度；河道断面氮素浓度的计算方式与步骤s4相同；溶解态有机氮、氨氮和硝态氮在坡面沟与河道中的迁移转化过程均使用对流弥散方程进行描述，除对流扩散过程之外还包括旁侧入流、上一个子流域的河道到下一个子流域的河道以及不同组分之间的相互转化作用。

优选的，若坡面沟或河道窄而直，则ex＝5.93·d·v，d表示河流深度；

若坡面沟或河道宽大，则ex＝0.11·v·w²/d，w表示坡面沟或河道宽度。

优选的，河道中有机氮、氨氮以及硝态氮的旁侧入流包括地下水与河道的交换量以及该子流域各等高带的坡面沟汇入河道的氮量。

本发明的有益效果是：1、考虑了地表、土壤、地下水、坡面沟与河道中氮素的迁移转化过程，溶解态氮随水移动，考虑了随各个阶段氮素的迁移。2、使用对流弥散方程计算冻融过程中土壤、坡面沟、河道污染物的迁移过程，提高了水文过程模拟的准确性，并且可以很好的模拟冻土融化期河道断面流量和污染物氮素的浓度。

附图说明

图1是本发明实施例中模拟方法的流程图；

图2是本发明实施例中地表氮迁移转化模拟流程图；

图3是本发明实施例中土壤中氮迁移转化模拟流程图；

图4是本发明实施例中地下水中氮迁移转化模拟流程图；

图5是本发明实施例中坡面沟中氮迁移转化模拟流程图；

图6是本发明实施例中河道中氮迁移转化模拟流程图；

图7是本发明实施例中流域硝态氮模拟与实测对比图。

具体实施方式

为了使本发明的目的、技术方案及优点更加清楚明白，以下结合附图，对本发明进行进一步详细说明。应当理解，此处所描述的具体实施方式仅仅用以解释本发明，并不用于限定本发明。

如图1至图7所示，本实施例中，提供了一种大流域尺度水氮迁移耦合模拟方法，包括如下步骤，

s1、根据地表水文过程，结合地表氮素总量，分别计算由于氮沉降和施肥沉降到地表的随地表径流流失、下渗到土壤以及储存在地表洼地的溶解态氮分量；

s2、以土壤水文过程为基础，计算源汇项，并根据源汇项利用对流弥散方程计算土壤各层氮素的浓度；

s3、根据地下水水文过程，利用质量守恒方程计算地下水氮素浓度；

s4、根据坡面沟水文过程，利用对流弥散方程计算坡面沟中氮素浓度；

s5、根据河道水文过程，利用对流弥散方程计算河道断面氮素浓度。

本实施例中，模拟方法主要包括五部分内容，具体为，地表氮迁移转化模拟、土壤中氮迁移转化模拟、地下水中氮迁移转化模拟、坡面沟中氮迁移转化模拟和河道中氮迁移转化模拟，下面分别针对这五部分内容进行详细的解释说明。

一、地表氮迁移转化模拟

该部分内容主要对应步骤s1，如图2所示，步骤s1具体为，根据地表水文过程，计算由于氮沉降和施肥沉积到地表的氮素随地表径流流失、下渗到土壤以及储存在地表洼地的比例，并利用比例乘以地表氮素总量分别计算随地表径流流失、下渗到土壤以及储存在地表洼地的溶解态氮量；计算公式分别如下，

min＝mfer+mwet

mrun＝raterun·min

minf＝rateinf·min

mstore＝ratestore·min

其中，rsuf_c为地表径流量，m；inf_c为地表下渗到土壤的水量，m；rstore_c为洼地储留量，m；raterun为地表径流流失的水量占地表总水量的比例；rateinf为下渗到土壤的水量占地表总水量的比例；ratestore为储存在地表的水量占地表总水量的比例；min为地表氮的输入总量，kg；mfer为地表施肥量，kg；mwet为地表氮沉降量(该地表氮沉降量主要是氨氮和硝态氮的沉降量)，kg；mrun为随地表径流流失的溶解态氮量，kg；minf为下渗到土壤的溶解态氮量，kg；mstore为储存在地表的溶解态氮量，kg。

二、土壤中氮迁移转化模拟

该部分内容对应步骤s2，如图3所示，步骤s2具体为，根据土壤水文过程，计算源汇相，源汇项包括由于土壤吸力引起的上下迁移量、使用一阶动力学反应方程计算氮素转化量以及通过壤中流流失量，在冻融期，施氮量以及植物根系吸收量为0；进而利用源汇项使用对流弥散方程计算土壤各层氮素的浓度；其中源汇项和土壤各层氮素的浓度分别计算如下：

1、利用下式计算在dt时间内的源汇项

其中，v为土壤体积，m³；mfer为施氮量，kg；mwet为每日氮沉降量，kg；mupdown为由于土壤吸引力引起的氮素在相邻两层的迁移量，kg；mtran为氮素不同成分之间的转化量(主要包括矿化作用、硝化作用和反硝化作用，可用一阶动力学反应方程进行计算)，kg，mtran＝-μici，μi为各种反应系数，ci表示土壤中污染物的浓度；muptake为植物根系吸收量，kg；mrunoff为壤中流损失量，kg；

土壤分为11层，只在第一层土壤考虑氮沉降与施肥，层与层之间考虑由于重力势引起的下渗，由于吸引力引起的土壤水分在上下层的迁移量，进而引起氮素在上下层的迁移。

2、利用对流弥散方程计算土壤各层氮素浓度(在非饱和带土壤，污染物的垂向运移主要采用对流-弥散方程进行计算)的计算公式如下，

其中，cdis为土壤水中溶解态氮的浓度，kg/m³；csorb为吸附态氮的浓度，kg/kg，csorb＝kd·cdis，kd表示固相、液相分配系数，m³/kg；ρ为土壤容重，kg/m³；v为土壤水流通量，m/d；z为表示各层土壤相对于基准面的深度，m；θ为t时刻位于深度z处的土壤体积含水率，m³/m³；d为水动力弥散系数，m²/d。

三、地下水中氮迁移转化模拟

该部分内容对应步骤s3，如图4所示，步骤s3具体为，根据地下水水文过程，统计土壤下渗到地下水的氮量、该等高带地下水与周围等高带地下水氮素交换量、地下水与河道的交换量以及地下水氮素转化量，并利用质量守恒方程计算地下水氮素浓度；

部分溶解态有机氮、氨氮和硝态氮会进一步淋溶到地下水中，其在地下水中的迁移转化过程采用质量守恒方程描述，

其中，p为饱和含水率，m³/m³；h为饱和带的厚度，m；qr为非饱和土壤带向饱和土壤带的渗漏速率，m/d；cri为非饱和带渗漏水中溶解态氮的浓度，kg/m³；qriver为地下水与河道交换水量，m/d；qdgd为地下水与周围等高带地下水的交换水量，m/d。

四、坡面沟中氮迁移转化模拟

该部分对应步骤s4的内容，如图5所示，s4具体为，根据坡面沟水文过程，统计由地表径流以及壤中流进入坡面沟的氮量以及坡面沟中氮素转化量，利用对流弥散方程计算坡面沟中氮素浓度；溶解态有机氮、氨氮和硝态氮在坡面沟中的迁移转化过程使用对流弥散方程进行描述，除对流扩散过程之外还包括旁侧入流、一个等高带的坡面沟到下一个等高带的坡面沟以及不同组分之间的相互转化作用，具体计算公式为，

其中，c′dis为河道中溶解态氮的浓度，kg/m³；v为水流的流速，m/d；ex为扩散系数，m²/d；s为有机氮、氨氮或硝态氮的源汇项，kg/m³/d，用于表示旁侧入流以及不同组分之间的相互转化作用，具体指代哪一种氮素的源汇项，需要根据计算的氮素种类进行确定。

若坡面沟窄而直，则ex＝5.93·d·v，d表示河流深度，m；

若坡面沟宽大，则ex＝0.11·v·w²/d，w表示坡面沟宽度，m。

对于有机氮，矿化作用会使有机氮转化为氨氮从而减少，有机氮在坡面沟中或河道中的源汇项可表示为，

其中，m1为有机氮的旁侧入流，kg；dx为河流的长度，m；a为河流的横截面积，m²；k1为有机氮的矿化速率，d^-1；c1为河流中有机氮的浓度，kg/m³；

对于氨氮，矿化作用会使有机氮转化为氨氮从而使河流中氨氮含量增加，硝化作用会使氨氮转化为硝态氮从而使河流中氨氮含量降低；则氨氮在坡面沟中或河道中的源汇项可表示为：

其中，m2表示氨氮的旁侧入流，kg；k2表示氨氮的硝化速率，d^-1；c2表示河流中氨氮的浓度，kg/m³；

对于硝态氮，硝化作用会使氨氮转化为硝态氮从而使河流中硝态氮含量增大，反硝化作用可使硝态氮转化为n2o，使河流中硝态氮含量降低；则硝态氮在坡面沟中或河道中的源汇项可表示为：

式中，m3表示硝态氮的旁侧入流，kg；k3表示硝态氮的反硝化速率，d^-1；c3表示河流中硝态氮的浓度，kg/m³。

五、河道中氮迁移转化模拟

该部分内容对应步骤s5的内容，如图6所示，s5具体为，根据河道水文过程，统计坡面沟流入河道的量以及地下水与河道交换量，计算河道中氮素转化量，利用对流弥散方程计算河道断面氮素浓度；其计算方法与步骤s4的计算方式相同。溶解态有机氮、氨氮和硝态氮在坡面沟和河道中的迁移转化过程均使用对流弥散方程进行描述，除对流扩散过程之外还包括旁侧入流、上一个子流域的河道到下一个子流域的河道以及不同组分之间的相互转化作用。

若河道窄而直，则ex＝5.93·d·v，d表示河流深度，m；

若河道宽大，则ex＝0.11·v·w²/d，w表示河道宽度，m。

河道中有机氮、氨氮以及硝态氮的旁侧入流包括地下水与河道的交换量、该子流域各等高带的坡面沟汇入河道的氮量。

通过采用本发明公开的上述技术方案，得到了如下有益的效果：

本发明提供了一种大流域尺度水氮迁移耦合模拟方法，考虑了地表、土壤、地下水、坡面沟与河道中氮素的迁移转化过程，溶解态氮随水移动，考虑了随各个阶段氮素的迁移。使用对流弥散方程计算冻融过程中土壤、坡面沟、河道污染物的迁移过程，提高了水文过程模拟的准确性，并且可以很好的模拟冻土融化期河道断面流量和污染物氮素的浓度。

以上所述仅是本发明的优选实施方式，应当指出，对于本技术领域的普通技术人员来说，在不脱离本发明原理的前提下，还可以做出若干改进和润饰，这些改进和润饰也应视本发明的保护范围。