2新安江流域水文模型

2024年1月12日发(作者：三国演义历史背景)

2新安江流域水文模型60年代初，河海大学(原华东水利学院)水文系赵人俊等开始研究蓄满产流模型，配合一定的汇流计算，将模型应用于水文预报和水文设计。1973年，他们在对新安江水库做入库流量预报的工作中，把他们的经验归纳成一个完整的降雨径流流域模型——新安江模型。模型可用于湿润地区和半湿润地区的湿润季节径流模拟和计算。最初的新安江模型为两水源模型，只能模拟地表径流和地下径流。80年代初期，模型研制者将萨克拉门托模型与水箱模型中，用线性水库函数划分水源的概念引入新安江模型，提出了三水源新安江模型，模型可以模拟地面径流、壤中流、地下径流。1984至1986年，又提出了四水源新安江模型，可以模拟地面径流、壤中流、快速地下径流和慢速地下径流。三水源新安江模型一般应用效果较好，但模拟地下水丰富地区的日径流过程精度不够理想。在新安江三模型中增加慢速地下水结构就成为四水源新安江模型。 当流域面积较小时，新安江模型采用集总模型，当面积较大时，采用分块模型。分块模型把流域分成许多块单元流域，对每个单元流域做产、汇计算，得到单元流域的出口流量过程。再进行出口以下的河道洪水演算，求得流域出口的流量过程。把每个单元流域的出流过程相加，就求得了流域出口的总出流过程。划分单元流域的主要目的是处理降雨分布的不均匀性，因此单元流域应当大小适当，使得每块面积上的降雨分布比较均匀。并有一定数目的雨量站。其次尽可能使单元流域与自然流域相一致，以便于分析与处理问题，并便于利用已有的小流域水文资料。如果流域内有大中型水库，则水库以上的集水面积即应作为一个单元流域。因为各单元流域的产汇、流计算方法基本相同，以下只讨论一个单元流域的情况。新安江模型包括4个计算环节：蒸散发计算；流域产流计算；径流

划分；汇流计算。4个计算环节分别概化了流域降雨径流的主要产、汇流物理过程。2.1流域蒸散发计算各种水源的蒸散发计算模型均可采用两层蒸发模型或两层蒸发模型，一般根据实际情况选用。原则是在模拟径流精度相同的情况下，尽量采用参数少的两层蒸散发模型。蒸散发模型不考虑面上分布的不均匀性，但可考虑土湿垂向分布的不均匀性。2.2.1两层蒸散发模型两层蒸散发模型将土层分为上、下两层，各层蓄水容量分别为、（）。流域土层下渗蓄水和蒸散发计算过程按下述原则进行：降雨先补充上层，上层蓄满后再补充下层；蒸发时先蒸发上层，上层蓄水量蒸发殆尽后再蒸发下层。两层蒸散发模型的计算式如下：⑴：， ，⑵ ： ，，式中：——是流域降雨量，；——是流域蒸散发能力，；——是上层土壤蓄水量，；——是下层土壤蓄水量，；——是下层土壤蓄水容量，；——上层蒸发量，。——下层蒸发量，。2.2.2三层蒸散发模型三层蒸散发模型将土层分为上层、下层、深层三个土层，各层蓄

