程 旭，馬細霞，肖 遙，2，王倩麗

（1.鄭州大學水利科學與工程學院，鄭州450001；2.中水東北勘測設計研究有限責任公司，長春130021）

0 引言

設計洪水是水利工程規(guī)劃設計、防洪策略和水資源管理的重要水文數據，由流量資料推求設計洪水是設計洪水計算的主要方法［1，2］。傳統上，通過擬合各歷史洪水要素觀測值概率分布函數以估計其不同重現期下的單變量設計洪水值，以其峰值、洪量和持續(xù)時間這3個內部存在物理聯系而并不相互獨立的隨機變量分析洪水過程，對于一個完整的洪水事件，僅依靠單變量概率分布分析存在一定的局限性，不能完整描述洪水過程［3，4］。因此要進行更完整的設計洪水頻率分析，就需要對各洪水要素進行多變量聯合統計分析［5，6］。設計洪水地區(qū)組成是流域開發(fā)方案設計以及上下游、干支流工程聯合調洪運行的重要依據［7，8］，目前，常用的地區(qū)組成方法有同頻率法和典型年法，該類方法主觀地考慮“對防洪不利”［9-11］。雖然便于使用和規(guī)劃，但更多是主觀臆斷而缺少客觀的判斷，不能充分反映各分區(qū)之間的空間相關性，也無法定量地描述各分區(qū)與設計斷面之間的相關程度［12］。Copula 函數是將變量聯合累積分布函數同變量邊緣累積分布函數連接起來的函數，能夠構建各洪水要素之間的聯合分布，也能夠構建各分區(qū)洪水之間的聯合分布［13-15］。

為求得同時考慮洪水要素相關性和空間相關性的設計洪水結果，本文以黃河中游河南段洪水地區(qū)組成復雜的花園口水文站上游為研究對象，選取12天洪量、洪峰流量兩個洪水要素，基于Copula 函數建立兩個洪水要素以及不同斷面之間洪水要素的聯合分布，并與遺傳算法相結合構建Copula-GA（Copula-Genetic Algorithms）算法，推求最有可能的設計洪水地區(qū)組成，以期為防洪預案制定提供參考依據。

1 研究區(qū)概況及數據來源

花園口水文站位于黃河中游，控制流域面積73 萬km2，占黃河流域總面積的92%［16］?；▓@口斷面上游干流河南段建有三門峽和小浪底兩座大型水庫，流域水系圖如圖1所示，研究區(qū)水系縱橫，呈不對稱分布，水患災害頻發(fā)，各量級洪水以上大洪水和上下同時大洪水為主占95.2%，即以全流域性洪水為主［17］。各防洪斷面洪水組成成分的確定對于水庫間防洪優(yōu)化聯調工作的開展起著至關重要的作用［18］?；▓@口斷面的洪水概率分布直接受到干支流水庫調節(jié)作用的影響［19］。所以，應與黃河中游各水系實際情況相結合，確定花園口洪水地區(qū)組成方案。本文將以花園口站過境洪量為控制標準，探求花園口以上，小浪底和相應區(qū)間洪水的相關關系及其地區(qū)組成。

圖1 研究區(qū)水系圖Fig.1 Water system diagram of the study area

本文基本資料為1950-2015年黃河干流河南段的花園口、小浪底以及小-花區(qū)間各站點的實測流量。1958年7月份洪水是花園口水文站記載以來的最大洪水，花園口站洪峰流量高達22 300 m3/s，該場洪水預見期短，對下游威脅較大，具有一定的代表性，對該場實測洪水過程進行插值，對洪峰所在時段取實測洪峰值，得到分辨率為6 h的洪水典型洪水過程。

2 研究方法

本文基于Copula 函數特性以及研究區(qū)洪水情況，建立花園口站洪量和洪峰流量的聯合分布，計算相應重現期的時段洪量和洪峰流量，并推求設計洪水過程線；采用Copula 函數構造各干支流洪水要素之間的聯合分布，并與遺傳算法結合構建Copula-GA算法，推算各干支流設計洪水的最可能地區(qū)組成。

