1. 引言

设计洪水是指水利水电工程规划、设计、施工中所指定的各种设计标准的洪水。20世纪50年代以来，我国开始大规模水利水电工程建设，并开展了水利水电工程设计洪水计算方法及应用研究。设计洪水标准和计算方法经历了从历史洪水资料加成法逐步过渡到频率分析计算方法，以及引入可能最大洪水的发展过程。1979年水利部和电力工业部颁发了《水利水电工程设计洪水计算规范SDJ22-79(试行)》，使我国设计洪水有了统一的标准；1993年编制了《水利水电工程设计洪水计算规范SL44-93》；2006年又修编了新的《水利水电工程设计洪水计算规范SL44-2006》 [1] 。新规范一方面保留了符合我国情况、切实可行的条文，另一方面以多年来我国行之有效的实用方法和科技成果为背景进行补充和发展，反映了科学性、先进性和实用性。但规范主要是针对天然流域和单一水库的设计洪水计算，是确保大坝防洪安全的设计值，可称为“水库建设期设计洪水”。

至2016年底，我国已建成大、中、小型水库98,460座，总库容8967亿m³。大部分河流已经或即将形成梯级水库群格局，明显改变了河川径流及洪水的时程分配过程。水库管理运用的主要任务是保证防洪安全并充分发挥水库兴利效益。防洪法规定水库在汛期不得擅自在汛期限制水位以上蓄水。我国水库调度管理沿用建设期的年最大设计洪水及特征值，要求水库每年整个汛期按年最大设计洪水设防，常造成汛期一出现洪水就被迫弃水，而到汛末又往往蓄不满。近年开展的分期设计洪水及分期汛限水位，汛期运行水位动态控制(中小洪水调度)，水库提前蓄水等研究和实践，取得了一些重大进展并产生了巨大的经济和社会效益。

设计洪水计算如何适应这种下垫面及河道汇流条件的改变，充分考虑上游水库工程的影响和水库群之间的相互补偿作用等，是水库运行管理实际工作中遇到的难题，也是亟待解决的水文科学问题。因此，研究探讨水库运用期设计洪水理论和方法，定量估算上游水库群调蓄对下游断面设计洪水和防洪安全标准的影响，推求水库运用期汛限水位，为水库科学调度方案提供理论依据与技术支撑，进一步补充修订我国《水利水电工程设计洪水计算规范》提供理论基础和科学依据，具有重大的理论价值和现实意义。

2. 设计洪水研究进展

郭生练等 [2] 系统总结了国内外洪水频率分析方法，包括抽样方法、分布线型、经验频率公式、参数估计方法、设计洪水过程线、历史洪水、区域洪水频率、PMP/PMF、分期和梯级水库设计洪水等方面的研究进展，重点阐述了中国近30年来的主要研究成果，并讨论了存在的主要问题和不足。下面重点综述讨论分期设计洪水、动态设计洪水，考虑上游水库影响设计洪水的研究进展及存在问题。

2.1. 分期设计洪水

分期设计洪水的核心内容包括汛期分期与分期设计洪水计算两个方面。我国地处亚洲季风区气候，河流年内汛期变化具有确定性、随机性、过渡性等特性。现有的汛期分期方法很多，可以概括分为成因分析法和数理统计法两种，分别属于定性和定量分析方法。数理统计法包括模糊分析法、分形分析法、变点分析法、系统聚类法、矢量统计法、相对频率法和圆形分布法等。成因分析法属于定性分析方法，早期应用较为广泛，但由于缺乏科学的数学计算依据，使得分期带有一定的不确定性；数理统计定量分析方法弥补了这方面的不足。虽然汛期分期的方法很多，各种方法都有各自的优缺点，得到的分期结论会有所差异，至今尚未形成公认的标准方法 [3] 。

我国现行的分期设计洪水计算方法基于分期最大值系列，采用单变量的P-III型分布进行描述，假定分期设计洪水的频率等于设计标准对应的重现期的倒数，求得各设计频率的分期设计洪水。大量研究表明：现行方法计算的汛期分期设计洪水系统偏小，通过应用它来提高水库兴利效益，是在冒降低水库防洪标准的风险。分期设计洪水既要符合防洪标准，又要能反映洪水的季节性规律。国内外学者在探索同时满足这两个条件的分期设计洪水计算技术方面做了一些尝试，主要包括全概率公式法和联合分布法两类 [2] 。宋松柏等 [4] 根据重现期的定义，运用数理统计法推导了独立同分布单变量和多变量水文事件重现期的计算公式。在此基础上，结合现有的3 种分期洪水概率模型，推导了年最大洪水重现期与分期最大洪水有关重现期计算模型和关系式。