水容量分别为、、（）。流域土层下渗蓄水和蒸散发的计算过程与两层蒸发模型类似，按下述原则进行：降雨首先补充上层，上层蓄满后再补充下层，下层蓄满后再补充深层；蒸发时先蒸发上层，上层蓄水量蒸发殆尽后再蒸发下层，下层蓄水量蒸发殆尽后再蒸发深层。由土壤蒸发规律知，在久旱时，流域下层土壤含水量趋于干涸，深层土壤水分将以薄膜水、气态水的形式向上输送，使蒸发维持在一个较小的稳定数值。因此，考虑来自了深层的水分蒸发。为了使流域的蒸散发不小于某个最小的稳定蒸散发量，模型中考虑了一个深层蒸散系数，并在结构上进行了控制，使计算过程中的最小蒸散发量不会小于。三层蒸散发模型的计算式如下：⑴当：， ， ，；⑵ ，时：，，，；⑶，时：，，，；⑷，时：，，，。式中：——是深层蒸散发系数。2.2新安江两水源模型模拟流域出口断面的径流过程时，如果水源划分只考虑直接径流和地下径流，而不考虑壤中流，这样的新安江流域模型称为新安江两水源模型。2.2.1模型的结构新安江两水源模型的4个计算环节中（蒸散发计算、流域产流计算、径流划分、汇流计算），蒸发计算可采用两层蒸发模型，也可采用三层蒸发模型。产流计算采用蓄满产流模型。径流划分计算利用稳定下渗率

将径流划分为地面径流（或称直接径流）和地下径流两种水源。地面径流采用单位线汇流，地下径流采用一次线性水库汇流，即地下汇流计算的坡面汇流和河网汇流用一个线性水库的出流过程概化。模型把流域面积划分为透水面积和不透水面积两部分，不透水面积上的降水在满足蒸发后将直接转化为地面径流。透水面积上将发生下渗，下渗的水量一部分存储于土层（降雨过后耗于蒸发）。满足了流域土壤蓄水容量后的下渗水量才能转化为径流。不透水面积用参数表示，它是用流域内不透水面积占全流域面积的百分比表示的。新安江模型的输出是流域出流过程和流域蒸散发过程，输入则为时段降雨量、蒸发皿观测蒸发量。新安江两水源模型结构如下图所示：图2.1-1 新安江两水源模型结构示意图图中方框内标注为状态变量，方框外标注的是模型参数。

2.2.2模型的产流计算新安江模型产流部分的计算是蓄满产流模式， 蓄满产流指在流域包气带土壤湿度达到田间持水量以前不产流，所有的降雨都被土壤吸收；而在土壤湿度达到田间持水量之后，所有的降雨（除去同期的蒸散发）都产流。在产流后，流域包气带土壤的下渗能力为稳定下渗率，下渗的水分成为地下径流和壤中流，超蓄的部分成为地面径流。考虑到流域内各点的蓄水容量并不相同，实际产流常常是在部分面积上产图2.2-1流域蓄水容量曲线示意图 流，新安江模型引入一条流域蓄水容量曲线来刻划流域内各点蓄水容量的不均匀性，把流域内各点的蓄水容量概化成如图2.3所示的一条抛物线（也可概化成其它函数形式），其方程为： （2.2-1）式中：为产流面积；——是点蓄水容量, ；——是蓄水容量小于等于蓄水量的全部点的面积,；——是单元流域面积, ；——是蓄水容量曲线指数。——是流域内点最大蓄水容量,；其与流域平均蓄水容量的关系为: （2.2-2）引入蓄水容量曲线后，在降雨过程中，只有在满足了包气带蓄水容量的面积上才可能产生径流，其余面积上不产生径流。扣除蒸发后的有效降雨量中，未转化为净雨产流量的水分通过下渗进入流域包气带土层，补充土层蓄水量。引入流域蓄水容量曲线后的产流计算式为：

当时为部分面积产流

（2.2-3）当时为全流域面积产流

（2.2-4） 式中：

——为有效降雨量， ；——为与流域土壤初始蓄水量相应的前期影响雨量，由下式计算。 （2.2-5）2.2.3模型的水源划分及汇流计算1.新安江两水源模型的水源划分新安江两水源模型的水源划分只考虑直接径流和地下径流，在计算出流域总产流量后，进一步利用稳定下渗率将总径流量划分为地面径流（或称直接径流）量和地下径流量两种水源。两水源划分的计算式如下： (2.2-6)

(2.2-7)2.新安江两水源模型的汇流计算直接径流采用一条单位线汇流直接汇流到流域出口断面，即：