2.1 兩變量重現期

選取時段洪量、洪峰流量兩個變量作為洪水特征要素，應用一組隨機變量（Qi，Wi）描述相互獨立的洪水事件系列｛I1，I2，...，In｝中的每一場洪水事件Ii，假設其聯合分布為F（q，w）。在設計洪水要素推求過程中，洪峰和洪量的表示方式為：

將以上單變量事件整合為聯合變量事件，表達方式為：

式中：Q和W分別表示洪水事件的洪峰流量和洪量；洪水事件Iq，w表示Q，W至少一個超過設定閾值。

若Q和W為年最大值，其單變量重現期為Tq、Tw［20］：

此時，Eq，w事件的重現期為Tq，w：

FQ（q）和FW（w）均大于或等于F（q，w），由（4）和（5）可得：

假定采用單變量方法得到的T年一遇設計洪峰qT和洪量wT，即變量Q和W分別取T年一遇的設計值，根據式（6）可得：

從上式可知，若要求Tq，w（qT，wT），則必定q和w分別大于或等于qT和wT。

2.2 建立洪峰流量和洪量的聯合分布

選取洪峰流量以及洪量兩個洪水要素，使用Copula 函數建立兩者的聯合分布，本文以GH-Copula函數為例［21，22］：

式中：FQ（q）和FW（w）分別為洪峰流量Q和洪量W的邊緣分布；θ為函數參數，變量之間的相關性越強，其值越大，θ的參數估計為：

式中：τ為Kendall 秩相關系數；（xi-xj）為觀測點據；sign（*）為符號函數，當（xi-xj）（yi-yj）＞0 時，sign（*）=1；當（xi-xj）（yi-yj）＜0 時，sign（*）=-1；當（xi-xj）（yi-yj）=0時，sign（*）=0。

對于特定重現期Tq，w，將式（8）帶入式（5）得［23］：

式中：Tu，v=Tq，w?？梢詽M足式（13）的q，w的無數種組合，在平面上繪制q，w關系圖，會出現相應重現期下的等值線，借助于式（8）～（10），式（13）可以表示為：

式中：v=L（u）為相應重現期下等值線的表達式。

根據式（14），可以得到重現期為Tq，w洪峰流量q與洪量w的等值線。該等值線上所有點的重現期均為Tq，w，當Tq，w=T時，該等值線上點（q，w）分別大于相應的qT和wT值。在設計洪水的洪峰流量頻率u與洪量頻率v相等時，有唯一解，即等值線v=L（u）與直線v=u的唯一交點，其解析表達式為：

式中：qTO和wTO為兩變量重現期Tq，w（Tq，w=Tu，v）對應的Q和W的設計值。

2.3 推求設計洪水過程線

根據所得到的設計洪峰流量qTO、設計洪量wTO以及典型洪水過程，采用變倍比放大法獲得設計洪水過程［24］：

式中：DF（t）、TF（t）分別為設計洪水和典型洪水在t時刻的流量；QTD和WTD分別為典型洪水過程的洪峰流量和歷時為TD的洪量。

2.4 設計洪水地區(qū)組成及遺傳算法的優(yōu)化

2.4.1 洪水地區(qū)組成的數學描述

設計洪水地區(qū)組成的核心內容是在匯流洪量一定的情況下，確定洪水來源的各分區(qū)洪水所占比例［25，26］。根據流域水系分布情況，本文重點研究花園口水文站、小浪底水庫及小-花區(qū)間的設計洪水地區(qū)組成。如圖2所示，A為小浪底水庫，X為小浪底入庫洪水；B為小浪底-花園口區(qū)間，Y小浪底-花園口區(qū)間洪水；Z為花園口水文站斷面，其中Z=X+Y。

圖2 研究區(qū)內洪水地區(qū)組成示意圖Fig.2 Schematic diagram of the flood region composition of study areas

本文各分區(qū)獨立選樣并采用P-III型分布作為頻率曲線，其概率密度函數為［27］：