2.2. 动态设计洪水

2010年，王善序提出了水库运用期动态设计洪水的理论和方法 [5] 。他采用随机过程理论描述径流年内不平稳随机过程，构建相应的水库运用动态设计洪水预报概率模型。他认为：传统的一个防洪标准对应一个设计洪水，这只有当设计涉及的是多年径流过程时才正确。在水库运用面临的不平穏过程，洪灾风险率、动态设计洪水都是时变的并与入库水情有关。这时，设计洪水也一定是时变的，并与入库水情有关。水库在汛期只有按动态设计洪水预留防洪库容(设置汛期运行水位)，才能使水库顺应自然径流演变规律、合理处理防洪与兴利矛盾、既保证防洪安全，又能提高兴利。

如何估算动态设计洪水，王善序研究了在汛期任一决策时刻估计未来时段设计(洪)水量，如三峡水库可取7~15 d的“时段水量”，根据历史观测资料推求时段设计水量，指导水库调度运行 [6] 。虽然该法具有一定的理论基础，但其实质是把整个汛期分成更多的未来时段(即更多的分期)，这样推求出来的动态设计洪水，无法满足现有防洪标准的要求，其风险要远大于分期设计洪水。

2.3. 考虑上游水库影响的设计洪水

梯级水库设计洪水，从暴雨上说，是要解决暴雨的地区分布即空间组合问题，并提出用同频率概念控制法和典型年法推求梯级水库的PMP；从洪水来说，是要解决洪水地区组成问题 [7] 。我国《水利水电工程设计洪水计算规范》推荐，洪水地区组成一般采用典型年和同频率组成两种方法。对于单库，一般多考虑同频率组成及典型年组成；梯级水库则大多采用典型年组成或通过自下而上逐级分析的方法拟定，即各级设计洪量可以采用不同的典型洪水进行分配，也可混合采用典型年法及同频率组成法分配洪量。

国外相关的研究很少，主要原因是西方发达国家的水利水电工程到20世纪70年代基本建成并投入运行。从查到的文献来看，采用的方法比国内现行方法简单。例如，美国陆军工程兵团，分析了不同的防洪规划方案和防洪措施对下游防洪区最大流量的影响，进而确定最优防洪控制策略时，通常采用下述两种方法推求受防洪工程影响后下游断面的设计洪水 [8] 。

① 当有流量资料时，先计算出防洪断面天然情况的设计洪水，根据防洪断面的防洪标准，选择一个大洪水典型，以防洪标准相应的洪峰流量之比为放大系数，线性放大各控制断面的洪水过程线。对防洪工程所在断面的洪水过程线，考虑其调洪效应后得到出流过程线，再与区间的洪水过程线组合(通过流量演算)，求得防洪断面调洪后的洪水过程线和最大流量，并假定天然情况最大流量与相应的调洪后的最大流量是同频率的，由此就可以求得调洪后不同频率的最大流量。

② 由暴雨资料推求防洪断面以上流域内各控制点的设计点暴雨，在暴雨与洪水同频的假定下，通过暴雨时面深关系和流域模型，将各种频率的暴雨转化成同频率的洪水过程线，然后采用与①相同的方法推求出防洪断面受调洪影响后不同频率的最大流量。不难看出，方法①只相当于国内现行方法中考虑典型年洪水地区组成一种形式。而方法②不仅假定暴雨与洪水同频，而且假定流域内各控制点的暴雨都同频，显然是不合理的。

李天元等 [9] 全面综述分析了近十几年来洪水地区组成方法：包括① 地区组成法；② 频率组合法；③ 随机模拟法；④ JC法；⑤ Copula函数法。是带有明显随机性的，科学地描述洪水地区组成规律，需要给出各地区洪水的联合概率分布函数。过去由于受到多维联合分布方法的限制，实际应用中只能寻求特定条件下的近似计算，如以上提到的同频率法和典型年法等。近年来随着Copula函数在水文领域的成功应用 [10] ，使得其推求设计洪水的地区组成成为可能。闫宝伟等 [11] 应用Copula函数构造了上游断面与区间洪水的联合分布，提出和推导了两种有代表性的设计洪水地区组成，即条件期望组成和最可能组成。推求的最可能组成考虑了地区间洪水的空间相关性，更加符合客观实际，为设计洪水地区组成分析计算提供了一条新的途径。李天元等 [12] 以Copula函数理论为基础，构造了水库断面洪量与区间洪量的联合分布，推求了条件概率函数的显式表达式，提出了基于Copula函数的改进离散求和法，通过直接对条件概率曲线进行离散，克服了《规范》中离散求和法需要进行变量独立性转换的问题。刘章君等 [13] 利用Copula函数建立了各分区洪水的联合分布，基于联合概率密度最大原则，推导得到最可能地区组成法的计算通式，并用来推求梯级水库下游断面的设计洪水。

