基于多源互济供水模式的分布式水循环模型的设计方法

1.本发明涉及水文模型的设计技术领域，具体是一种基于多源互济供水模式的分布式水循环模型的设计方法。


背景技术：

2.随着人类生产力的提高和技术的进步，人类活动对水循环的影响明显加大，流域水循环演变呈现出复杂的自然
‑
人工复合特性。传统的水文模型多是以产汇流机制为核心开发的，一般以地表径流模拟为主，对资源开发利用的处理较为粗略，且不能对人工取用水活动及调控决策进行考虑，难以适应这种基于强人类扰动过程的精细化模拟需求。众多来自水文及水文相关学科的专家学者开始尝试在传统的水文模型框架上寻求新的突破，研究人类扰动下的水循环演变规律，分布式自然
‑
人工水循环模型应运而生。
3.将分布式水文模型与集总式水资源配置模型耦合是构建自然
‑
人工分布式水循环模型的最常用方式之一。分布式水文模型以精确模拟水循环过程为目标，能够补充水资源配置模型无法给出的水循环过程，提供其所需要的实时水资源边界情况；水资源配置模型能够处理人类调控下的水资源开发利用过程，给出取、用、耗、排水数据，有效提高水文模型的模拟精度。二者的耦合可以取长补短，充分发挥各自优势，实现自然
‑
人工水循环精确模拟。如文献一（赵勇.广义水资源合理配置研究[d]. 中国水利水电科学研究院, 2006.）中，将平原区分布式水循环模型(pdwcm)与广义水资源合理配置模型(gwram)耦合，采用分解和聚合的方式进行信息交互，实现了自然
‑
人工水循环全过程模拟；在文献二（张洪刚,熊莹,邴建平,李明新.nam模型与水资源配置模型耦合研究[j].人民长江,2008(17):15
‑
17.）中，利用水文模型nam与水资源配置模型的耦合，研究强人类活动影响下的汉江流域水资源响应；在文献三（张守平. 基于流域水循环及其伴生过程的水量水质联合配置研究[d].中国水利水电科学研究院,2015.）中，采用集总式水资源配置模型模拟人工“取水
‑
用水
‑
耗水
‑
排水”过程，将以三级区套地市为单位的水资源配置分区数据传输到以子流域高程带划分的计算单元上，实现了与分布式水循环模型wep
‑
l的耦合，进行了渭河流域自然
‑
人工水循环及伴生过程模拟。
[0004]
但是，这些模型的构建多以分布式水文模型与集总式水资源配置模型松散耦合的方式进行，也即是采用公用参数或者输出文件单向传递的方式实现模型的耦合，忽略自然水循环与社会水循环过程之间的动态互馈作用，对人类活动用水过程与天然水文过程之间的相互影响反映不足，不利于水资源的精细化模拟与调控。


技术实现要素：