式中：α、β、a0分別為P-III 分布的形狀、尺度以及位置參數，α＞0，β＞0；Γ(α)為α的伽馬函數。

α、β、a0與總體的統計參數Xˉ、Cv、Cs的關系為：

2.4.2 同頻率地區(qū)組成法

根據防洪要求，上游分區(qū)的洪量與下游設計斷面洪量重現期（設計頻率）相同。對于本研究，同頻率法可以有以下的兩種組合形式，分別為：當花園口斷面發(fā)生重現期為T（頻率為P）的洪水zp時，小浪底水庫也發(fā)生重現期為T的洪水xp，小-花區(qū)間發(fā)生相應洪水y，有y=zp-xp；當花園口斷面發(fā)生重現期為T（頻率為P）的洪水zp時，小-花區(qū)間也發(fā)生重現期為T的洪水yp，小浪底水庫發(fā)生相應洪水x，即x=zp-yp。

2.4.3 最可能洪水地區(qū)組成

針對傳統洪水地區(qū)組成算法缺乏對空間相關性的考慮，本文利用GH-Copula 函數構建研究區(qū)內小浪底洪水X與小-花區(qū)間洪水Y之間的聯合分布，其中，各個斷面獨立選樣并計算邊緣分布，以此推求在花園口站發(fā)生一定量級的洪水時，小浪底和小-花區(qū)間最可能的設計洪水過程。

利用GH-Copula 函數建立小浪底水庫洪水X和小-花區(qū)間洪水Y的聯合分布如下式：

式中：FX（x）和FY（y）分別為小浪底水庫洪水與小-花區(qū)間洪水的P-III型分布函數。

由式（19）可得聯合密度函數如下：

當下游花園口斷面發(fā)生設計洪水zp時，尋求上游小浪底和小-花區(qū)間洪水過程的最可能組合，即在滿足x+y=zp的條件下，求解概率密度函數f（x，y）的最大值點。將y=zp-xp代入f（x，y）中，可式（19）轉化為單變量x的函數，如下式：

上式對x求導可得：

令df/dx，整理得：

式中：αx、βx、a0.x和αy、βy、a0，y分別為小浪底與小花區(qū)間洪水要素的P-III分布FX（x）和FY（y）對應參數。

對上述的方程求解，可以得到在給定的花園口站設計洪水值zp時，小-花區(qū)間洪水與小浪底的最可能的洪水地區(qū)組合（x，y）。

2.4.4 構建Copula-GA算法

對于已經構建的基于Copula 洪水地區(qū)組成模型，求得最有可能的組合情況，即在整個變量序列里求得Copula 密度函數的極值所對應的設計洪量值。但是由于Copula 函數特性，該模型求解過程繁瑣，不能有效和準確的推求其全局最有可能的極值組合。因此，本文將遺傳算法（Genetic Algorithms）對與構建的基于Copula 洪水地區(qū)組成模型相耦合，以探尋既客觀又高效的地區(qū)組成結果［28-30］。

此處所構建的Copula-GA（Copula-Genetic Algorithms）模型，以式（23）作為遺傳算法的目標函數，即適應度函數，計算個體適應度：

Copula-GA模型的運行過程可概述如下：

（1）建立Copula洪水地區(qū)組成模型；

（2）確定目標函數f（x），即率密度函數式（23）；

（3）初始化種群；

（4）隨機生成初始種群；

（5）利用式（23）計算個體適應度，并保留適應度值較高的染色體；

（6）對種群交叉、變異得到下一代群體；

（7）重復（5）～（6），滿足終止條件后輸出適應度最大的染色體。

3 結果分析

3.1 花園口站設計洪水

本文的研究區(qū)屬于半濕潤區(qū)，暴雨強度大，一次暴雨歷時2～3 d，最長達5 d，洪水過程歷時長，可達10～12 d，因此本文構建年最大12 d 洪量W12與年最大洪峰流量Q的聯合分布。Q與W12的邊緣分布函數分別為FQ（q）和FW（w12）的P-III 型分布。使用線性矩法計算FQ（q）、FW（w12）的參數，成果如表1所示。