上述方法在工程实践中都有应用，但都存在较明显的缺陷。规范中的地区组成法虽然简便易行，但人为不确定性较大；频率组合法需要对分区的频率曲线进行离散求和，在独立性转换中难免出现数据失真；对于复杂的梯级水库群，随机模拟法难以保持各分区的洪水涨落特性；JC法在分配洪量计算的同时，所得的风险或概率值可作为决策者判定其成果可靠性的依据，但是目前国内外在可接受风险方面并无相关参考，这将给决策带来困难。基于Copula函数的梯级水库设计洪水方法具有较强的统计基础，且所得设计结果客观合理，但其应用也存在一些问题，如当梯级水库较多时，基于Copula函数的最可能地区组成法需要应用高维Copula函数，计算难度较大。此外，当上下游水库距离较远时，现有的马斯京根洪水演进方法的精度难以满足要求。

3. 水库运用期设计洪水的科学问题

随着一大批大型水电站开工建设和建成投运，以三峡工程为骨干的长江上游干支流梯级水库群已形成规模(如图1所示)，已建或在建总调节库容1000亿m³以上，总防洪库容超过500亿m³ [14] 。长江上游水库群的防洪与兴利矛盾表现十分突出，尤其是在汛期和汛末蓄水期，水库之间以及同中下游地区水资源利用矛盾等问题更加严峻，常常出现汛期有水不能蓄、不敢蓄、汛后蓄不满的被动局面。水库汛期运行水位动态控制和提前蓄水，虽可增加发电量，提高水库蓄满率和汛后水资源供给量，但汛期运行水位太高会占用一部分的防洪库容，需要定量估算可能增加的防洪风险。

水库群建成投入运行后，洪峰经过水库调蓄后会变小变缓，不同时段(特别是短时段，如1、3、7、15 d的设计洪量)的设计洪水也会变小。由于下游水库设计洪水变小，在确保防洪安全(包括枢纽本身和防洪保护对象)的基础上，水库运用期的汛限水位可以适当上浮。关键是如何科学分析计算受上游水库调蓄影响下的下游水库设计洪水，有了科学合理的设计洪水成果，才可能分析确定水库汛期运行水位(因为设计洪水减小)的浮动上限值。因此，如何推求水库运用期的设计洪水及其汛限水位，如何协调各个水库防洪、发电、航运、供水等调度目标以及水库群之间的蓄泄关系，是水库防洪调度和洪水资源化所面临的重大科学问题和关键技术难题。下面以三峡水库为例，水库运用期的设计洪水的科学问题包括：

1) 三峡水库建设期千年一遇设计洪水为98,000 m³/s，在同样的来水条件下，经过长江上游水库群的调蓄，三峡水库运用期千年一遇设计洪水(洪峰、不同时段洪量)是多少？

2) 根据长江流域大洪水年(如1954年、1998年)观测资料和各水库调度规程，经过现有上游水库群调蓄和演算后的三峡入库洪水过程及洪峰洪量是多少？

3) 如何寻找复杂梯级水库群的洪水遭遇和地区组成规律？如何计算三峡水库运用期的设计洪水？

4) 如何推求水库运用期的汛限水位？能否在实际调度运行中替代原有的汛限水位？在确保防洪安全标准不变的前提下，可增加的综合利用效益有多大？

4. 水库运用期的设计洪水研究内容和方法

我国现有《规范》制定了水库建设期的设计洪水，无法回答上述科学问题。因此有必要在综述分析比较现有方法的基础上，探讨理论依据充分、方案切实可行、易于推广应用的水库运用期的设计洪水计算理论和方法，指导水库调度运行和科学管理。

4.1. 研究内容

梯级水库群调度运行期的设计洪水计算，一般都要涉及到推求受上游水库调蓄影响后的设计洪水。由上、下两个水库组成的梯级水库是最常见的，具有一定的代表性，因为多级水库可以看成是两级水库的各种组合。当干支流串联水库并联到一起，即形成了混联梯级水库群，因此，本项目的研究对象可归纳为以下四种类型：