[0005]
针对当前分布式自然
‑
人工水循环模型在自然
‑
人工水循环耦合模拟及水资源开发利用与调控方面的不足，本发明基于自然
‑
人工水循环理论，以swat模型为工具，开发一种分布式自然
‑
人工水循环动态互馈模拟模型，通过添加多源互济供水模块，使swat模型具有了自然
‑
人工水循环之间的动态互馈模拟、水资源开发利用与调控等功能，在模型运行
中，时刻保持自然水循环与人工侧枝水循环之间的动态互馈关系，不仅反映水循环过程对人工取用水的影响，同时还反映经济社会活动对水循环过程的实时干预作用。
[0006]
本发明采用如下技术方案实现：一种基于多源互济供水模式的分布式水循环模型的设计方法，包括如下步骤：步骤s1：嵌套式hru划分：采用基于“流域
‑
水资源区
‑
行政区
‑
灌区
‑
土地利用
‑
土壤
‑
坡度”多元属性的嵌套式坡面离散方法进行hru划分，划分后的hru具有相应的属性；步骤s2：hru属性设计：构建hru属性识别模块，所述hru属性识别模块用于识别hru的属性；步骤s3：设计多源互济供水模块，所述多源互济供水模块用于调用hru属性识别模块，以识别hru的属性，根据识别的属性判定土地利用类型、对应的供水源及供水源的供水优先次序，根据供水源的供水优先次序调用相应的水源模块进行取水；步骤s4：swat模型改进：将多源互济供水模块与swat模型中修改的模块进行衔接以实现数据实时交互，所述hru根据输入的需水信息、供水水源类型信息、供用水优先序规则信息、水利工程信息和由swat模型提供的天然水文条件信息进行逐日水资源分配与调控，输出逐日“供水
‑
用水
‑
耗水
‑
排水”的人工侧枝水循环信息，将人工侧枝水循环信息传递给swat模型，swat模型根据人工侧枝水循环信息模拟刻画天然水循环过程，实时模拟水文变化对水资源开发利用的影响，以及模拟人工用水对下一阶段的水资源和供水的影响，为人工侧枝水循环提供实时天然水文边界条件信息，实现“自然
‑
人工”水循环动态互馈式模拟。
[0007]
进一步的，所述步骤s1包括：（1）采用arcgis对dem提取流域河网，划分天然子流域；（2）在天然子流域上叠加土地利用信息、土壤类型信息和坡度信息，划分出天然hru；（3）将天然hru套上行政区、水资源区边界，将天然hru进一步切割；（4）根据灌区的分布，将灌区与天然hru进行叠加，最终完成hru划分。这样，每一个hru都具有子流域属性、水资源区属性、行政区属性、灌区属性、土地利用类型属性和土壤类型属性。
[0008]
进一步的，所述步骤s2包括：（1）hru属性识别模块构建构建hru属性识别模块，用于读取指定的hru的属性，所述指定的hru的属性包括子流域属性、水资源区属性、行政区属性、灌区属性；（2）hru属性识别模块的调用将构建的hru属性识别模块放在swat模型中的主模块中，便于对hru属性识别模块的调用。
[0009]
进一步的，所述步骤s3包括：（1）水源代码信息文件设计：所述水源代码信息文件用于读取指定的供水源信息，设置的供水源类型包括外调水、水库水、河道水、浅层地下水、深层地下水、坑塘水6种水源类型，水源代码信息文件通过程序指令读入；
（2）供水优先次序信息文件设计：所述供水优先次序信息文件用于读取水源供水优先次序信息，指定供水源供水优先次序，供水优先次序信息文件通过程序指令读入；（3）取水控制量信息文件设计：所述取水控制量信息文件用于读取地表水供水控制总量和地下水开采控制总量信息，识别hru所属行政区的年度地表水供水控制总量和地下水开采控制总量，用于后续的水源的取水量计算。
[0010]
（4）多源互济供水计算设计，计算流程如下：1）hru类型识别首先识别hru的土地利用类型，如果是居工地，则程序进入城市供水计算流程；如果是灌溉用地，则程序进入灌溉用水计算流程；如果是其他土地利用类型，程序结束；2）供水源及其优先序识别通过识别hru的取水源识别代码，调用相应的水源模块，通过水源代码信息文件和供水优先次序文件的读取，识别每一个hru的供水源个数、类型、取水先后顺序，根据取水源识别代码依次调用各水源模块；3）多源互济供水计算根据hru的供水次序寻求各水源进行取水，直至满足该hru的当日的生活、工业用水需求、农业灌溉用水需求、或者到最后一个供水源供水完毕为止。
[0011]
进一步的，所述多源互济供水计算的步骤包括：
①
明确目标hru设定的当日需水量wd；
②
明确目标hru的供水水源个数k、水源代码及其供水优先次序，其中，k≤30；
③
按次序调用水源模块，计算水源的取水量，水源模块包括河道取水模块（rchuse模块）、水库取水模块（res模块）、watuse模块（包含地下水、坑塘水和外调水的取水计算），水源取水量取决于当日hru的需水量和水源的可供水量，可供水量则取决于水源的可利用水量、取水工程（如引水渠道、供水管道、电井等）供水能力和取水控制量，计算公式如下：式中，i表示hru序号；j表示水源的供水优先次序号；wsp表示水源的当日实际取水量，m3；wd表示hru的当日需水量，m3；wsc表示水源的当日可供水量，m3；wf表示取水工程的供水能力，m3；wa表示水源的当日可利用水量，m3；wmx表示年度取水控制量，m3，wum表示年用水控制总量，m3；wsm表示年度地表取水控制总量，m3；wgm表示年度地下水开采控制总量，m3；
④
对供水优先次序为1的水源，则优先从该水源取水，如果水源的可供水量wsc1＞wd，那么该水源供水量为wsp1=wd，供水程序结束，hru的水源总供水量为wsp= wsp1；否则，
wsp1=wsc1，并且hru的需水量变为wf =wd
‑
wsc1，程序将继续寻找下一级供水源；
⑤
对供水优先次序为j（j=2,
…
,k
‑
1；k≤30）的水源，则优先从该水源取水。如果该水源的当日可供水量wsc
j
＞wf，那么该水源供水量为wsp
j
=wf，程序结束，hru的水源总供水量为wsp=wsp+wsp
j
；否则，wsp
j
=wsc
j
，并且hru的需水量变成wf=wf
‑
wsc
j
，程序将继续寻找下一级供水源；
⑥
对供水优先次序为k（k≤30）的水源，如果该水源的当日可供水量wsc
k
＞wf，那么该水源供水量为wsp
k
=wf，程序结束，hru的水源总供水量为wsp=∑wsp
i
；否则，wsp
k
=wsc
k
，并且hru的需水量变成wf=wf
‑
wsc
k
，程序结束。
[0012]
进一步的，所述步骤s4中修改swat模型具体包括（1）模块屏蔽屏蔽各水源模块，即rchuse 模块、res模块、watuse模块、irr_rch 模块、irr_res模块和irrsub模块，不再采取单一水源取水模式；将以上模块放在多源互济供水模块中进行调用；（2）修改rchuse 模块、res模块分别添加代码，用参数wspi（i=1,2）分别替换掉rchuse 模块、res模块中的参数waterrch、wuresn，实现多源互济供水模块multi_sc与rchuse 模块、res模块的衔接和调用；（3）修改watuse模块
①
由于该模块中没有外调水供水功能，通过修改程序增加外调水取水以及外调水量限制约束，供水量不能超过调水限制总量，计算公式表达如下：其中，waterout (i,j)表示第i天、第j个hru的调水利用量，m3；wx5表示调水限制总量，m3；
②
在watuse模块中添加代码，用参数wspi（i=3,4,5,6）分别替换掉watuse模块中的参数watershal、waterdeep、waterout和waterpnd，实现多源互济供水模块multi_sc与watuse模块的衔接与调用；（4）添加管网输水模拟功能城市供水管网在向用水户输送水量的过程中，会因管网材质老化或破损因素产生“跑、冒、滴、漏”输水损失，因此有必要考虑城市供水损失，在rchuse 模块、res模块、watuse模块中添加如下公式的计算程序：式中，pip表示供水管网漏损率；（5）修改irrsub模块由于模块中没有坑塘灌溉，添加坑塘灌溉功能；完善外调水灌溉功能，通过增加供水量限制约束，灌溉取水量不能超过调水限制总量：
其中，wirrout (i,j)表示第i天、第j个hru的调水灌溉利用量，m3；wx5表示调水限制总量，m3，wirrpnt (i,j)表示第i天、第j个hru的坑塘灌溉利用量，m3；wx6表示坑塘可供水量，m3；（6）添加灌溉渠道模拟功能
①
swat模型并不考虑灌溉渠道的输水损失对水循环的影响，灌溉渗漏视为系统损失，通过修改irr_rch模块、irr_res模块、irrsub模块的源代码，增加渠系输水过程模拟，包括渠道水量损失和渠道退水，其中，渠道水量损失包括渠道内的水分蒸发损失和渠道渗漏损失两部分，主要计算公式如下：式中，et
can
表示渠系蒸发损失量，mm；irr
can
表示进入渠道的灌溉水量，mm；ls
can
表示渠系渗漏损失量，mm；surp
can
表示渠系退水量，mm；表示渠系水有效利用系数；表示渠系蒸发系数；表示渠系渗漏系数；
②
渗漏损失水量进入上层土壤补给土壤水，参与土壤水分循环，通过修改percmain模块相关代码，添加渗漏损失的计算程序，计算公式如下：式中，表示第t天第一层土壤的土壤含水量，mm；表示第t
‑
1天第一层土壤的土壤含水量，mm；inf
pcp
表示降雨入渗量，mm；inf
irr
表示灌溉入渗量，mm；inf
wet
表示湖库湿地渗漏量，mm；（7）修改gwmod模块代码供水管网漏失水量进入地下含水层，补给潜水，通过修改gwmod模块中的地下水补给代码，实现对管网漏损水量的模拟，计算公式如下：式中，表示第t天的潜水补给量，mm；表示第t
‑
1天的潜水补给量，mm；prc表示补给地下水的土壤渗漏水量，mm；gw_delay表示地下水补给延迟系数，mm；area表示hru的面积，m2；（8）subbasin模块修改在subbasin模块中添加多源互济供水模块，便于swat模型运行时的调用，通过读取指定的水源类型、水源个数、水源代码、取水量、取水时间，进行年内的水源动态互补供水运算，实现swat模型运行时的多水源联合供水模拟；（9）surface模块修改渠系退水是指渠道排泄超引灌溉水，渠道退水直接进入河道，参与河道汇流计算，
通过修改surface模块相关代码，将叠加地表产流量，计算公式如下：式中，表示渠道退水后的径流量，mm；表示渠道退水前的径流量，mm；（10）点源模块修改城市生活和工业污水产生后，通过排水管网系统直接排入河流，或者输送至污水处理厂进行处理，处理达标后的污水一部分直接排放进入河道，一部分通过深度处理达到再生水利用标准后回用于绿化、生活杂用以及生产；swat模型采用其点源模块进行模拟，所述点源模块包括recday模块和recmon模块，在recday模块和recmon模块中修改相关代码，将污染排放参数wdr分别替换参数floday、flomon，相关计算公式如下：其中，wdr为城市污水产生量，m3；wp为污水排放量，m3；r表示用水消耗率；v表示污水处理厂的污水处理率；re表示再生水利用率；（11）main主模块修改将构建的readattr模块放在swat模型中的main主模块中，便于对readattr模块的调用。
[0013]
与现有技术相比，本发明具有如下有益效果：（1）模型采用嵌套式坡面离散方法，很好地处理了天然hru边界、行政区边界、灌区边界互不重合的问题，既能体现传统分布式水循环模型单元划分特点，也满足了水资源区管理和行政区管理相结合的需求；（2）嵌入的多源互济供水模块实现社会水循环的多水源、多行业联合配置，具有刻画多水源（包括河道、水库、地下水、外调水、坑塘水等）、多工程（包括蓄水工程、引水工程、提水工程、调水工程等）的功能，将水资源系统中的各种水源与用水户之间的拓扑关系、水分在降雨
‑
产流过程以及各社会生产部门中的传递关系进行客观、清晰的描述，真实体现人类活动对水循环过程的影响；（3）通过添加多源互济供水模块和修改swat模型的相关模块，构建了分布式自然
‑
人工水循环动态互馈模拟模型，弥补了以往的水循环模型不能很好地反映人类活动用水过程与天然水文过程之间的相互影响的不足，并且模型在进行长序列模拟时可以根据生活生产取用水的需要以及水源的水量情况，灵活取用不同的水源，更加适用于强人类活动地区。该模型具有了自然
‑
人工水循环之间的动态互馈模拟、水资源开发利用与调控等功能，实现了自然水循环与人工侧枝水循环之间的双向耦合模拟，充分反映复杂水循环系统的动态性、互馈性特征，可以作为区域自然
‑
人工水循环互馈模拟以及水资源系统精细化管理的有力支撑工具。
附图说明
[0014]
图1为本发明实施例基于多源互济供水模式的分布式水循环模型的设计方法的流程图；图2为本发明实施例hru划分流程图；
图3为本发明实施例多源互济供水模块的运行流程图；图4为本发明实施例水源供水优先计算流程图；图5为本发明实施例新店铺站的实测与模拟月径流过程对比；图6为2016水平年白河流域自然
‑
人工水循环转化关系图；图7为2016水平年白河流域各水源年内供水过程图。
具体实施方式
[0015]
为使本发明实施例的目的、技术方案和优点更加清楚，下面将结合本发明实施例中的附图，对本发明实施例中的技术方案进行清楚、完整地描述，显然，所描述的实施例是本发明的一部分实施例，而不是全部的实施例。基于本发明中的实施例，本领域普通技术人员在没有做出创造性劳动的前提下所获得的所有其他实施例，都属于本发明保护的范围。
[0016]
本发明实施例提供一种基于多源互济供水模式的分布式水循环模型的设计方法，主要步骤包括嵌套式hru划分、hru属性设计、多源互济供水模块设计和swat模型改进四部分内容。多源互济供水模块设计就是设计并添加相关模块，增加swat模型的城乡多水源供水和多水源灌溉模拟功能，完善人工侧枝水循环，进而构建自然
‑
人工分布式水循环模型。swat模型改进，即是修改swat模型相关模块，便于水源类型、水源个数、水源代码、取水优先序信息的读取以及多源互济模块的调用。所涉及的模块主要包括rchuse模块、res模块、watuse模块、irr_rch模块、irr_res模块、irrsub模块、subbasin模块、percmain模块、gwmod模块、surface模块、recday模块、recmon模块、main模块等。
[0017]
本发明所述方法包括以下步骤：步骤s1：嵌套式hru划分；步骤s2：hru属性设计；步骤s3：多源互济供水模块设计；步骤s4：swat模型改进。
[0018]
在所述的步骤s1中，采取基于“流域
‑
水资源区
‑
行政区
‑
灌区
‑
土地利用
‑
土壤
‑
坡度”多元属性的嵌套式坡面离散方法进行hru划分。其主要步骤（图2）主要包括：（1）天然子流域划分；（2）天然hru划分；（3）天然hru叠加管理单元；（4）天然hru叠加灌区。
[0019]
（1）采用arcgis中的 arc hydro tools 工具进行天然子流域划分。首先对 dem 底图进行填洼等修正操作，依据设定的生成子流域临界集水面积csa（critical source area）阈值要求，确定格网水流流向，识别流域分水线，并依据流域坡度、坡向、坡长等河网特征参数划分天然子流域；（2）按照子流域内的土地利用类型、土壤类型和坡度类型进行叠加，划分出若干个天然hru。 hru是指下垫面特征相对单一和均匀的区域，在这个区域中的下垫面具有相同的水文特性，每个hru内只有一种植被、一种土壤和一种坡度类型的组合。
[0020]
（3）采用选取hru叠加地市/县行政区和水资源区的嵌套方式。采用gis的“editor toolbar”功能，将管理单元边界gis图直接叠加到上一步所划分的天然hru图上，在天然hru与管理区划交界地区，依照边界线将天然hru一分为二，并将划分的hru附上管理单元属性；（4）采用gis的“editor toolbar”功能，在天然hru与灌区交界地区，依照边界线将天然hru一分为二，并将划分的hru附上灌区属性，最终完成hru的划分。最后，每一个hru都
具有子流域属性、水资源区属性、行政区属性、灌区属性、土地利用类型属性和土壤类型属性。
[0021]
在所述步骤s2中，hru属性的设计包括如下内容：（1）hru属性识别模块readattr构建：构建hru属性识别模块readattr，用于读取指定的hru属性。设置的hru属性分别是子流域属性、水资源区属性、行政区属性、灌区属性。输入文件格式见表1。
[0022]
（2）hru属性识别模块的调用将构建的readattr模块（hru属性识别模块）放在swat模型中的main主模块中，便于对readattr模块的调用。
[0023]
在所述步骤s3中，多源互济供水模块multi_sc的设计包括如下内容：（1）水源代码信息文件设计；（2）供水优先次序信息文件设计；（3）多源互济供水计算设计。
[0024]
（1）水源代码信息文件设计：用于读取指定的供水源信息。设置的供水源类型包括外调水、水库水、河道水、浅层地下水、深层地下水、坑塘水6种水源类型，并且设定单个hru的同类型水源最多为5个，即每个hru最多可以设置30个供水源。水源代码以输入文件读入的方式进行指定，并针对不同的hru类型，居工地hru以行政区为单位进行指定，同一行政区内的hru供水源一致；灌溉用地hru以灌区为单位进行指定，同一灌区内的hru供水源一致。输入文件格式分别见表2和表3。
[0025]
（2）供水优先次序信息文件设计：用于读取指定的水源供水优先次序信息。每一个水源可以为多个hru利用，每一个hru也可以从多个供水源取水。根据每一个hru的供水源类型和个数，设定供水次序。指定供水源的供水优先次序，供水优先次序以输入文件读入的方式进行指定，并针对不同的hru类型，居工地hru以行政区为单位进行指定，同一行政区内的hru水源供水顺序一致；灌溉用地hru以灌区为单位进行指定，同一灌区内的hru水源供水顺序一致。优先次序以1到k（k为hru的供水水源个数，k≤30）进行排序，1为最高优先序，对未指定的水源，序号设置为0。输入文件格式见表4和表5。
[0026]
（3）取水控制量信息文件设计：用于读取各行政区各年度的供水控制总量信息。每个hru根据自身的行政区属性，识别所属行政区的年度地表水供水控制总量和地下水开采控制总量，用于后续的水源的取水量计算。输入文件格式见表6。
[0027]
（4）多源互济供水计算设计：多源互济供水模块multi_sc的运行流程（图3）主要包括：1）hru类型识别；2）供水源及其优先序识别；3）多源互济供水计算。
[0028]
1）hru类型识别首先识别hru的土地利用类型。如果是建设用地（城镇用地、农村用地），则程序进入城乡供水计算流程；如果是灌溉用地，则程序进入灌溉用水计算流程；如果是其他土地利用类型，则该hru无供水要求，程序退出。
[0029]
2）供水源及其优先序识别多源互济供水模块通过识别hru的取水源识别代码，调用相应的水源模块，河道、
水库、浅层地下水、深层地下水、外调水和坑塘水的取水源识别代码分别是1、2、3、4、5、6。通过水源代码信息文件和供水优先次序文件的读取，识别每一个hru的供水源个数、类型、取水先后顺序，根据取水源识别代码依次调用各水源模块。
[0030]
3）多源互济供水计算多源互济供水模块根据hru的供水次序寻求各水源进行取水，直至满足该hru的当日的生活、工业用水需求（建设用地hru）或者农业灌溉用水需求（农业用地hru）或者到最后一个供水源供水完毕为止。供水优先计算的步骤（图4）如下：
①
明确目标hru设定的当日需水量wd；
②
明确目标hru的供水水源个数k（k≤30）、水源代码及其供水优先次序；
③
按次序调用水源模块，计算水源的取水量。水源模块包括rchuse模块（河道取水模块）、res模块（水库取水模块）、watuse模块（地下水、坑塘水的取水计算）。水源取水量取决于当日hru的需水量和水源的可供水量。可供水量则取决于水源的可利用水量、取水工程（如引水渠道、供水管道、电井等）供水能力和取水控制量。计算公式如下：式中，i表示hru序号；j表示水源的供水优先次序号；wsp表示水源的当日实际取水量，m3；wd表示hru的当日需水量，m3；wsc表示水源的当日可供水量，m3；wf表示取水工程的供水能力，m3；wa表示水源的当日可利用水量，m3；wmx表示年度取水控制量，m3。wum表示年用水控制总量，m3；wsm表示年度地表取水控制总量，m3；wgm表示年度地下水开采控制总量，m3。
[0031]
④
对供水优先次序为1的水源，则优先从该水源取水。如果水源的可供水量wsc1＞wd，那么该水源供水量为wsp1=wd，供水程序结束，hru的水源总供水量为wsp= wsp1；否则，wsp1=wsc1，并且hru的需水量变为wf =wd
‑
wsc1，程序将继续寻找下一级供水源；
⑤
对供水优先次序为j（j=2,
…
,k
‑
1；k≤30）的水源，如该水源的当日可供水量wsc
j
＞wf，那么该水源供水量为wsp
j
=wf，程序结束，hru的水源总供水量为wsp=wsp+wsp
j
；否则，wsp
j
=wsc
j
，并且hru的需水量变成wf=wf
‑
wsc
j
，程序将继续寻找下一级供水源；
⑥
对供水优先次序为k（k≤30）的水源，如果该水源的当日可供水量wsc
k
＞wf，那么该水源供水量为wsp
k
=wf，程序结束，hru的水源总供水量为wsp=∑wsp
i
；否则，wsp
k
=wsc
k
，并且hru的需水量变成wf=wf
‑
wsc
k
，程序结束。
[0032]
在所述步骤s4中，swat模型改进主要是将多源互济供水模块与swat模型中修改的模块进行衔接以实现数据实时交互，所述hru根据输入的需水信息、供水水源类型信息、供用水优先序规则信息、水利工程信息和由swat模型提供的天然水文条件信息进行逐日水资源分配与调控，输出逐日“供水
‑
用水
‑
耗水
‑
排水”的人工侧枝水循环信息，将人工侧枝水循环信息传递给swat模型。其中修改swat模型中的相关模块，便于水源类型、水源个数、水源
代码、取水优先序信息的读取以及多源互济模块的调用。修改的模块主要有rchuse模块、res模块、watuse模块、irr_rch模块、irr_res模块、irrsub模块、subbasin模块、surface模块、gwmod模块、recday模块、recmon模块和main模块等。各模块主要功能如下：rchuse 模块：河道供水模块，从指定的河道中取水，用于生活工业用水。
[0033]
res模块：水库供水模块，从指定的水库中取水，用于生活工业用水。
[0034]
watuse模块：地下水、坑塘供水模块，从指定的子流域内的浅层地下含水层、深层地下含水层和坑塘中取水，用于生活工业用水。
[0035]
irr_rch 模块：河道灌溉模块，从指定的河道中取水，用于农田灌溉用水。
[0036]
irr_res模块：水库灌溉模块，从指定的水库中取水，用于农田灌溉用水。
[0037]
irrsub模块：地下水、外调水灌溉模块，从指定的子流域内的浅层地下含水层、深层地下含水层和外源水源中取水，用于农田灌溉。
[0038]
percmain模块：渗漏损失模块，根据降雨入渗、灌溉入渗、坑塘入渗计算结果，模拟土壤水变化。
[0039]
gwmod模块：地下水模块，根据地表入渗结果，模拟浅层、深层地下水变化以及地下径流。
[0040]
subbasin模块：子流域模块，主要包括蒸发模拟、产流模拟、地下水模拟、植被生长模拟、农田管理模拟、营养物质迁移转化模拟等内容，通过对水文、水质过程相关模块的调用实现。
[0041]
surface模块：地表径流模拟模块，模拟冠层截流、积雪融雪、地表产流等主要地表水文过程。
[0042]
recday模块：日点源模块，通过输入逐日污染物排放量，模拟逐日点源变化过程。
[0043]
recmon模块：月点源模拟，通过输入逐月污染物排放量，模拟逐月点源变化过程。
[0044]
main模块：swat模型主模块，主要包括读入输入文件，分配数组大小、初始化参数、水文过程模拟等内容。
[0045]
修改swat模型中的相关模块，具体包括：（1）模块屏蔽屏蔽水源模块，即rchuse 模块、res模块、watuse模块、irr_rch 模块、irr_res模块和irrsub模块，不再采取单一水源取水模式。将以上模块放在多源互济供水模块中进行调用。
[0046]
（2）修改rchuse 模块、res模块添加代码，将参数wsp
i
（i=1,2）分别替换掉rchuse 模块、res模块中的参数waterrch、wuresn，实现多源互济供水模块multi_sc与rchuse 模块、res模块的衔接和调用；（3）修改watuse模块
①
由于该模块中没有外调水供水功能，通过修改程序增加外调水取水以及外调水量限制约束，供水量不能超过调水限制总量，计算公式表达如下：
其中，waterout (i,j)表示第i天、第j个hru的调水利用量，m3；wx5表示调水限制总量，m3。
[0047]
②
在watuse模块中添加代码，用参数wsp
i
（i=3,4,5,6）分别替换掉watuse模块中的参数watershal、waterdeep、waterout和waterpnd，实现多源互济供水模块multi_sc与watuse模块的衔接与调用。
[0048]
（4）添加管网输水模拟功能城市供水管网在向用水户输送水量的过程中，会因管网材质老化或破损等因素产生“跑、冒、滴、漏”等输水损失。因此有必要考虑城市供水损失。在rchuse 模块、res模块、watuse模块中添加如下公式的计算程序代码：式中，pip表示供水管网漏损率。
[0049]
（5）修改irrsub模块由于模块中没有坑塘灌溉，添加坑塘灌溉功能；完善外调水灌溉功能，通过增加供水量限制约束，灌溉取水量不能超过调水限制总量。
[0050]
其中，wirrout (i,j)表示第i天、第j个hru的调水灌溉利用量，m3；wx5表示调水限制总量，m3。wirrpnt (i,j)表示第i天、第j个hru的坑塘灌溉利用量，m3；wx6表示坑塘可供水量，m3。
[0051]
（6）添加灌溉渠道模拟功能
①
swat模型并不考虑灌溉渠道的输水损失对水循环的影响，灌溉渗漏视为系统损失。通过修改irr_rch模块、irr_res模块、irrsub模块的源代码，增加渠系输水过程模拟，包括渠道水量损失和渠道退水。其中，渠道水量损失包括渠道内的水分蒸发损失和渠道渗漏损失两部分。主要计算公式如下：式中，et
can
表示渠系蒸发损失量，mm；irr
can
表示进入渠道的灌溉水量，mm；ls
can
表示渠系渗漏损失量，mm；surp
can
表示渠系退水量，mm；表示渠系水有效利用系数；表示渠系蒸发系数；表示渠系渗漏系数。
[0052]
②
渗漏损失水量进入上层土壤补给土壤水，参与土壤水分循环。通过修改percmain模块相关代码，添加渗漏损失计算程序。计算公式如下：式中，表示第t天第一层土壤的土壤含水量，mm；表示第t
‑
1天第一层土壤的土壤含水量，mm；inf
pcp
表示降雨入渗量，mm；inf
irr
表示灌溉入渗量，mm；inf
wet
表示湖库湿地渗漏量，mm；
（7）修改gwmod模块代码供水管网漏失水量进入地下含水层，补给潜水。通过修改gwmod模块中的地下水补给代码，实现对管网漏损水量的模拟。计算公式如下：式中，表示第t天的潜水补给量，mm；表示第t
‑
1天的潜水补给量，mm；prc表示补给地下水的土壤渗漏水量，mm；gw_delay表示地下水补给延迟系数，mm；area表示hru的面积，m2；（8）subbasin模块修改在subbasin模块中添加多源互济供水模块，便于swat模型运行时的调用。通过读取指定的水源类型、水源个数、水源代码、取水量、取水时间，进行年内的水源动态互补供水运算，实现swat模型运行时的多水源联合供水模拟。
[0053]
（9）surface模块修改渠系退水是指渠道排泄超引灌溉水，渠道退水直接进入河道，参与河道汇流计算。通过修改surface模块相关代码，将叠加地表产流量。计算公式如下：式中，表示渠道退水后的径流量，mm；表示渠道退水前的径流量，mm；（10）点源模块修改城市生活和工业污水产生后，通过排水管网系统直接排入河流，或者输送至污水处理厂进行处理。处理达标后的污水一部分直接排放进入河道，一部分通过深度处理达到再生水利用标准后回用于绿化、生活杂用以及生产等用途。swat模型采用其点源模块（recday模块和recmon模块）进行模拟。在recday模块和recmon模块中修改相关代码，将污染排放参数wdr分别替换参数floday、flomon。相关计算公式如下：其中，wdr为城市污水产生量，m3；wp为污水排放量，m3；r表示用水消耗率；v表示污水处理厂的污水处理率；re表示再生水利用率。
[0054]
（11）main主模块修改将构建的readattr模块放在swat模型中的main主模块中，便于对readattr模块的调用。
[0055]
用fortran语言在windows平台应用程序开发环境visual studio 2012中根据以上结构对swat源代码进行修改，实现了对多源互济供水模块的调用，具有了城市多水源供水和多水源农业灌溉模拟功能，从而构建了以swat模型为基础的、具有多源互济功能的分布式自然
‑
人工水循环双向互馈模拟模型。
[0056]
本发明选取白河流域作为实施案例，白河流域自然
‑
人工水循环模拟过程介绍如下：1、研究区概况白河流域地处汉江流域中部，北起伏牛山地区，西与老鹳河为邻，东靠唐河流域，
南抵汉江，大部分地区位于南阳盆地内，流域总面积为1.23万km2。地形呈北高南低态势，水系主要由白河及其支流湍河、刁河、严陵河等，多年平均水资源量20.8亿m3。流域内有河南省南阳市卧龙区、宛城区、新野县、湖北省襄阳市襄州区等3个地级市、共12个县级行政单元。白河流域境内有鸭河口灌区、引丹灌区等大型灌区4个，赵湾水库灌区、高丘水库灌区、陡坡水库灌区等中型灌区18个。为防洪抗旱、保障农业灌溉用水，流域内建有大型水库（鸭河口水库）1座，赵湾水库、虎山水库、陡坡水库、官沟水库等中型水库十几座。
[0057]
2、基础数据收集模型构建所需要的数据有空间数据，包括dem（90m
×
90m），土地利用图（1:10万）、土壤分布图（1:100万）、行政区分布图、灌区分布图、以及水系图等；白河流域内及邻近的4个气象站的1990~2016年气象数据，包括日降水量、日最高及最低气温、相对湿度、日照时数风速等气象要素；水利工程数据，主要包括水库位置、死库容、调节库容、总库容等数据，灌溉渠道过水能力、电机井日提水能力等数据；新店铺水文站1991~2016年的逐月径流数据，用于模型的率定和验证；农业灌溉用水数据，通过查询2006~2016年南阳市水资源公报、襄阳市水资源公报等资料，获得历年农业灌溉用水数据；灌溉面积与种植结构数据，包括小麦、水稻、花生、油菜、芝麻、棉花、蔬菜、瓜果等各类作物种植面积。此外，还有播种、灌溉、施肥、收割等作物管理措施有关的时间、次数、单次（灌溉、施肥）数量等基本信息等。
[0058]
3、建模过程针对白河流域水资源分布、水利工程特点及各行业取用水情况，模型构建步骤如下：（1）hru划分：基于dem数据、土地利用类型图、土壤类型图划分天然子流域和hru，共获得34个天然子流域和507个天然hru；在此基础上，根据行政区划及灌区分布情况，将天然hru细分成1027个hru，每个hru均具有子流域属性、水资源区属性、行政区属性、灌区属性、土地利用类型属性和土壤类型属性。
[0059]
（2）农业种植与管理信息输入：包括各行政区县的种植作物类型、作物种植面积与灌溉面积、轮作制度与灌溉制度等。
[0060]
（3）供水规则信息输入：主要包括水源供水对象、分区供水原则、水源供水优先级、用水户分水原则、供水水源类型与个数、行业用水优先级等信息；此外，需要输入的信息还有水库、渠道、电机井等水利工程信息，渠系水有效利用系数、田间水有效利用系数等、耗水率等信息。
[0061]
（4）气象数据信息输入：选择建立好的降水、气温、风速、辐射和相对湿度数据输入，再把所有属性数据及水库数据输入，模型构建完成后即可开始运行。
[0062]
4、参数率定与模型校验通过上述过程和分析，确定模型的主要参数值，参数调整的最终值见表7所示。通过对模型参数的调整，模型的径流模拟与实测径流过程对比结果见表8和图5所示。在率定期（1995~2005年），新店铺站月径流模拟值与实测值的相关系数为0.792，纳什效率系数为0.756；在验证期（2006~2016年），新店铺站月径流模拟值与实测值的相关系数为0.643，纳什效率系数为0.635。可以看出，月径流模拟值与实测值流量过程线拟合程度较好，模型的模拟精度达到了要求值。
[0063]
对全部位于流域内的行政区（卧龙区、宛城区、镇平县）的模拟结果的偏差百分比
见表9所示。可以看出，2010~2016年三个县（区）的供用水量偏差很小，其生活用水量、工业用水量、农业用水量、地表水供水量以及地下水供水量与实际水量的误差均在10%以内。可见，模型的模拟结果很好地反映了白河流域实际用水和供水情况。
[0064]
5、模拟结果分析通过对白河流域2016年水循环模拟汇总，获得流域水资源循环转化关系见图6所示。2016年全流域总降水量95.47亿m3，地表产流量22.73亿m3，年末土壤水蓄变量8.62亿m3、地表水蓄变量1.25亿m3、地下水蓄变量为0.65亿m3。全流域耗水总量为70.61亿m3，其中土壤蒸发29.90亿m3，植被散发33.09亿m3，截留蒸发6.25亿m3，积雪升华0.003亿m3，水面蒸发1.17亿m3，生活工业消耗0.45亿m3。考虑外调水利用，全流域的经济社会总用水量为9.53亿m3，其中，地表水利用量6.21亿m3（含外调水4.47亿m3），地下水利用量3.32亿m3；人工总耗水量为5.34亿m3，人工排水量为1.71亿m3，总出境水量20.30亿m3。
[0065]
汇总2016年白河流域各水源供水情况如图7所示。2016年，白河流域水资源利用以地下水和外调水（丹江口水库供水）为主，全年地下水开发利用量占流域用水总量的34.83%（3.32亿m3），外调水利用量占流域用水总量的46.92%（4.47亿m3），河道水利用量最小，仅有0.18亿m3，占流域用水总量的1.95%。全年的用水峰值集中在3月和8月，3月份是冬小麦用水的关键的时期，外调水利用量为0.86亿m3，其次为蓄水供水（0.42亿m3），由于地表水资源量相对充足，地下水仅用于补充供水，补充供水量为0.26亿m3。8月份为玉米生长用水的关键的时期，外调水利用量为0.77亿m3，由于蓄水资源量不足，供水量仅为0.2亿m3，地下水起到了重要的补充供水作用，补充供水量为0.75亿m3，在数量上几乎与外调水供水量相等。此外，在10
‑
11月、1月出现用水小高峰，主要是由于种植冬小麦后的灌溉用水。由于这段时间处于枯水季节，地表水资源量有限，供水主要由丹江口水库供水（外调水）和地下水补充灌溉。
[0066]
表1
hru代码子流域属性水资源区属性行政区属性灌区属性字段格式整形整形整形整形说明子流域代码水资源区代码行政区代码灌区代码内容hru所属子流域代码hru所属水资源区代码hru所属的行政区代码hru所属的灌区代码
表2参数名称数据类型说明内容cntyid整形行政区代码行政区序号riv1整形1号河道水源河道所在的子流域代码riv2整形2号河道水源河道所在的子流域代码riv3整形3号河道水源河道所在的子流域代码riv4整形4号河道水源河道所在的子流域代码riv5整形5号河道水源河道所在的子流域代码res1整形1号水库水源水库代码res2整形2号水库水源水库代码res3整形3号水库水源水库代码res4整形4号水库水源水库代码res5整形5号水库水源水库代码
shal1整形1号浅水层水源浅水层所在的子流域代码shal2整形2号浅水层水源浅水层所在的子流域代码shal3整形3号浅层水水源浅水层所在的子流域代码shal4整形4号浅水层水源浅水层所在的子流域代码shal5整形5号浅水层水源浅水层所在的子流域代码deep1整形1号深水层水源深水层所在的子流域代码deep2整形2号深水层水源深水层所在的子流域代码deep3整形3号深层水水源深水层所在的子流域代码deep4整形4号深水层水源深水层所在的子流域代码deep5整形5号深水层水源深水层所在的子流域代码out1整形1号坑塘水源坑塘所在的子流域代码out2整形2号坑塘水源坑塘所在的子流域代码out3整形3号坑塘水源坑塘所在的子流域代码out4整形4号坑塘水源坑塘所在的子流域代码out5整形5号坑塘水源坑塘所在的子流域代码pnd1整形1号外流域水源外源水代码pnd2整形2号外流域水源外源水代码pnd3整形3号外流域水源外源水代码pnd4整形4号外流域水源外源水代码pnd5整形5号外流域水源外源水代码表2中若无供水，代码以0表示。
[0067]
表3参数名称数据类型说明内容irrid整形灌区代码灌区序号riv1整形1号河道水源河道所在的子流域代码riv2整形2号河道水源河道所在的子流域代码riv3整形3号河道水源河道所在的子流域代码riv4整形4号河道水源河道所在的子流域代码riv5整形5号河道水源河道所在的子流域代码res1整形1号水库水源水库代码res2整形2号水库水源水库代码res3整形3号水库水源水库代码res4整形4号水库水源水库代码res5整形5号水库水源水库代码shal1整形1号浅水层水源浅水层所在的子流域代码shal2整形2号浅水层水源浅水层所在的子流域代码shal3整形3号浅层水水源浅水层所在的子流域代码shal4整形4号浅水层水源浅水层所在的子流域代码shal5整形5号浅水层水源浅水层所在的子流域代码
deep1整形1号深水层水源深水层所在的子流域代码deep2整形2号深水层水源深水层所在的子流域代码deep3整形3号深层水水源深水层所在的子流域代码deep4整形4号深水层水源深水层所在的子流域代码deep5整形5号深水层水源深水层所在的子流域代码out1整形1号坑塘水源坑塘所在的子流域代码out2整形2号坑塘水源坑塘所在的子流域代码out3整形3号坑塘水源坑塘所在的子流域代码out4整形4号坑塘水源坑塘所在的子流域代码out5整形5号坑塘水源坑塘所在的子流域代码pnd1整形1号外流域水源外源水代码pnd2整形2号外流域水源外源水代码pnd3整形3号外流域水源外源水代码pnd4整形4号外流域水源外源水代码pnd5整形5号外流域水源外源水代码表3中若无供水，代码以0表示。
[0068]
表4参数名称数据类型说明内容cntyid整形行政区代码行政区序号sup_riv1整形1号河道水源的取水序号供水序号代码sup_riv2整形2号河道水源的取水序号供水序号代码sup_riv3整形3号河道水源的取水序号供水序号代码sup_riv4整形4号河道水源的取水序号供水序号代码sup_riv5整形5号河道水源的取水序号供水序号代码sup_res1整形1号水库水源的取水序号供水序号代码sup_res2整形2号水库水源的取水序号供水序号代码sup_res3整形3号水库水源的取水序号供水序号代码sup_res4整形4号水库水源的取水序号供水序号代码sup_res5整形5号水库水源的取水序号供水序号代码sup_shal1整形1号浅水层水源的取水序号供水序号代码sup_shal2整形2号浅水层水源的取水序号供水序号代码sup_shal3整形3号浅层水水源的取水序号供水序号代码sup_shal4整形4号浅水层水源的取水序号供水序号代码sup_shal5整形5号浅水层水源的取水序号供水序号代码sup_deep1整形1号深水层水源的取水序号供水序号代码sup_deep2整形2号深水层水源的取水序号供水序号代码sup_deep3整形3号深层水水源的取水序号供水序号代码sup_deep4整形4号深水层水源的取水序号供水序号代码sup_deep5整形5号深水层水源的取水序号供水序号代码
sup_out1整形1号坑塘水源的取水序号供水序号代码sup_out2整形2号坑塘水源的取水序号供水序号代码sup_out3整形3号坑塘水源的取水序号供水序号代码sup_out4整形4号坑塘水源的取水序号供水序号代码sup_out5整形5号坑塘水源的取水序号供水序号代码sup_pnd1整形1号外流域水源的取水序号供水序号代码sup_pnd2整形2号外流域水源的取水序号供水序号代码sup_pnd3整形3号外流域水源的取水序号供水序号代码sup_pnd4整形4号外流域水源的取水序号供水序号代码sup_pnd5整形5号外流域水源的取水序号供水序号代码表4中若无供水，代码以0表示。
[0069]
表5参数名称数据类型说明内容irrid整形灌区代码灌区序号irr_riv1整形1号河道水源的取水序号供水序号代码irr_riv2整形2号河道水源的取水序号供水序号代码irr_riv3整形3号河道水源的取水序号供水序号代码irr_riv4整形4号河道水源的取水序号供水序号代码irr_riv5整形5号河道水源的取水序号供水序号代码irr_res1整形1号水库水源的取水序号供水序号代码irr_res2整形2号水库水源的取水序号供水序号代码irr_res3整形3号水库水源的取水序号供水序号代码irr_res4整形4号水库水源的取水序号供水序号代码irr_res5整形5号水库水源的取水序号供水序号代码irr_shal1整形1号浅水层水源的取水序号供水序号代码irr_shal2整形2号浅水层水源的取水序号供水序号代码irr_shal3整形3号浅层水水源的取水序号供水序号代码irr_shal4整形4号浅水层水源的取水序号供水序号代码irr_shal5整形5号浅水层水源的取水序号供水序号代码irr_deep1整形1号深水层水源的取水序号供水序号代码irr_deep2整形2号深水层水源的取水序号供水序号代码irr_deep3整形3号深层水水源的取水序号供水序号代码irr_deep4整形4号深水层水源的取水序号供水序号代码irr_deep5整形5号深水层水源的取水序号供水序号代码irr_out1整形1号坑塘水源的取水序号供水序号代码irr_out2整形2号坑塘水源的取水序号供水序号代码irr_out3整形3号坑塘水源的取水序号供水序号代码irr_out4整形4号坑塘水源的取水序号供水序号代码irr_out5整形5号坑塘水源的取水序号供水序号代码
irr_pnt1整形1号外流域水源的取水序号供水序号代码irr_pnt2整形2号外流域水源的取水序号供水序号代码irr_pnt3整形3号外流域水源的取水序号供水序号代码irr_pnt4整形4号外流域水源的取水序号供水序号代码irr_pnt5整形5号外流域水源的取水序号供水序号代码表5中若无供水，代码以0表示。
[0070]
表6参数名称数据类型说明内容cnty_surf(1,1)实型1号行政区第1年的地表水取水量取水量cnty_surf(1,2)实型1号行政区第2年的地表水取水量取水量::::cnty_surf(m,n)实型m号行政区第n年的地表水取水量取水量cnty_gw(1,1)实型1号行政区第1年的地下水取水量取水量cnty_gw(1,2)实型1号行政区第2年的地下水取水量取水量::::cnty_gw(m,n)实型m号行政区第n年的地下水取水量取水量表7参数名称调参方式物理意义调参值cn2r湿润条件ii下的初始scs径流曲线数1.37gwqmnv浅层含水层产生“基流”的阈值深度1098gw_delayv地下水补给延迟系数10.4alpha_bfv基流alpha因子0.38escov土壤蒸发补偿系数0.49epcov植物吸收补偿因子0.29gw_revapv浅层地下水再蒸发系数0.10revapmnv浅层地下水再蒸发阈值522αv渠系输水损失蒸发比0.09sol_kr土壤饱和水力传导度0.52rchrg_dpv深含水层渗透比0.15βv渠系输水损失入渗比0.50sol_awcr土壤有效含水量1.11pipv管网漏损率0.10φv渠系水有效利用系数0.56~0.69ωv田间水有效利用系数0.95rv居工地耗水率0.2~0.45vv污水处理率1.0rev再生水利用率0表7中，v表示调参值替代原参数值；r表示原参数值乘以调参值。
[0071]
表8
表9本发明以swat模型为基础，开发一种分布式自然
‑
人工水循环模型，使其具有了自然
‑
人工水循环动态互馈及水资源开发利用与调控一体化模拟功能，可完成基于城市多水源供水和多水源灌溉供水模式的流域自然
‑
人工水循环过程模拟，在模型运行中，能够时刻保持自然水循环与人工侧枝水循环之间的动态互馈关系，不仅反映水循环过程对人工取用水的影响，同时还反映水资源开发利用与调控对水循环过程的实时干预作用，从而为深入了解强人类活动影响下的流域水循环机制，合理开发和利用水资源提供科学的参考依据。
[0072]
以上所述，仅为本发明的具体实施方式，但本发明的保护范围并不局限于此，任何属于本技术领域的技术人员在本发明揭露的技术范围内，可轻易想到的变化或替换，都应涵盖在本发明的保护范围之内。因此，本发明的保护范围应该以权利要求的保护范围为准。