(2.2-8)式中：——第时刻的直接径流流量，；——第时段的直接径流量，；——单位线第时段末的流量，。地下径流采用一个线性水库汇流到流域出口断面，即： (2.2-9)

式中：——单位转换系数，，是流域面积，。流域出口断面时刻的流量为： (2.2-10)2.2.4模型的参数及调试 采用三层蒸发计算的新安江两水源模型共有9个参数，一条单位线。 ——流域蒸发折算系数，是流域蒸散发能力与蒸发皿蒸发量之比； ——深层蒸散发系数； ——不透水面积占全流域面积的百分比，％； ——流域平均蓄水容量(指张力水)，； ——流域平均上层蓄水容量, ； ——流域平均下层蓄水容量, ；——蓄水容量曲线指数；——稳定下渗率；——地下水消退系数；——单位线。当模型初定后就可上机调试，在调试的过程中，应注意各参数的径流响应特征，以及参数之间相依性的影响，以便调试能有目的地顺利进行。以下简介的调试技术可供参考。流域蒸发折算系数的调试，对一个具体流域来说，此参数完全靠优选。调试此参数主要考察年径流模拟误差，当年径流模拟误差达最小时，此值最优，但个别参数的最优并不表示模型总体最优。根据现有经验，年径流模拟误差控制在5％～8％左右就可以了。流域平均蓄水容量的调试，一般根据实测的降雨径流资料分析或经验选取后，适当微调即可。经验选定后，此值调试范围不大，比较容易确定。通常在南方湿润地区可选80～150左右，江淮一带约在110左右,

燕山东北东部地区，选150左右。流域平均上层蓄水容量，一般取10～20，流域平均下层蓄水容量一般取60～90。蓄水容量曲线指数的调试，表示流域蓄水容量分布的不均匀性，当全流域蓄满后，这个参数就不起作用了，因此，应该选取流域没有达到蓄满的那些洪水点据作为调试的依据。此值一般在0.2～0.5之间取值。注意值与值之间有相依性，二者对径流模拟的结果相互有影响。稳定下渗率的调试，对水源划分起决定性作用，但至今对的研究还不充分，目前可对作一些简单处理，如果取为一个常数使模拟效果不好时，可考虑作为变动参数。考察是否合适主要观察地下径流的模拟精度。地下水消退系数的调试，主要观察洪水退水段的拟合精度。不透水面积参数对模型拟合精度影响很小，一般取0～0.05。2.3新安江三水源模型模拟流域出口断面的径流过程时，如果水源划分考虑直接径流、壤中流和地下径流三种水源，这样的新安江流域模型称为新安江三水源模型。2.3.1新安江三水源模型结构新安江三水源的4个计算环节中（蒸散发计算、流域产流计算、径流划分、汇流计算）。流域的蒸散发计算和产流计算与新安江二水源模型相同，水源划分则完全不同，因此汇流计算也不相同。新安江三水源模型结构示意图如下图所示。图2.3-1 新安江三水源模型结构示意图2.3.2新安江三水源模型的水源划分

1.不考虑自由水蓄水容量分布不均匀的划分方法求得的产流量包括地面径流 、壤中流和地下径流 三部分。新安江三水源模型用一个自由水蓄水库解决水源划分，自由水蓄水库有两个出流孔，底孔为地下径流出流孔，边孔为壤中流出流孔。新安江模型考虑了产流面积（）的变化，自由水蓄水库实际只发生在产流面积上，其底宽为产流面积，显然它是随时间变化的。产流量进入水库即在产流面积上产生的径流深，也就是自由水蓄水库所增加的蓄水深，当自由水蓄水深超过其最大值时，超过部分成为地面径流。壤中流和地下径流按线性水库出流，其出流系数分别为和。底孔出流量和边孔出流量分别进入各自的水库，并按线性水库的退水规律流出，分别成为地下水出流和壤中流出流。模型认为蒸散发在张力水中消耗，自由水蓄水库的水量全部为径流。水源划分示意图见图2.3-2。水源划分计算式为：图2.3-2 不考虑自由水蓄水容量分布不均匀的划分示意图地面径流计算式：当 (2.3-1)当 (2.3-2)壤中流计算式：