1) 两串联水库均不承担下游防护对象的防洪任务。如图2(a)所示，B水库的洪水是经A水库调洪后的下泄洪水与区间Y的洪水组合而形成的，所以在进行B水库的防洪设计时，就需要推求B受A水库调洪影响后的设计洪水。

2) 两串联水库下游有防洪对象。如图2(b)所示，如果所要设计的工程是A水库工程，为研究A、B两个梯级水库对防护对象C的防洪效果，就需要推求C断面受上游A、B两水库调洪综合影响后的设计洪水；如果所要设计的工程是B水库工程，就需要推求B受A水库调洪影响的设计洪水，同时还要推求C断面受A、B两水库调洪共同影响的设计洪水。

3) 两串联水库之间有防洪对象。如图2(c)所示，在设计A水库时，为研究A水库对防护对象C的防洪作用，需要推求C断面受A水库调洪影响后的设计洪水。在设计B水库时，就需要推求B受A水库调洪影响后的设计洪水。

4) 梯级水库群下游有防洪对象。如图2(d)所示，当两条河流上的梯级水库共同承担C断面的防洪任务时，需要推求C断面受上游干流A~B梯级水库和支流A^’~B^’梯级水库组成的混联水库群联合调蓄影响后的设计洪水。

结合具体的研究对象，项目的主要研究内容可分为以下几个部分。

1) 研究流域洪水遭遇特性与地区组成规律。从气候成因角度，结合统计学的方法，分析流域暴雨洪水特性、以及上下游范围内暴雨洪水的各种成因、地域分布、强度、时间的分布；分析上下游洪水遭遇特性，以及遭遇时洪峰、洪量组成规律，全面揭示流域洪水的地区组成规律。

2) 研究梯级水库设计洪水地区组成方案。要分析梯级水库调蓄作用对下游断面洪水的影响效果，首先要拟定下游防洪断面以上水库断面及区间的洪水地区组成方案。研究采用规范规定的地区组成方法确定梯级水库各断面洪量分配方案，重点探讨研究具有统计基础、切实可行的地区组成新方案，多方案相互比较分析。

3) 研究受梯级水库群影响的设计洪水计算新方法。结合拟定的设计洪水地区组成新方案，根据梯级水库的防洪调度规则，分析得到下游防洪断面受梯级调蓄影响后的设计洪水，分析梯级水库对流域洪水的削减作用，以及对下游断面防洪标准和防洪风险的影响。

4) 研究梯级水库调蓄对下游已建水库运行管理的影响。根据计算得到的受梯级水库调蓄影响后的设计洪水，进一步分析梯级水库对河川径流及洪水时空分配过程的影响，为下游已建水库的防洪调度和运行管理提供更全面准确的信息，为流域梯级规划提供技术支撑。

5) 根据水库运用期的设计洪水，通过调洪演算推求水库运用期的汛限水位，分析防洪风险和兴利效益的增量，指导梯级水库群调度方案编制和运行。

4.2. 研究方法

基于还原后的天然流量系列，采用最可能地区组成法构建天然状态下梯级水库群设计洪水计算模块(见图3)，分别推求各水库控制断面的设计洪水及特征值。

1) 应用基于Copula函数的多变量联合分布理论和方法，用来处理、分析、研究流域梯级水库断面与各区间洪水的时空相关性及组成规律。

2) 通过对各种方法计算结果的分析比较，提出理论依据充分、切实可行的水库运用期设计洪水的新理论和新方法，并论证方法的合理性和实用价值。

4.3. 关键技术

研究的关键是寻找洪水遭遇特性与地区组成规律。如何从气候成因角度，结合统计学的方法，分析流域暴雨洪水特性、以及上下游范围内暴雨洪水的各种成因、地域分布、强度、时间的分布；分析上下游洪水遭遇特性，以及遭遇时洪峰、洪量组成规律，全面揭示流域洪水的地区组成规律。探讨梯级水库设计洪水地区组成的不同方案，分析梯级水库调蓄作用对下游断面洪水的影响效果，重点探讨研究具有统计基础、切实可行的地区组成新方案。下面给出最可能地区组成法的研究思路和方法：

1) 梯级水库设计洪水的最可能地区组成法