表1 花園口站洪水統計特征值成果Tab.1 The results of flood statistics characteristic value of Huayuankou Hydrological Station

由花園口站實測流量資料，計算得花園口站洪峰流量和12天洪量的Kendall 秩相關系數τ= 0.776，對應θ≈4.464。使用聯合變量同頻率法以及單變量同頻率法對花園口站設計洪峰流量以及12天洪量進行推求，成果見圖3、4以及表2所示。

圖3 花園口站洪峰流量和洪量聯合頻率分布圖Fig.3 Joint probability diagram of flood peak and flood volume at Huayuankou Station

由表2可知，采用聯合變量同頻率法，計算得到的各個重現期的設計洪峰流量Q和設計洪量W12大于通過單變量同頻率法計算獲得的設計洪峰流量Q和設計洪量W12。即使傳統單變量法能夠單獨控制變量的重現期，但是缺乏對變量之間相關性的考慮，然而聯合變量重現期所得的設計值不會低于設計標準。鑒于對防洪尤為不利的情況，取聯合變量重現期作為本次推求花園口站的設計洪水特征值，將會加大花園口洪水威脅程度，而由此設計結果優(yōu)化得到的調度方案會提高花園口站整體抵御洪水的能力，從而降低花園口站遭受洪災的風險，按變倍比放大獲得花園口站設計洪水過程線，結果如圖5所示。

圖5 花園口站設計洪水過程線Fig.5 Design flood process of Huayuankou Hydrological Station

表2 花園口水文站設計洪峰和洪量設計成果Tab.2 Design results of flood peak and flood volume of Huayuankou Hydrological Station

3.2 設計洪水地區(qū)組成

小浪底壩址、小-花區(qū)間采用1955-2015年連續(xù)還原天然狀態(tài)下洪水系列資料，通過矩法估計獲得小浪底水庫以及小-花區(qū)間洪水統計參數，見表3。

表3 小浪底和小-花區(qū)間洪水統計參數Tab.3 Results of flood statistics in Xiaolangdi and Xiao-hua interval

小浪底壩址和小-花區(qū)間及花園口的12天洪量邊緣分布分別采用P-III型分布，參數估計結果如表4所示。

表4 花園口、小浪底以及小-花區(qū)間P-III型分布參數估計結果Tab.4 Estimation results of P-III distribution parameters in Huayuankou Hydrological Station，Xiaolangdi and Xiao-hua interval

根據現有資料系列，分別采用同頻率地區(qū)組成法和基于GH-Copula函數和遺傳算法構建的Copula-GA 算法對小浪底與小-花區(qū)間洪水的地區(qū)組成進行計算。

（1）同頻率地區(qū)組成法通過假設小浪底洪水重現期與花園口相同推求小浪底設計洪水特征值，再利用干支流關系計算小-花區(qū)間洪水特征值，計算結果如表5所示。

表5 小浪底與花園口同頻率地區(qū)組成計算結果Tab.5 Flood region composition results of Xiaolangdi and Huayuankou with the same frequency

（2）基于Copula 函數建立小浪底和小-花區(qū)間聯合分布，并推求其Copula 密度函數，利用干支流洪量組合關系，將在特定重現期下的Copula 密度函數轉變?yōu)閱巫兞棵芏群瘮?。然后采用GA 算法尋優(yōu)求解，以推求其密度函數最大極值的可能組合，并推求小浪底與小-花區(qū)間12 天洪量的對應值，得到花園口站設計洪水最可能地區(qū)組成，成果見圖6以及表6。

圖4 花園口站洪峰和洪量聯合重現期概率圖Fig.4 Probability Diagram of Joint Return Period of Flood Peak and Flood Volume at Huayuankou Hydrological Station

圖6 小浪底和小-花區(qū)間聯合分布函數圖及小浪底洪量概率密度圖（P=2%）Fig.6 Xiaolangdi and Xiao-Hua interval joint distribution function diagram and Xiaolangdi flood probability density diagram（P=2%）