(2.3-3)地下径流计算式：

(2.3-4)注意每步计算时应将自由水蓄水库的蓄量折算为当前时段产流面积上的深度。2.考虑自由水蓄水容量分布不均匀的划分方法流域自由水蓄水容量实际上也是不均匀的，作为图2.3-1表示的水源划分方法的改进，人们引入流域自由水蓄水容量分布曲线来刻化自由水蓄水量的不均匀性，其线型与流域蓄水容量曲线类似，曲线见图2.3-3，曲线方程如下式。 （2.3-5）

式中：为产流面积；——是点自由水蓄水容量；——是自由水蓄水容量小于等于蓄水量的全部点的面积；——是单元流域面积；——是自由水蓄水容量曲线指数。——是流域内点最大自由水蓄水容量；与流域平均自由水蓄水容量的关系为： （2.3-6）由于认为在产流面积上才有自由水，因此，产流面积上的点最大自由水蓄水容量是随产流面积变化而变化的，见图2.3-3。对当前时段，已知，则将其代人（2.3-5）得 （2.3-7）1.0图2.3-3 考虑自由水蓄水容量分布不均匀的划分示意图从（2.3-7）解出得： （2.3-8）同样，在产流面积上的自由水平均蓄水容量与该面积上的点最大自由水蓄水容量的关系为： （2.3-9）在划分水源时，产流面积上只有在满足了蓄水容量的水库面积上才可能产生地面径流；其余面积上不产生地面径流。进入水库的水分扣

除地面径流量后，剩余的补充水库蓄水量S，这时的水源划分计算式为：⑴地面径流出流当时为部分面积产流

（2.3-10）当时为全流域面积产流

（2.3-11）式中：——为计算时段初的自由水蓄水库的蓄水量（用产流面积上的平均深度表示）；——为计算相应于的自由水蓄水库前期影响入流量，由下式计算； （2.3-12）进入水库的水分扣除地面径流量后，剩余的补充水库蓄水量： （2.3-13）自由水蓄水库的时段入流量为时段产流量，但由（2.3-3）和（2.2-4）计算的产流量是用流域平均水深表示的水量，应将折算为产流面积上表示的水深。可以证明折算为产流面积上表示的水深值恰为流域有效降雨值，证明如下：为折算为产流面积上表示的水深，根据降雨～径流关系的性质，水文学原理中已经证明相对产流面积与产流量之间有关系：

（2.3-14）从此式可看出⑵壤中流出流：

（2.3-15）⑶地下径流出流： （2.3-16）由于,都是产流面积上的深度表示的水量，计算所得、、应折算回

全流域深度表示的水量，即：、、在当前时段计算结束后应乘。2.3.3流域汇流计算汇流计算可以采用两种方法：1.传统流域汇流这种汇流方法是：地面径流采用单位线汇流，壤中流和地下径流采用线性水库汇流至出口断面后与单位线汇流结果叠加，最后得到流域出口断面的模拟出流。即图2.3-2中下面的两个水库的出流就是壤中流和地下径流在流域出口断面的出流流量，该两个水库分别称为壤中流调蓄水库和地下径流调蓄水库。⑴壤中流汇流计算 （2.3-17）⑵地下径流汇流计算 （2.3-18）式中： ——时段壤中流出流流量，；——壤中流消退系数，可假定为；——时段地下径流出流流量，；——地下径流消退系数，可假定为；——是径流深转换为流量时的单位转换系数，.上式中假定了调蓄水库与自由水水库的出流系数相同。在实践中，为延缓壤中流和地下径流出流，有人把下面的两个水库的出流看作是坡地调蓄水库，而在它们后面再分别串联一个水库反映壤中流和地下径流受河槽调蓄的作用。2.河网汇流单位线汇流河网汇流单位线也称时变河网汇流单位线，这种汇流方法忽略地面径流的坡面汇流阶段，直接将地面径流产流量转换为流量作为河网地面径流入流（看作集中在计算时段末入流）。地下径流和壤中流汇流均采