梯级水库如图4所示，C为设计断面，其上游有n个水库 ${\text{A}}_{\text{1}},{\text{A}}_{\text{2}},\cdots ,{\text{A}}_n{}_{-\text{1}},{\text{A}}_n$ ；n个区间流域 ${\text{B}}_{\text{1}},{\text{B}}_{\text{2}},\cdots ,{\text{B}}_n{}_{-\text{1}},{\text{B}}_n$ 。随机变量X、Y_i和Z分别表示水库A₁、区间流域B_i和断面C的天然来水量，取值依次为x、y_i和z ( $i=1,2,\cdots ,n$ )。

受上游 ${\text{A}}_{\text{1}}{\text{-A}}_{\text{2}}\text{-}\cdots {\text{-A}}_n{}_{-\text{1}}{\text{-A}}_n$ 梯级水库的影响，分析断面C设计洪水的地区组成需要研究天然情况下水库A₁断面和n个区间 ${\text{B}}_{\text{1}},{\text{B}}_{\text{2}},\cdots ,{\text{B}}_n{}_{-\text{1}},{\text{B}}_n$ 共(n + 1)个部分洪水的组合。由于河网调节等因素的影响，往往难以推求设计洪峰流量的地区组成，且对调洪能力大的水库，洪量起主要作用。因此，通常将断面C某一设计频率p的时段洪量z_p分配给上游(n + 1)个组成部分，以研究梯级水库的调洪作用。由水量平衡原理得：

$x+{\displaystyle \underset{i=1}{\overset{n}{\sum }}{y}_i}=z_p$ (1)

式中：z_p为断面C的天然设计洪量，x、y_i分别为水库A₁、区间B_i相应的天然洪量。

设计洪水的地区组成本质上是给定断面C的设计洪量z_p，在满足式(1)约束条件下分配z_p，得到组合

$\left(x,y_1,y_2,\cdots ,y_n{}_{-1},y_n\right)$ 。得到洪量分配结果后，可以从实际系列中选择有代表性的典型年，放大该典型年各分

区的洪水过程线可得各分区相应的设计洪水过程线，然后输入到 ${\text{A}}_{\text{1}}{\text{-A}}_{\text{2}}\text{-}\cdots {\text{-A}}_n{}_{-\text{1}}{\text{-A}}_n$ 梯级水库系统进行调洪演算，就可以推求出同一频率p断面C受上游梯级水库调蓄影响的设计洪水值。

理论上，满足式(1)约束的洪量组合 $\left(x,y_1,y_2,\cdots ,y_n{}_{-1},y_n\right)$ 有无数种可能，但不同组合推求的断面C受上游梯

级水库影响的设计洪水值不同。因此，如何选择合理的地区组成方法来计算洪量组合至关重要。同频率地区组成法是目前使用最为广泛的方法，但该法只是一种人为假设，是否符合流域洪水地区组成规律，要视分区与设计断面洪水的相关性密切程度而定。此外，随水库数量的增加，需要拟定的地区组成方案数呈指数倍增长，计算工作量将急剧增加，难以满足复杂梯级水库设计洪水计算的实际需要。

在工程设计中，人们通过对实际发生洪水的时空特性规律分析，通常关心最可能发生且对下游防洪较不利的地区组成形式。因此，为了避免组合的任意性和盲目性，从发生可能性最大的角度，推导梯级水库设计洪水的最可能地区组成计算通式。

不同洪量组合发生的相对可能性大小，可以用X和Y_i ( $i=1,2,\cdots ,n$ )的联合概率密度函数值 $f\left(x,y_1,y_2,\cdots ,y_n{}_{-1},y_n\right)$ 大小来度量。联合概率密度函数值越大，表明该种地区组成发生的可能性越大。欲得到最可能地区组成，即为求解 $f\left(x,y_1,y_2,\cdots ,y_n{}_{-1},y_n\right)$ 在满足式(1)条件下的最大值，即

$\begin{array}{l}\max f\left(x,y_1,y_2,\cdots ,y_{n-1},y_n\right)\\ \text{s}\text{.t}\text{.}x+{\displaystyle \underset{i=1}{\overset{n}{\sum }}{y}_i}=z_p\end{array}\}$ (2)

X和Y_i ( $i=1,2,\cdots ,n$ )的联合分布函数用 $F\left(x,y_1,y_2,\cdots ,y_{n-1},y_n\right)$ 表示，表达式如下：

$F\left(x,y_1,y_2,\cdots ,y_{n-1},y_n\right)=P\left(X\le x,Y_1\le y_1,\cdots ,Y_n\le y_n\right)$ (3)

