1. 引言
防洪设计中，设计断面上游有调蓄作用较大的梯级水库或水库群，对洪水具有调节作用，改变了天然洪水状况，直接影响下游设计断面的设计洪水。这些水库的调洪对设计断面洪水的影响，除了与水库本身的库容、泄洪能力等特性有关外，还与设计断面以上各部分面积上的洪水组成情况有关。因此，分析设计断面设计洪水的地区组成，对研究流域开发方案，计算水库对下游的防洪作用，进行梯级水库或水库群的联合洪水调度等均有重要意义[1] 。
由于流域内暴雨分布不均，各地区洪水来量不同但有一定的相关性，且各组成部分洪水的相互遭遇组合具有明显的随机性，给实际分析带来了很大的困难。我国《水利水电工程设计洪水计算规范》[2] (SL44- 2006)(以下简称《规范》)规定采用地区组成法(包括典型年法和同频率组成法)来拟定设计洪水的地区组成，尤以同频率组成法在工程实践中应用最为广泛。同频率组成法的基本出发点是指定某一局部地区的洪量与下游设计断面的洪量同为设计频率，其余洪量再根据水量平衡原则分配到流域的其他地区[3] 。然而，这只是人为假定的一种情形，并没有充分考虑各分区洪水的空间相关性，统计理论基础不足。此外，其设计成果尽管主观上考虑“对防洪不利”的原则，而客观上却往往使防洪断面达不到预定的防洪标准，并且对于复杂的梯级水库防洪设计也只能采用简单的逐级同频率法，方法较为单一，成果的不确定性较大[4] [5] 。
洪水地区组成的实质是多个相关随机变量的组合，科学地描述这一问题的方法是给出各地区洪水的联合分布。Copula函数可以描述水文变量的相关性结构，能够灵活地构造边缘分布为任意分布的水文变量联合分布，在描述具有多个属性的水文事件方面及多个水文事件遭遇规律方面具有很大的潜力，近年来在水文领域得到了广泛应用[6] -[10]。有学者将其引入到设计洪水地区组成的分析与计算[11] [12] ，构造了上游水库断面和区间洪水的联合分布，推求了设计断面上游有单个水库时两种具有统计意义的地区组成，即条件期望组成和最可能组成。然而，基于Copula函数的设计洪水地区组成仅限于单个水库的研究，对于梯级水库尚未涉及，难以满足工程实际需求。
本文将其推广到设计断面上游有梯级水库的情形，利用Copula函数推导得到梯级水库防洪系统设计洪水的条件期望组成和最可能组成。以清江流域为例，计算高坝洲断面的设计洪水地区组成，定量分析上游水布垭–隔河岩梯级水库的调蓄作用对其设计洪水的影响，并与《规范》规定的同频率组成进行比较，以期为梯级水库设计洪水的地区组成分析计算提供一条新的途径。
2. 梯级水库设计洪水地区组成的定义
由上、下两个水库组成的梯级水库是最常见的，具有一定的代表性，多级水库可以看成是两级水库的各种组合，因此本文以两个串联水库组成的梯级水库为例，探讨其设计洪水的地区组成。如图1，A1、A2分别为上、下水库，C为设计断面，B1、B2分别为水库A1、A2之间和水库A2、断面C之间的区间；随机变量X、Y1、Y2、Z1和Z分别表示水库A1、区间B1、区间B2、水库A2和断面C的天然来水量，取值依次为x、y1、y2、z1和z。
受上游A1-A2梯级水库的影响，分析断面C设计洪水的地区组成需要研究上水库A1、上区间B1和下区间B2三部分洪水的组合。当断面C出现某一频率p的洪量zp时，由水量平衡原理得：
(1)
设计洪水地区组成的实质是给定断面C的设计洪量zp，按照“可能发生”和“对防洪不利”的原则，在满足式(1)约束条件下分配zp，得到组合(x，y1，y2)。
从实际系列中选择有代表性，对防洪不利的大洪水年作为典型年；按设计洪水地区组成方案得到的洪量分配结果(x，y1，y2)放大该典型年各分区的洪水过程线；将水库A1断面洪水过程线经调洪演算得到下泄洪水过程，与区间B1洪水组合后，输入到水库A2进行调洪演算所得的下泄洪水过程，再与区间B2洪水组合，就可以推求出设计断面C受上游A1-A2梯级水库调蓄影响的洪水过程线及洪峰设计值。
3. 同频率地区组成
当断面C出现某一频率p的洪量zp时，水库A1断面、区间B1和区间B2来水有无数种可能组合。工程实际中，为了分析和比较设计洪水不同地区组成的防洪效果，常需拟定若干种有代表性的地区组成方案，经调洪演算和综合分析，从中选取能满足工程设计要求的，作为设计的依据。这些方案必须是可能发生的，而且对防洪来说是偏于不利的。
Figure 1. Sketch diagram of cascade reservoir flood control system
图1. 梯级水库防洪系统示意图
计算梯级水库设计洪水的地区组成时，同频率组成法一般自下而上，分两步分配断面C的设计洪量zp。第一步确定水库A2断面或区间B2出现的洪量与设计断面C同频率，第二步确定水库A1断面或区间B1出现的洪量与水库A2断面同频率。由乘法原理知，完成zp的分配分两步，每步有两种可能，共有22 = 4中不同的洪水组成方案，列于表1。
上述四种同频率组成方案，最后选用哪种组成作设计，要视分析的结果并结合工程性质和需要综合确定。这四种组成只是有一定代表性的地区组成，是设计考虑的四种特殊情况，实际上它们既不是最可能出现的地区组成，也不一定是最恶劣的地区组成[13] 。此外，这四种组合出现的可能性也不一样。一般来说，当某部分地区的洪水与下游断面的相关关系比较密切时，二者同频率的可能性比较大，反之若某部分地区的洪水与下游断面的相关关系较差时，则不宜采用下游与该部分地区同频率地区组成的方式。
同频率组成指定某一分区出现与断面C出现同频率的洪水，其余两分区分别出现相应的洪水，没有充分考虑各分区洪水的空间相关性，统计基础不足。下面，我们将利用Copula函数推导考虑空间相关性的两种梯级水库设计洪水地区组成，即条件期望组成和最可能组成，对传统的同频率组成进行补充完善。
4. 基于Copula函数的条件期望组成和最可能组成
4.1. Copula函数
Copula函数是定义域为[0,1]均匀分布的多维联合分布函数，令
为一个n-维分布函数，其边缘分布为
。由Sklar定理，则存在一个n-Copula函数C，使得对任意
[14] ：
(2)
式中：
为Copula函数的参数。
Sklar定理是Copula函数构造联合分布的理论根据，它可以将联合分布分为边缘分布和相关性结构两部分，形式灵活，构造简单。大量研究表明，Archimedean Copula函数族中的Gumbel-Hougaard Copula函数适用于变量存在正相关的情形，比较适合构造河流上下游站点的洪峰之间或洪量之间的联合分布[6] -[8] [15] 。因此，本文将分别采用二维和三维非对称Gumbel-Hougaard Copula函数来构造X、Y1和Y2之间的二维和三维联合分布。
二维Gumbel-Hougaard Copula函数表达式为[6] ：
(3)
式中：u和v分别为边缘分布函数，参数θ与Kendall秩相关系数
的关系为
(4)
Table 1. Four equivalent frequency composition schemes of cascade reservoirs
表1. 四种同频率洪水地区组成方案
注：p'为
对应于随机变量Z1频率曲线的设计频率。
三维非对称Gumbel-Hougaard Copula函数表达式为[16] ：
(5)
式中：u、v和w分别为边缘分布函数，参数θ1和θ2通过极大似然法进行估计。
4.2. 条件期望组成
当水库A1断面出现洪量x时，区间B1、B2所对应的洪量y1、y2并非是唯一的，只是出现不同取值的概率不同而已，即存在条件概率分布函数
和
。引入Copula函数后条件概率函数的表达式可写成如下形式：
(6)
(7)
式中：
、
和
分别为X、Y1和Y2的概率分布函数；
和
分别为X、Y1和X、Y2的联合分布函数。
当水库A1断面出现洪量x时，区间B1和区间B2的洪量存在条件期望值，表达式如下：
(8)
(9)
式中：
和
分别为
和
的密度函数；
和
分别为
和
的密度函数；
和
分别为
和
的反函数；
和
分别为Y1和Y2的概率密度函数。
假设
的设计值
(
为试算得到的水库A1断面设计频率，与断面C设计频率
一般不同)，此时，若
、
的条件期望值
和
满足：
(10)
则
为断面C设计洪水的条件期望组成。
条件期望组成是一种平均情况，具有一定代表性，如同频率地区组成一样，既不是最可能出现的地区组成，也不一定是最不利的地区组成。由于它考虑了水库A1断面、区间B1和区间B2洪水的空间相关性，存在一定的统计基础，因此比同频率法更具说服力。
4.3. 最可能组成
不同的地区组成发生的概率是不同的，设计上通常只关心可能发生且对下游防洪不利的组成。用X、Y1和Y2的联合概率密度函数值大小来衡量该组合发生的相对可能性大小，其表达式为[14] ：
(11)
式中：
为
的密度函数；
为X的概率密度函数。
当断面C出现设计洪量
时，为了寻求水库A1断面、区间B1和区间B2的最可能组成，即为求解
、
和
在满足
的条件下
的最大值。将
代入式(11)得：
(12)
欲使
取最大值，需满足：
(13)
将式(13)化简得：
(14)
式中：
，
，
，
；
、
和
分别为相应概率密度函数的导函数。
若X、Y1和Y2均服从P-III分布，则式(14)可表述为：
(15)
式中：
、
、
，
、
、
和
、
、
分别为X、Y1、Y2分布(P-III)的参数。利用牛顿迭代法求解非线性方程组式(15)可以得到水库A1断面、区间B1和区间B2洪水的最可能组成
。
5. 实例研究
清江流域面积17,000 km2，干流中下游建有水布垭–隔河岩–高坝洲梯级水库。本实例中，图1为A1为水布垭水库，坝址以上流域面积10,860 km2，下游92 km处建有隔河岩水库A2，控制面积14,430 km2，两水库为长江错峰各预留5.0亿m3防洪库容。C为高坝洲断面，位于隔河岩水库下游50 km处，控制面积15,650 km2，为径流式电站。B1、B2分别为水布垭–隔河岩区间流域(以下简称水–隔区间)和隔河岩–高坝洲区间流域(以下简称隔–高区间)，面积分别为3570 km2和1220 km2。以高坝洲断面设计洪水地区组成为例，对三种地区组成进行比较研究，分析上游水布垭–隔河岩梯级水库的调洪对高坝洲断面设计洪水的影响。
5.1. 边缘分布及联合分布的确定
根据水布垭–隔河岩梯级水库对洪水的调节特性，选取3天为有效控制时段。令Qm和W3分别表示各分区天然洪峰流量和最大3天洪量。我国设计洪水理论频率曲线线型一般采用P-III型分布，其概率密度函数为[2] ：
(16)
式中：
、
和
分别为形状、尺度和位置参数。
参数估计结果见表2，采用χ2检验法[17] 对6个随机变量进行分布拟合检验，结果见表2。可知，在5%的显著性水平下各变量均通过了检验。
采用二维或三维Gumbel-Hougaard Copula函数分别建立了4个联合分布：a) 水布垭断面W3和水–隔区间W3；b) 水布垭断面W3和隔–高区间W3；c) 水–隔区间W3和隔–高区间W3；和d) 水布垭断面W3、水–隔区间W3和隔–高区间W3。分别采用水布垭断面、水–隔区间和隔–高区间1965~2010年的资料，得到的3个二维联合分布(a)、(b)和(c)的秩相关系数和采用式(5)估计的参数结果均列于表3。三维联合分布(d)的参数采用极大似然法估计，结果也列于表3。经验分布与理论分布分别见图2(a)、(b)、(c)和(d)。可知，经验频率与理论频率值的拟合情况都较好。
Table 2. Estimated parameters of marginal distributions in each region
表2. 各分区边缘分布参数估计结果
(a) 水布垭断面和水–隔区间(b) 水布垭断面和隔–高区间 (c) 水–隔区间和隔–高区间(d) 水布垭断面、水–隔区间和隔–高区间
Figure 2. Theoretical and observed joint probability distributions
图2. 联合理论分布和经验分布
Table 3. Estimated parameters of joint distributions
表3. 联合分布分布参数估计结果
采用离差平方和准则(OLS)评价联合分布的拟合情况，计算公式如下[16] ：
(17)
式中：
和
分别为经验频率与理论频率；
为资料系列长度。
表3给出了4个联合分布的OLS值，各联合分布的OLS值均较小，不超过0.05，可知拟合情况较好。图2以及表3的结果表明所建立的4个联合分布是合理的。
5.2. 高坝洲断面设计洪水的地区组成
利用以上基于Copula函数的条件期望组成和最可能组成，同时结合《规范》推荐的同频率组成，分析计算高坝洲断面的设计洪水地区组成。
1) 同频率组成
采用同频率组成拟定高坝洲断面的设计洪水地区组成方案时，需要分为两步。第一步确定隔河岩断面或隔–高区间与高坝洲断面发生同频率洪水，第二步确定水布垭断面或水–隔区间与隔河岩断面发生同频率洪水，共有四种不同的方案：
①同频率Ⅰ。高坝洲断面发生设计频率洪水时，隔河岩断面发生同频率洪水，隔–高区间发生相应的洪水；水布垭断面发生与隔河岩断面同频率洪水，水–隔区间发生相应的洪水。
②同频率Ⅱ。高坝洲断面发生设计频率洪水时，隔河岩断面发生同频率洪水，隔–高区间发生相应的洪水；水–隔区间发生与隔河岩断面同频率洪水，水布垭断面发生相应的洪水。
③同频率Ⅲ。高坝洲断面发生设计频率洪水时，隔–高区间发生同频率洪水，隔河岩断面发生相应的洪水；水布垭断面发生与隔河岩断面同频率洪水，水–隔区间发生相应的洪水。
④同频率Ⅳ。高坝洲断面发生设计频率洪水时，隔–高区间发生同频率洪水，隔河岩断面发生相应的洪水；水–隔区间发生与隔河岩断面同频率洪水，水布垭断面发生相应的洪水。
图3是根据1965~2010的各分区实测资料点绘的各分区与下游设计断面最大3天洪量的设计频率相关关系图，直线为同频率线。点据越靠近同频率线，证明该分区与设计断面发生同频率洪水的可能性越大。
分析两步的同频率情况，第一步：从图(a)可以看出隔河岩断面极有可能与高坝洲断面发生同频率洪水；图(b)的结果表明当高坝洲断面发生5年一遇(设计频率20%)以上洪水时，隔–高区间洪水频率系统偏大，不宜机械地采用该分区与高坝洲断面同频率的组成方式。因此，第一步选取隔河岩断面与高坝洲断面同频率。第二步：从图(c)和(d)中不难发现，水–隔区间与隔河岩断面的洪水相关关系很差，水布垭断面则相对较好，故第二步选取水布垭断面与隔河岩断面同频率相对较为合理。两步相比之下，隔河岩断面与高坝洲断面同频率的可能性又要显著大于水布垭断面与隔河岩断面同频率。这是因为隔–高区间面积仅占高坝洲断面控制面积的7.8%，而水–隔区间面积占隔河岩断面控制面积的比例高达24.7%。区间面积比例越小，两水库坝址断面的洪水特性越密切。综上所述，四种同频率方案中同频率Ⅰ方案发生的可能性最大，较为合理，但仍然具有一定的人为性。为了进行比较，本文将选取该方案分析高坝洲断面的设计洪水地区组成，计算结果见表4。
2) 条件期望组成
利用Copula函数分别构造水布垭断面和水–隔区间，水布垭断面和隔–高区间最大3天洪量的二维联
(a)
(b)
(c)
(d)
Figure 3. Occurrence probability of equivalent frequency floods for each region and design site
图3. 各分区洪水同频率的可能性分析
Table 4. Regional 3-day flood volume proportions for three methods at Gaobazhou site (%)
表4. 三种高坝洲断面设计洪水地区组成的各分区3天洪量比例 (%)
合分布，通过条件分布函数，拟定条件期望地区组成。即水布垭断面为设计值，水–隔区间和隔–高区间为相应的条件期望值，结果列于表4。可以看出，同频率组成中相应的隔–高区间洪量要小于条件期望组成中的期望洪量，且在稀遇洪水的设计中，差别更大。然而，隔–高区间洪水又是不可控的，一定程度上增加了高坝洲断面的防洪风险。事实上，洪水地区组成本质上是一种洪水分配方案，在总量一定的条件下，某一区域洪水所占比重增大，势必导致另一些区域所占比例减小。同频率地区组成的主观假定可能会造成一方偏大，另一方偏小的结果，对防洪安全不利。本文所提出的条件期望组成考虑了洪水的空间相关性，一定程度上避免了同频率组成的这种任意性。
3) 最可能组成
基于Copula函数构造水布垭断面、水–隔区间和隔–高区间最大3天洪量的三维联合分布，以联合概率密度函数最大为目标函数，从最可能发生的角度拟定最可能组成。表4列出了高坝洲断面设计洪水的最可能组成的计算结果。结果显示，当高坝洲坝址发生频率5%以上洪水时，水布垭断面所占的洪量比例为[54.0%，59.2%]，小于其面积比69.4%；水–隔区间所占洪量比例为[30.0%，33.5%]，大于其面积比22.8%；隔–高区间所占洪量比例为[10.8%，12.5%]，大于其面积比7.8%。此外，高坝洲断面发生较大洪水主要是由于水布垭断面发生洪水所致，且随着设计洪水量级的增大，水布垭断面所占比例均逐渐降低，水–隔区间和隔–高区间所占洪量比例逐渐增加。对于高坝洲断面千年一遇的设计洪水，最可能的组成是水布垭断面、水–隔区间和隔–高区间洪水重现期分别为438年、457年和813年，而非同频率组成对应的1000年、384年和77年。可见，同频率组成并不是最可能出现的地区组成。
5.3. 水布垭–隔河岩梯级水库对高坝洲断面设计洪水的影响
得到高坝洲断面设计洪水的地区组成结果后，就可以定量估算水布垭–隔河岩梯级水库对高坝洲断面设计洪水的影响。以水布垭断面、水–隔区间和隔–高区间分配到的相应洪量为控制，按1997年典型洪水过程线同倍比放大得到各自的设计洪水过程线，输入到水布垭水库和隔河岩水库按各自的调洪规则[10] 进行调洪演算。水布垭–隔河岩区间考虑了河道洪水演算，由于隔河岩–高坝洲区间洪水传播时间较短，隔河岩水库下泄洪水过程与隔–高区间洪水过程线叠加后，未进行河道演算，直接得到高坝洲断面受水布垭–隔河岩梯级水库影响的洪水过程线。不同设计频率的高坝洲断面天然洪峰流量和三种地区组成推求的受水布垭–隔河岩梯级水库影响的洪峰流量见表5。
Table 5. Design peak discharges at Gaobazhou site under the impact of cascade reservoirs (m3/s)
表5. 梯级水库调蓄影响的高坝洲断面洪峰流量对比 (m3/s)
由表5可知，与高坝洲断面天然洪峰流量相比，受水布垭–隔河岩梯级水库调蓄影响的洪峰流量有较大程度的削减。对于千年一遇的高坝洲断面设计洪水，高坝洲断面洪峰流量由天然情况的24,300 m3/s减小为[13,980 m3/s，14,870 m3/s]，削减[38.8%，42.5%]，平均40.2%；对于百年一遇而言，由天然情况的18,800 m3/s减小为[11,800 m3/s，12,240 m3/s]，削减[34.9%，36.8%]，平均35.6%。三种地区组成计算的结果相差不大，与基于Copula函数的两种组成相比，同频率组成的结果系统偏小，可能会增加高坝洲断面的防洪风险。造成这一结果的原因可能是同频率组成指定水布垭断面出现与高坝洲断面同频率的洪水，水–隔区间和隔–高区间则分别出现相应的洪水，而真实情况并非如此。基于Copula函数的条件期望组成和最可能组成则充分考虑了水布垭断面、水–隔区间和隔–高区间洪水的随机性和空间相关性，具有较强的统计基础，更加符合客观实际。
6. 结论
针对目前工程设计中常用的同频率地区组成具有一定的人为性，统计理论基础不足的问题，利用Copula函数推求了两种具有统计基础的梯级水库设计洪水地区组成，即条件期望组成和最可能组成。以清江流域高坝洲坝址断面的设计洪水地区组成为例，对以上三种地区组成进行比较研究，分析了上游水布垭–隔河岩梯级水库的调洪对高坝洲断面设计洪水的影响，可以得到如下结论：
1) 受上游梯级水库的调蓄影响，高坝洲断面洪峰流量有较大程度削减，千年一遇设计洪水的削峰率为40%左右。
2) 同频率组成有一定的代表性，但它既不是最可能出现的，也不是最恶劣的地区组成，其设计洪水成果系统偏小，可能会增加设计断面的防洪风险。
3) 基于Copula函数的条件期望组成和最可能组成，充分考虑了各分区洪水的随机性和空间相关性，具有较强的统计基础，更加合理。
基金项目
国家自然科学基金(51190094, 51079100)。
NOTES
作者简介：刘章君(1991-)，男，江西吉安人，博士研究生，主要从事水文分析与计算等方面的研究。