用线性水库作坡地调蓄后得到相应河网入流、，即：把图2.3-2或图2.3-3中的转换为流量作为河网地面径流入流，下面的两个水库的出流看作壤中流和地下径流经坡地调蓄后注入河网的入流量。三种水源的河网入流流量相加就是所谓的河网总入流。河网总入流经流域河网调蓄后得到流域出口断面的流量过程。河网汇流单位线的计算式如下：

（2.3-19）式中： ——时段河网的总入流，，；——河网消退系数；——河网汇流时间，；——时刻的河网出流（即时刻的流域出口断面流量），；——时刻的河网总入流经河网调蓄后，延迟时间的出流，。是随时间而变化的时变参数，它取决于时刻的河网的总入流量和流域的河槽特性，用下式计算： （2.3-20）式中： ——反映流域河槽特性的系数。2.3.4线性水库不同单位时段的出流系数转换 从(2.3-17)和(2.3-18)，知道新安江模型的自由水水库出流是按线性水库出流处理的。在客观世界中，流域的降雨、蒸散发、径流等水文物理量都是连续量，流域内的蓄水量也是连续量，观测和计算时都只能处理为一个一个的时段平均量；从数学观点来看，这实际上是将连续量处理为离散量，将表达物理过程的连续方程处理为表达离散量的差分方程，这样就必然带来差分误差，特别是水文模型计算通常采用的时段较大，如果在模型的算法设计上不注意，差分误差可能较大。为了减小计算的差分误差，本可以缩小计算时段，但是，在水文上常常受到限制，这是由于原始水文观测数据时段本身较大，如果在所有的计算环节上都缩小计算时段并不能提高成果的精度。尽管如此，在模型的算法上，至

少应要求算法不会带来更大的差分误差。仔细考察三水源新安江流域水文模型计算的各个环节，就会发现(2.3-17)和(2.3-18)的自由水水库出流计算可能会出现较大的差分误差，为此，有的模型程序设计者在算法上作了一定处理，以设法减小差分误差。基本的想法是在这个计算环节上减小计算步长，即将模型总体的计算时段在这个环节上进一步细分；为此，需要将水库长时段的出流系数转换为短时段的出流系数。设水库在一个计算时段内的入流量为，为减小差分误差，将等分为份入流，这实际上是将原计算时段细分为个更小的计算时段，每个时段的水库入流为。可以证明水库计算时段为的出流系数和与计算时段为的出流系数和之间存如下关系： ； （2.3-21）证明如下：1.只有一个出流孔的线性水库不同时段出流系数的转换先对只有一个出流孔的线性水库进行证明。设线性水库的出流方程和退水时的水量平衡方程如下：

（2.3-22）

（2.3-23）式中： ——为线性水库的瞬时出流系数；——为线性水库时刻的蓄量；——为线性水库时刻的出流流量。合并（2.3-22）和（2.3-23）得到: （2.3-24）在上积分上式得 （2.3-25）此式称为线性水库退水方程。此式的另一个表达形式是：

（2.3-26）这是（2-27）在区间上积分的结果。令，即，取为一个单位计算时段，则上式记为： （2.3-27）称为线性水库消退系数（又称退水系数）。在区间上积分，得到水库在时段内流出的水量。

（2.3-28）令，即，取为一个单位计算时段，则上式记为 （2.3-29）（2.3-29）中的就是线性水库单位时段为的出流系数。比较（2.3-29）和（2.3-22）可知和的意义是完全不同的，总是与一个单位时段相联系的，单位时段不同，其值也不同。令：