通過表6 可以看出，隨著花園口洪量的增大，小浪底與小-花區(qū)間洪量量級都相應的隨之增大，并且小浪底洪量所占的比例明顯加大。由此可以得出：小浪底的洪水對花園口洪水的發(fā)生起到了決定性作用，如果僅發(fā)生小-花區(qū)間為主的洪水，花園口一般不會出現災害性的洪水；只有發(fā)生全流域型洪水時，即小-花區(qū)間與小浪底同時發(fā)生洪水時花園口才會出現災害性大洪水，這與黃河中下游洪水成因及來源相印證［17，31］。

表6 花園口站設計洪水最可能地區(qū)組成成果Tab.6 The most likely flood region composition results of the design flood at Huayuankou Station

由花園口站設計洪水最可能地區(qū)組成成果可知，將小浪底在各頻率地區(qū)組成所分配的洪量值，由小浪底洪量頻率曲線推求出其相對應的頻率值；然后再根據洪峰、洪量同頻率的這一假設，分別推求出與小浪底各洪量頻率對應相同的洪峰值，成果見表7。

表7 小浪底壩址設計洪水洪峰、洪量成果Tab.7 Design flood peak and flood volume results of Xiaolangdi dam site

對同頻率地區(qū)組成法和Copula-GA 法計算的小浪底的洪峰流量和洪量結果分析可知：當花園口發(fā)生的洪水重現期為200年（P=1%）、100年（P=1%）以及50年（P=2%）時，小浪底最可能發(fā)生的洪水重現期約為52 a（P=1.92%）、29 a（P=3.44%）和18 a（P=5.76%）；通過假設小浪底洪水重現期與花園口相同計算得到的小浪底設計洪水，相較Copula-GA 法計算得到的最可能的設計洪水地區(qū)組成，前者所推求的小浪底洪水特征值較大，進而導致小-花區(qū)間洪水洪水特征值較小，不能達到相應的設計頻率。

3.3 小浪底及小-花區(qū)間設計洪水過程

根據求得的小浪底對應花園口各設計頻率的洪峰、洪量設計值，放大所選取的典型洪水，推求小浪底52年一遇、29年一遇、18年一遇設計洪水過程，成果見圖7。

圖7 小浪底設計洪水過程線Fig.7 Xiaolangdi design flood process

根據小浪底設計洪水過程，由馬斯京根方法將小浪底洪水演進至花園口斷面。通過花園口各重現期設計洪水過程線逐時刻減去小浪底演進至花園口過程線，即可得到小-花區(qū)間設計洪水過程線。

4 結論

本文使用Copula 函數構建了花園口站洪峰流量與洪量的聯合分布，并計算得到了特定重現期的洪峰流量以及洪量，充分地考慮了洪峰流量與洪量之間的相關性，通過對比分析現有設計成果發(fā)現，兼顧洪峰和洪量相關性計算得到的各重現期下的洪水峰值流量和時段洪量設計值略高于現有的設計成果值［13］，更全面的考慮洪水要素之間的相互作用，也將更適用于受洪峰和洪量共同控制的區(qū)域。

在計算花園口站各重現期下的設計洪水地區(qū)組成時，針對傳統的設計洪水地區(qū)組成計算方法中存在問題，本文采用Copula 函數構建了小浪底與小-花區(qū)間12 天洪量的聯合分布，并與遺傳算法相耦合構建了Copula-GA 模型，以快速準確的推求Copula 密度函數的最大值組合，即設計洪水最可能地區(qū)組成，得到的結果更具備客觀性；當下游防洪斷面將流經特定重現期洪水時，通過本研究最可能的地區(qū)組成方法所推得上游各斷面洪水重現期，相比傳統同頻率方法，也更符合研究區(qū)洪水的實際情況。

洪水歷時作為另一個洪水特征量，并與洪峰流量、洪量存在一定的相關性［21］。本文僅考慮了洪峰流量、洪量的相關性，建立了洪峰流量、洪量的聯合分布，在后續(xù)研究中還應建立洪峰流量、洪量以及洪水歷時間的聯合分布，推求不同重現期下的洪峰流量、洪量、洪水歷時組合?！?/p>