大伙房流域降雨径流模型

大伙房流域降雨径流模型大伙房流域降雨径流预报模型又简称“DHF ”模型，该模型于1973年由辽宁省大伙房水库管理局刘爱杰、王本德等人提出，至今已使用30余年，为水库洪水调度做出了很大贡献。

“DHF ”模型是适用于我国湿润地区的超渗产流模型，目前已在辽宁省多个水库的水情自动测报系统中使用，效果较好。

建立在“DHF ”模型基础上的降雨径流预报方案，在大伙房流域经过调试和精度验证后进行使用，在使用中平均精度令人满意。

尤其在“957”特大洪水调度中，发挥了显著作用，准确预报出了第一非常溢洪道溢流时间和水库最高库水位，为省防制定调度决策提供了科学依据，使水库工程发挥了强大的调蓄作用，最大限度地配合了下游抢险，共减免下游直接经济损失74.89亿元。

“DHF ”模型由两部分组成，一是八参数超渗产流计算模型，引用双层入渗曲线进行扣损计算，并以抛物线描述表层蓄水量和下层渗率的分布状况；二是八参数变强度、变速度的经验单位线汇流计算模型，参数随降雨分布而变，采用“前期影响净雨”描述汇流速度的变化。

这是一个集总的概念模型，模型的参数多半在满足其物理意义的前提下确定，只有6个需要优选法选定或试错法确定。

1 大伙房模型概化流程流域下垫面分为表层、下层和地下水蓄存三部分，计算流程如图10-1所示。

2 大伙房模型产流计算产流模型将下垫面分为表层，下层和深层三部分。

表层土壤中的张力水蓄量与植物截流、填洼储存合称表层蓄水量aS ，其极值为表层蓄水容量S ；下层土壤中的张力水蓄量称为下层蓄水量aU ，其极值为下层蓄水容量U ；地下水储水层的蓄水量以a V 表示，其极值为地下水库蓄水容量V 。

图10-1 大伙房产流模型概化流程图其中：S 0-表层蓄水容量； U 0-下层蓄水容量；D 0-下层蓄水容量（U 0）与地下水库蓄水容量（V 0）之和； g-不透水面积占全流域面积之比值； K c -流域蒸散发能力与大水体蒸发量的比值； A-表层蓄水容量曲线形状参数； B-下渗率抛物线分配曲线形状参数； K-下层下渗曲线曲率；2.1 不透水面积上的径流计算根据上面的流程图，不透水面积上直接产生径流，通过下式来计算：c Ey g P =* （1）/E D C R P P E P K E C=-=- （2）式中P 为时段降雨量，KC 是流域蒸散发能力与ER 大水体蒸发量的比值，C 为日降雨观测的时段数，g 为流域不透水面积占总流域面积的比值。

2.2 地面壤中流的计算地面壤中流通过下式计算u y R f=- （3）2.2.1 表层扣损R 的计算当m C S P aS +<时，000(1()/())aC a m C R P S S S S P aS =+-+-+ （4）当m C S P aS +≥时，0()C a R P S S =-- （5）式中Sa 是表层流域平均蓄水量，Sm 是相应于Sa 的点蓄水量；a 为曲线形状参数。

Sm 和Sa 由下式计算得到：1/001(1/)a m a S aS S S ⎡⎤=--⎣⎦ （6）[]/1(1/)a M m M S S a S S =-- （7）式中MS 为表层点最大蓄水量，当m M S S =时，则a S S =，并有M S aS = （8）其中CP 通过下式来计算：C D CP P E y =-- （9）2.2.2 下渗计算时段平均总下渗率f 与a D的关系:10()a f Z D D =- (10)时段平均下层下渗率uf 与aD 的关系:20()u a f Z u u =- (11)其中2002/121;1k tu D k tZ e Z e --=-=-，1k 和2k 的关系如下:1200/k k u D = (12)1k 和2k 分别是反映总下渗和下层下渗曲线曲率系数；a D 是au 的函数，二者的关系式如下:00/001(1/)u D a a D D u u ⎡⎤=--⎣⎦(13)D 是下层蓄水容量u 与地下水库蓄水容量V 之和。

2.2.3 径流深计算公式总径流y 的计算公式：2202020()1Bn a Z u R y R Z u u Z u Z Bu ⎡⎤+=+-+-⎢⎥⎣⎦ (14)当220n R Z u Z Bu +≥时，得20()a uy R Z u u R f =--=- (15)地面壤中流uy 的计算公式：当110n Z D R Z BD +<时，00/11001110(1/)(1)u D Bn u a Z D R y R Z D u u Z D Z BD +=--+-(16)当110n Z D R Z BD +≥时，00/100(1/)u D u a y R Z D u u =-- (17)nD 和nu 通过下列两式计算00/()001(1/)u BD n a D BD u u ⎡⎤=--⎣⎦ (18)1/001(1/)B n a u Bu u u ⎡⎤=--⎣⎦(19)地下径流Ly 的计算：L W L y K R = (20)其中WK 为Ly 与LR 之比值，在干旱地区LR 全部渗入深层，不形成地下径流，则W K =；在湿润地区LR 全部形成地下径流流入河道，所以1W K =。