记单位时段为的水库消退系数为。由（2.3-27）的退水方程，可推导得 （2.3-30）此式中的1代表一个时段，表示个时段。显然。对单位时段的退水方程，有； 或 （2.3-31）这里，1代表一个时段。比较（2.3-30）最后一个式子与（2.3-31）式可知，线性水库单位计算时段为的出流系数与单位计算时段为的出流系数之间有如下关系：

（2.3-32）已知对单位计算时段，有；记单位计算时段的出流系数为，则可推得 （2.3-33）此式就是所需要得结果。2.有二个出流孔的线性水库不同时段出流系数的转换记单位计算时段为的出流系数分别为和，单位计算时段为的出流系数分别为和，不同计算时段的两孔总出流系数分别为和，代人（2.3-

33）式得到： （2.3-34）要求解和还须增加一个条件，可假定不同计算时段的两孔出流系数分配比例相同，即成立： （2.3-35）利用（2.3-34）和（2.3-35）消去得到： （2.3-36）利用（2.3-34）和（2.3-35）消去得到： （2.3-37）利用2.3-36）和（2.3-37）就可进行线性水库不同时段的出流系数转换。2.3.5新安江三水源模型的参数调试采用河网汇流单位线的新安江三水源模型共有14个参数，根据参数的功能可把参数分为4类： ⑴与蒸散发计算有关的参数：蒸散发折算系数，上层土壤蓄水容量下层土壤蓄水容量，深层蒸散发系数。⑵与产流计算有关的参数：流域蓄水容量，流域蓄水容量曲线指数，流域不透水面积比值。⑶与水源划分有关的参数：流域自由水蓄水容量，流域自由水蓄水容量曲线指数，自由水蓄水库壤中流出流系数与自由水蓄水库地下径流出流系数。⑷与汇流有关的参数：壤中流消退系数和地下径流消退系数；河网汇流单位线参数和。壤中流汇流调蓄水库和地下径流汇流调蓄水库与自由水水库是不同的水库，但如果假定二调蓄水库的出流系数与自由水水库对应的径流的出流系数相同，即假定，，则参数可减为12个。4类参数相对独立，同类参数有一定相依性。1.参数初定根据流域自然地理资料、水文气象资料和前人以积累的经验初定参

数。⑴蒸散发折算系数可采用汛期连续大雨，流域蓄水量处于蓄满状态的累计降雨、累计径流、累计蒸发观测资料，用下式估计: （2.3-38）⑵流域蓄水容量、上层土壤蓄水容量、下层土壤蓄水容量值一般根据实测的降雨径流资料分析或经验选取后，适当微调即可。经验选定后，此值调试范围不大，比较容易确定。通常在南方湿润地区可选80～150左右，北方半湿润地区可达170，江淮一带约在110左右,燕山东北东部地区，选150左右。 中包含有植物截留，少林地5～10，多林地10～20左右。在60～90范围。 ⑶流域蓄水容量曲线指数此值一般在0.2～0.5之间取值，通常流域越小，该值越小，特小流域可能小于0.2。1000平方公里以上流域一般大于0.3。 ⑷不透水面积比值值对模型拟合精度影响很小，一般取0～0.05。⑸深层蒸散发系数 南方多林地区左右 ；北方半干旱地区左右。 ⑹流域自由水蓄水容量 一般由经验初定，土层浅薄的山区，土厚多林地区，甚至更大。 ⑺流域自由水蓄水容量曲线指数 根据经验，常在1～1.5左右。⑻自由水蓄水库壤中流出流系数与自由水蓄水库地下径流出流系数一般来说，约为0.7左右，二者分配比与具体流域有关，只能调试决定。⑼壤中流消退系数和地下径流消退系数：一般趋近0.9；约为0.98～0.998，具体取值与具体流域有关，可