相应的概率密度函数为：

$f\left(x,y_1,y_2,\cdots ,y_{n-1},y_n\right)={\partial }^{n+1}F\left(x,y_1,y_2,\cdots ,y_{n-1},y_n\right)/\partial x\partial y_1\cdots \partial y_n$ (4)

我国一般假定X和Y_i ( $i=1,2,\cdots ,n$ )均服从P-III型分布 [1] ，且各变量之间还存在一定的相关性。若用传统办法来确定(n + 1)个随机变量的联合概率分布函数和相应的密度函数解析表达式，是非常困难的。近年来Copula函数理论的应用和日趋完善使之成为可能，它能够灵活地构造具有任意边缘分布的多变量联合分布 [10] 。

2) 基于Copula函数推导最可能地区组成计算通式

假设随机变量 $X_i\left(i=1,2,\dots ,n\right)$ 的边缘分布函数分别为 $F_{X_i}\left(x_i\right)=P_{X_i}\left(X_i\le x_i\right)$ ，其中n为随机变量的个数， $x_i$ 为随机变量 $X_i$ 的取值， $H\left(x_1,x_2,\cdots ,x_n\right)$ 为随机变量 $X_1,X_2,\cdots ,X_n$ 的联合分布函数，则存在唯一的Copula函数使得 [15] ：

$H\left(x_1,x_2,\cdots ,x_n\right)=C\left(F_{\text{1}}\left(x_1\right),F_{\text{2}}\left(x_2\right),\cdots ,F_n\left(x_n\right)\right)=C\left(u_1,u_2,\cdots ,u_n\right)$ (5)

式中： $u_i=F_{X_i}\left(x_i\right)$ ，在区间[0,1]上连续。

应用Copula函数建模时，可分为两步：确定随机变量各自的边缘分布；根据各变量的边缘分布，选取适当的Copula函数，准确地描述边缘分布间的相关结构。

借助Copula函数，联合分布函数 $F\left(x,y_1,y_2,\cdots ,y_{n-1},y_n\right)$ 可以表示为：

$F\left(x,y_1,y_2,\cdots ,y_{n-1},y_n\right)=C\left(u,v_1,v_2,\cdots ,v_{n-1},v_n\right)$ (6)

式中： $C\left(u,v_1,v_2,\cdots ,v_{n-1},v_n\right)$ 为Copula函数； $u=F_X\left(x\right)$ 、 $v_i=F_{Y_i}\left(y_i\right)$ 分别为水库A₁、区间B_i洪水的边缘分布函数。

相应的联合概率密度函数：

$f\left(x,y_1,y_2,\cdots ,y_{n-1},y_n\right)=c\left(u,v_1,v_2,\cdots ,v_{n-1},v_n\right)f_X\left(x\right){\displaystyle \underset{i=1}{\overset{n}{\prod }}{f}_{Y_i}\left(y_i\right)}$ (7)

式中： $c\left(u,v_1,v_2,\cdots ,v_{n-1},v_n\right)={\partial }^{n+1}C\left(u,v_1,v_2,\cdots ,v_{n-1},v_n\right)/\partial u\partial v_1\partial v_2\cdots \partial v_{n-1}\partial v_n$ 为Copula函数的密度函数； $f_X$ 、 $f_{Y_i}$ 分别为X和Y_i ( $i=1,2,\cdots ,n$ )的概率密度函数。

将 $x=z_p-{\displaystyle \underset{i=1}{\overset{n}{\sum }}{y}_i}$ 代入 $f\left(x,y_1,y_2,\cdots ,y_{n-1},y_n\right)$ 得到：

$\begin{array}{l}f\left(z_p-{\displaystyle \underset{i=1}{\overset{n}{\sum }}{y}_i},y_1,y_2,\cdots ,y_{n-1},y_n\right)\\ =c\left(F_X\left(z_p-{\displaystyle \underset{i=1}{\overset{n}{\sum }}{y}_i}\right),F_{Y_1}\left(y_1\right),F_{Y_2}\left(y_2\right),\cdots ,F_{Y_{n-1}}\left(y_{n-1}\right),F_{Y_n}\left(y_n\right)\right)f_X\left(z_p-{\displaystyle \underset{i=1}{\overset{n}{\sum }}{y}_i}\right){\displaystyle \underset{i=1}{\overset{n}{\prod }}{f}_{Y_i}\left(y_i\right)}\end{array}$ (8)