LR 为地下水库下渗强度，为总径流与壤中流之差，即：L uR y y =- (21)2.3 大伙房模型蒸散发计算（1）当降雨量超过流域蒸发能力时，即净渗雨强0C P >，产生下渗，土壤蓄水量增加，此时：表层蓄水量计算公式(1)001(1()/)a a t mt ct S S S P aS +⎡⎤=--+⎣⎦(22)式中t 为单位日，以8:00为日分界；mS 以式（6）计算当mt ct S P aS +≥时，上式可写成(1)0a t S S += (23)依照流程图的下渗程序，逐时段计算式，(1)a i ai ci iS S P R +=+- (24)式中i 为时段序；cP 、R 分别用（9）、（4）、（5）式进行计算。

下层土壤蓄水计算公式 据流程图，逐时段计算au 的公式(1)a i ai i iu u R y +=+- (25)若逐日计算，可用下式计算 若220nt t Z u R Z Bu +<时2(1)202020()(1)Bnt t a t at at Z u R u u Z u u Z u Z Bu ++=+---(26)若220nt t Z u R Z Bu +≥时(1)20()a t at at u u Z u u +=+- (27)式中的nu 与2Z 分别用式（19）及221k tZ e -=-求出。

（2）当日流域平均降雨量小于流域蒸发能力，即DP E <时，土壤蓄水量将减少进行蒸散发计算表层蒸散发公式假定：a.表层蒸发等于流域平均蒸散发能力；b.每次下渗后的初次蒸散发皆从地表开始。

uE 为表层可蒸发量，计算公式为：,1,,000(1)(1)B t B t a a u t E E E S aS aS -⎡⎤=---⎢⎥⎣⎦ (28) 式中：,,1,1,,1,c t u t t t B t cN A B t B t c tE E P E e E E E E ---=-===+∑ (29)因此，1t +日表层蓄水量计算式为,1,,a t a t u tS S E +=- (30)下层蓄水量au 计算公式假定：a.表层蓄水量蒸发尽后，该点的下层开始蒸发；b.下层蒸发量LE 与下层蓄水量au 成比例关系。

因此，()/L c u a E E E u u =- (31)当()0u S E -<时，则()/L c a a E E S u u =- (32)这样，1t +日下层土壤蓄水量以下式计算,1,a t a t Lu u E +=- (33)2.4 大伙房模型汇流计算大伙房汇流模型是根据多年预报经验，总结出的8参数变强度、变汇流速度的经验单位线，其线型公式为CCM T n t AA MnT t e nT K t DDM ⎥⎦⎤⎢⎣⎡=∂∂∙∙-πωπsin )/(3式中 t ∂∂ω——汇流曲线，是M T 与t 的函数；M T ——汇流曲线底宽；t ——汇流时间（时刻）； K3 ——汇流曲线的比例系数；n ——地下径流与地面壤中流汇流曲线底宽的比例系数； AA 、DD 、CC ——反映汇流曲线形状的参数；地表径流与地下径流具有不同的汇流曲线形状，可用形状参数区分；不同的降雨分布地表径流也具有不同的汇流曲线，其参数也不同。

M T 是时段净雨（径流）与前期影响净雨（已产生的径流）的函数，模型中假定： 0)(K T a B M y y L T -+= )(,,1,t T t a t a y y K y +=+式中 B L —特征河长/B L L B =；L —流域最大河长；0K —拟合线的斜率；ay —前期影响净雨量，表征河槽蓄水状况，制约汇流速度等；K —前期净雨影响程度衰减系数，反映已产生的径流在河槽中消退速度； 2.5 大伙房模型产流参数优选大伙房模型产流参数的最优化问题是多目标的有约束的非线性规划问题。

目标函数取为N 次洪水产流的拟合精度和N 次洪水计算产流与实际产流误差总和的线性加权，各取0.5的权重进行求解。

将问题变为：{}12min max max{*1*}..i i i w ACU w TYTC s t x x x -≤≤式中 ACU1-N 次洪水产流的拟和精度TYTC-N 次洪水计算产流与实际产流误差总和1/()/ACU N NTYTC TOTAL TYTC TOTAL ==-合格其中N 为洪水的次数，TOTAL 为N 次洪水的实际净雨的总和，TYTC 为N 次洪水的预报净雨的总和， N 合格为拟合的合格洪水的次数。

若某次洪水满足以下条件之一即为合格：计算产流与实际产流的绝对误差小于3mm ；计算产流与实际产流的相对误差小于20%且绝对误差小于20mm ；计算时判断合格的绝对误差和相对误差的标准可以根据具体情况来设定；本程序采用可变容差法求解问题，把多个约束求极小值问题变为一个单约束求极小值问题。

设所求解的非线性规划问题为：{}min ()()0(1,2,...,)..()0(1,...,)nf x x E hi x i m s t gi x j m p ∈==⎧⎨≥=+⎩本算法将上述两问题化为单约束问题进行求解:{}()min ()..()0n k f x x E s t T x ∈Φ-≥式中T(x)是约束破坏的估计量，表达式为：1/22211()()()pm i j j i j m T x h x U g x ==+⎡⎤=+⎢⎥⎣⎦∑∑其中，Uj 为Heaviside 算子。

()k Φ是第K 步搜索中给出的关于可行性的可变容差准则值，表达式为：()(1)()()2(0)1min ,12(1)k k k k i r m x x r m t -++⎧⎫Φ=Φ-⎨⎬+⎩⎭Φ=+∑其中 t-初始多面体的大小; m-等式约束数目; r=n-m ，()f x 的自由度;()k i x -nE 多面体中第K 次搜索时的第i 个顶点值; ()2k r x +-多面体（除hx ）外的形心顶点;由于各参数的初值相差较大，而且变化范围也各不相同，所以将参数归一化，即初值都取为1.0。