进一步调试决定。⑽河网汇流单位线参数和在0.02左右，与具体流域的河网特性有关，可进一步调试决定。与具体流域的河网特性有关，为计算方便，通常取计算单位时段的整倍数。2.参数率定首先率定第一类参数，主要观察率定期的总水量是否与实测水量相近。其次率定第二类参数，主要观察次洪水的水量是否与实测水量相近。第三类参数的率定，主要观察次洪水的退水段过程的水量是否与实测水量相近。第四类参数的率定，主要观察次洪水的过程线是否拟合良好。参数率定过程中要反复对参数进行多遍调试，直至满足精度要求。2.3.6新安江三水源模型的VB程序介绍新安江三水源，两层蒸发模型（时变单位线）VB源程序Function XingAnJang4() As Integer Dim Fno As Integer, I As Integer，Nian(2) As Integer，AX(20)As Single, Qo As Single Dim Wm As Single, Wum As Single, Wlm As Single, k As Single,Sm As Single, b As Single Dim Bx As Single, FR0 As Single, FAs Single, U As Single, X As Single, XX As Single Dim Kss As Single, Kg As Single, CRss As Single, CRg AsSingle, Immp As SingleDim P As Single, E As Single, Pe As Single, Eu As Single, El

As Single, EE As Single

Dim R As Single, Rs As Single, Rss As Single, Rg As Single, QAs Single, Qz As Single,Dim W As Single, Wu As Single, Wl As Single, S As Single, WmmAs Single, Smm As Single

Dim QRss0 As Single, QRg0 As Single, Cr As Single, Cs AsSingle, Cx As SingleDim Qa() As Single Qo = 5 'Qo是基流。 Open "e:test流域模型" For Input As #1 Open "e:test流域模型" For Input As

#5 Open "e:test流域模型" For Output As

#6 Fno = 1ReadRecord Fno:TurnRecord Nian, AX, 11 + 1, 1: Clo (Fno)Wm = AX(1): Wum = AX(2): k = AX(3): b = AX(4)Sm = AX(5): Bx = AX(6): Kss = AX(7) Kg = AX(8)CRss = AX(9): CRg = AX(10): Immp = AX(11) W = 85: Wu = 1: S = 0: FR0 = 0: Dt = 3: F = 2928: U = F /(3.6 * Dt): QRss0 = 0 QRg0 = 10.1 Wmm = (1 + b) * Wm: Smm = (1 + Bx) * Sm: Wlm = Wm - WumWl = W - Wu: X = 1 / (1 + b): XX = 1 / (1 + Bx) Cr = 0.025: Cx = 6

Fno = 5: n = 0: Do Until EOF(Fno): Line Input #Fno, TextLien:n = n + 1: Loop: Seek #Fno, 1 ReDim Qa(n + 1000) For I = 1 To n ReadRecord Fno: TurnRecord Nian, AX, 4 + 1, 1 n1 = Nian(1): m1 = Int(AX(1)): P = AX(2): E = AX(3): qb =AX(4)

Pe = P - k * E If Pe > 0 Then Eu = k * E: El = 0 XuManChanLiu1 Pe, R, Wm, Wmm, Wum, Wlm, W, Wu, Wl, b, X El R = 0 EuEl P, k, E, Wlm, Wu, Wl, Eu, El End IfRsRssRg1 Pe, R, S, Sm, Smm, Rs, Rss, Rg, Bx, XX, FR0, Kss,Kg, Immp EE = Eu + El: W = Wu + Wl HeWangZongRuLiu Rs, Rss, Rg, QRss0, QRg0, Qz, CRg, CRss, Qo,

U Tao = 1Cs = 1 - Cr * Qz ^ 0.4 Qa(I + Tao) = Cs * Qa(I) + (1 - Cs) * QzQ = Qa(I + Tao)Write #6 I, P, EE, R, Rs, Rss, Rg, W, Q, qb End IfNext IClo (6) = CStr(666)End Function

Public Sp(20) As VariantFunction ReadRecord(Fno) As Integer Dim TextLien As String, Stext As String, sss As String Dim I As Integer, j As Integer, LineLen As Integer j = 1