欲使 $f\left(z_p-{\displaystyle \underset{i=1}{\overset{n}{\sum }}{y}_i},y_1,y_2,\cdots ,y_{n-1},y_n\right)$ 取最大值，需满足：

$\partial f\left(z_p-{\displaystyle \underset{i=1}{\overset{n}{\sum }}{y}_i},y_1,y_2,\cdots ,y_{n-1},y_n\right)/\partial y_i=0$ $\left(i=1,2,\cdots ,n\right)$ (9)

刘章君等 [13] 通过推导归纳得出，由n个水库组成的梯级防洪系统设计洪水最可能地区组成的方程组通式为：

$\{\begin{array}{l}{c}_1{f}_X\left(x\right)-c_2{f}_{Y_1}\left(y_1\right)+c\cdot \left[\frac{{f^{\prime }}_X\left(x\right)}{f_X\left(x\right)}-\frac{{f^{\prime }}_{Y_1}\left(y_1\right)}{f_{Y_1}\left(y_1\right)}\right]=0\\ c_1{f}_X\left(x\right)-c_3{f}_{Y_2}\left(y_2\right)+c\cdot \left[\frac{{f^{\prime }}_X\left(x\right)}{f_X\left(x\right)}-\frac{{f^{\prime }}_{Y_2}\left(y_2\right)}{f_{Y_2}\left(y_2\right)}\right]=0\\ ⋮\\ c_1{f}_X\left(x\right)-c_n{f}_{Y_{n-1}}\left(y_{n-1}\right)+c\cdot \left[\frac{{f^{\prime }}_X\left(x\right)}{f_X\left(x\right)}-\frac{{f^{\prime }}_{Y_{n-1}}\left(y_{n-1}\right)}{f_{Y_{n-1}}\left(y_{n-1}\right)}\right]=0\\ c_1{f}_X\left(x\right)-c_{n+1}{f}_{Y_n}\left(y_n\right)+c\cdot \left[\frac{{f^{\prime }}_X\left(x\right)}{f_X\left(x\right)}-\frac{{f^{\prime }}_{Y_n}\left(y_n\right)}{f_{Y_n}\left(y_n\right)}\right]=0\\ x+{\displaystyle \underset{i=1}{\overset{n}{\sum }}{y}_i}-z_p=0\end{array}$ (10)

式中： $c=c\left(u,v_1,v_2,\cdots ,v_{n-1},v_n\right)$ ， $c_1=\partial c/\partial u$ ， $c_{i+1}=\partial c/\partial v_i$ ； ${f^{\prime }}_X,{f^{\prime }}_{Y_i}\left(i=1,2,\cdots ,n\right)$ 分别为相应概率密度函数的导函数。

若X和Y_i ( $i=1,2,\cdots ,n$ )均服从P-III分布，则

${f^{\prime }}_X\left(x\right)=\frac{{\beta }_x^{{\alpha }_x}}{\Gamma \left({\alpha }_x\right)}{\left(x-{\gamma }_x\right)}^{{\alpha }_x-2}{\text{e}}^{-{\beta }_x\left(x-{\gamma }_x\right)}\left[\left({\alpha }_x-1\right)-{\beta }_x\left(x-{\gamma }_x\right)\right]$ (11)

${f^{\prime }}_{Y_i}\left(y_i\right)=\frac{{\beta }_{y_i}^{{\alpha }_{y_{{}_i}}}}{\Gamma \left({\alpha }_{y_i}\right)}{\left(y_i-{\gamma }_{y_i}\right)}^{{\alpha }_{y_{{}_i}}-2}{\text{e}}^{-{\beta }_{y_{{}_i}}\left(y_i-{\gamma }_{y_i}\right)}\left[\left({\alpha }_{y_i}-1\right)-{\beta }_{y_i}\left(y_i-{\gamma }_{y_i}\right)\right]$ (12)

故

$\frac{{f^{\prime }}_X\left(x\right)}{f_X\left(x\right)}=\frac{{\alpha }_x-1}{x-{\gamma }_x}-{\beta }_x$ ， $\frac{{f^{\prime }}_{Y_i}\left(y_i\right)}{f_{Y_i}\left(y_i\right)}=\frac{{\alpha }_{y_i}-1}{y_i-{\gamma }_{y_i}}-{\beta }_{y_i}$ (13)

式中： ${\alpha }_x$ 、 ${\beta }_x$ 、 ${\gamma }_x$ 和 ${\alpha }_{y_i}$ 、 ${\beta }_{y_i}$ 、 ${\gamma }_{y_i}$ 分别为变量X和Y_i形状、尺度和位置参数。将式(13)代入式(10)，则可表述为：

$\{\begin{array}{l}{c}_1{f}_X\left(x\right)-c_2{f}_{Y_1}\left(y_1\right)+c\cdot \left[\frac{{\alpha }_x-1}{x-{\gamma }_x}-\frac{{\alpha }_{y_1}-1}{y_1-{\gamma }_{y_1}}+{\beta }_{y_1}-{\beta }_x\right]=0\\ c_1{f}_X\left(x\right)-c_3{f}_{Y_2}\left(y_2\right)+c\cdot \left[\frac{{\alpha }_x-1}{x-{\gamma }_x}-\frac{{\alpha }_{y_2}-1}{y_2-{\gamma }_{y_2}}+{\beta }_{y_2}-{\beta }_x\right]=0\\ ⋮\\ c_1{f}_X\left(x\right)-c_n{f}_{Y_{n-1}}\left(y_{n-1}\right)+c\cdot \left[\frac{{\alpha }_x-1}{x-{\gamma }_x}-\frac{{\alpha }_{y_{n-1}}-1}{y_{n-1}-{\gamma }_{y_{n-1}}}+{\beta }_{y_{n-1}}-{\beta }_x\right]=0\\ c_1{f}_X\left(x\right)-c_{n+1}{f}_{Y_n}\left(y_n\right)+c\cdot \left[\frac{{\alpha }_x-1}{x-{\gamma }_x}-\frac{{\alpha }_{y_n}-1}{y_n-{\gamma }_{y_n}}+{\beta }_{y_n}-{\beta }_x\right]=0\\ x+{\displaystyle \underset{i=1}{\overset{n}{\sum }}{y}_i}-z_p=0\end{array}$ (14)

非线性方程组式(14)即为基于Copula函数推求的最可能地区组成法应满足的计算通式。式(14)分别有(n+1)个未知数 $x,y_1,y_2,\cdots ,y_n{}_{-1},y_n$ 和(n+1)个方程，根据问题的实际意义， $f\left(x,y_1,y_2,\cdots ,y_n{}_{-1},y_n\right)$ 的最大值客观上存在且唯一，因此该方程组必定有唯一的解。显然，其最大值不会在边界取得，而是在定义域内部，即 $0<x,y_1,y_2,\cdots ,y_n{}_{-1},y_n<z_p$ 取得，因此由偏导数为零求解方程组得到的驻点即为最大值点。

由于问题的复杂性，采用数值方法计算非线性方程组(14)的近似解。牛顿法是目前求解非线性方程组最为行之有效的一种算法，具有收敛速度快、算法稳定性强等优点。可将以各分区控制面积占设计断面的比例分配洪量z_p得到的结果作为初始解，采用牛顿迭代法进行迭代求解，得到水库A₁断面、区间流域B₁、B₂…B_n-1、B_n

洪水的最可能地区组成 $\left(x^{\ast },y_1^{\ast },y_2^{\ast },\cdots ,y_{n-1}^{\ast },y_n^{\ast }\right)$ 。

5. 结语

水利水电工程兴建和气候环境变化，流域下垫面和天然河流均发生了巨大的变化，导致水文过程的非一致性。现有《规范》以天然流量系列推求的水库设计洪水及特征值，主要是为保证水库长期运用大坝及防洪安全、需要确定水库防洪能力，但无法满足已建水库运用期洪水资源化的调控需求。因此，考虑人类活动影响和适应变化环境的水文分析计算，推求水库运用期的设计洪水及相应的水位特征值，指导水库防洪调度和洪水资源优化，是水文学与水资源学科领域研究的前沿课题，国内外类似的研究报道甚少。因此研究受上游水库运行调蓄影响条件下水库运用期的设计洪水理论和方法，可望填补国内外现有设计洪水计算规范的空白。

现有《规范》的设计洪水，即“建设期设计洪水”，是为水库防洪能力和大坝安全设计服务的；而水库运用期的设计洪水，是为水库科学调度决策服务的。作为水库设计洪水的组成部分，两者皆不可或缺、不可替代。当预测水库流域可能发生稀遇大洪水或存在其他重大安全隐患时，按水库建设期年最大设计洪水及其汛限水位调度，确保大坝和防洪安全；其他时间可按水库运用期的设计洪水及其汛限水位调度。

基金项目

国家重点研发计划(2016YFC0402206)，国家自然科学重点基金(51539009)资助项目。

参考文献

参考文献