Line Input #Fno, TextLien LineLen = Len(TextLien) For I = 1 To LineLen Stext = Mid(TextLien, I, 1) Select Ca Asc(Stext) Ca 45, 46, 48, 49, 50, 51, 52, 53, 54, 55, 56, 57 sss = sss + Stext If I = LineLen Then Sp(j) = sss End If Ca El If sss <> "" Then Sp(j) = sss: j = j + 1: sss = "" End If End Select Next I ReadRecord = jEnd Function

'蓄满产流函数一Function XuManChanLiu1(Pe, R, Wm, Wmm, Wum, Wlm, W, Wu, Wl, b,X) As Integer Dim Dr As Single, A As Single A = Wmm * (1 - (1 - W / Wm) ^ X) If Pe + A < Wmm Then R = Pe - Wm + W + Wm * (1 - (Pe + A) / Wmm) ^ (1 + b) El R = Pe + W - Wm

End If If Wu + Pe - R < Wum Then Wu = Wu + Pe - R El Wl = Wl + Wu + Pe - R - Wum Wu = Wum If Wl > Wlm Then Dr = Wl - Wlm Wl = Wlm R = R + Dr End If End IfEnd Function

'两层蒸发函数Function EuEl(P, k, E, Wlm, Wu, Wl, Eu, El) As Integer If Wu + P > k * E Then Eu = k * E Wu = Wu + P - Eu El = 0 El Eu = Wu + P El = (k * E - Eu) * Wl / Wlm Wu = 0 Wl = Wl - El End IfEnd Function

'三水源划分函数Function RsRssRg1(Pe, R, S, Sm, Smm, Rs, Rss, Rg, Bx, XX, FR0,Kss, Kg, Immp) As Integer Dim j As Integer, nn As Integer Dim FR As Single, Au As Single, Smmf As Single, Smf AsSingle, RR As Single Dim Rs0 As Single, Rss0 As Single, Rg0 As Single, Kss0 As

Single, Kg0 As Single Rs = 0: Rss = 0: Rg = 0 If R > 0 Then FR = R / Pe S = S * FR0 / FR If FR < 1 Then Smmf = Smm * (1 - (1 - FR) ^ (1 / Bx)) El Smmf = Smm End If Smf = Smmf / (1 + Bx) nn = 5 RR = Pe / nn Kss0 = (1 - (1 - Kss - Kg) ^ (1 / nn)) / (1 + Kg / Kss) Kg0 = Kss0 * Kg / Kss For j = 1 To nn If S < Smf Then Au = Smmf * (1 - (1 - S / Smf) ^ XX) El Au = Smmf S = Smf

End If If (RR + Au > Smmf) Then Rs0 = (RR + S - Smf): S = Smf El Rs0 = (RR - Smf + S + Smf * (1 - (RR + Au) / Smmf)^ (1 + Bx)) S = S + RR - Rs0 End If Rs = Rs + Rs0 Rss0 = S * Kss0: Rss = Rss + Rss0 Rg0 = S * Kg0: Rg = Rg + Rg0 S = S - (Rss0 + Rg0) Next j FR0 = FR El Rs = 0 Rss = S * Kss Rg = S * Kg S = S - (Rss + Rg) End If If(Pe<=0) Then Pe=0Rs = Rs * FR0 * (1 - Immp) + Pe * Immp Rss = Rss * FR0 * (1 - Immp) ‘当Rs＝0时，自由水水库面积仍前一时段的面积。 Rg = Rg * FR0 * (1 - Immp)End Function

'三水源河网总入流函数（地面径流直接入流）

Function HeWangZongRuLiu(Rs, Rss, Rg, QRss0, QRg0, Qz, CRg,CRss, Qo, U) As Integer Dim Qs As Single, QRss As Single, QRg As Single Qs = Rs * U QRss = CRss * QRss0 + (1 - CRss) * Rss * U QRg = CRg * QRg0 + (1 - CRg) * Rg * U Qz = Qs + QRss + QRg + Qo QRss0 = QRss: QRg0 = QRgEnd Function