不同重現期標準雙變量設計洪水計算方法

晉 恬聞 昕方國華黃顯峰袁 玉

(1.河海大學水利水電學院,江蘇南京 210098;2.中國水利水電科學研究院流域水循環模擬與調控國家重點實驗室,北京 100038)

設計洪水是我國水利工程在規劃、設計、運行、管理等階段的重要防洪依據,常通過同倍比、同頻率放大典型洪水的方法得到。同倍比法僅考慮單個洪水特征變量因素,同頻率法雖然考慮了洪峰流量和洪量兩個因素,但未考慮二者之間的聯系,對于感潮河流等洪水成因復雜的流域存在較大的防洪隱患和風險。雙變量設計洪水則考慮了洪水成因的復雜性以及洪水事件本身的多屬性特征,可為工程和區域防洪安全提供新的視角和科學依據。目前,對于雙變量設計洪水的計算方法已有較多探索和研究,但這些研究多針對單一重現期標準或單一計算方法,不同重現期標準和不同計算方法,特別是不同雙變量重現期標準適用的設計洪水值計算方法尚不明晰。

Copula函數[1]是研究變量相關性測度的工具,可以刻畫變量之間的非線性關系,具有邊緣分布靈活以及其形式不受邊緣分布限制的特點,可用于洪水、暴雨、干旱等多變量水文事件的聯合概率分布分析[2鄄4]。 Archimedean Copula 函數由于具有構造簡單、僅含有一個參數的優點,水文事件聯合概率分布多選用該類函數進行擬合分析。目前常用于描述雙變量洪水頻率的重現期標準有聯合重現期和同現重現期兩種,Salvadori等[5鄄7]根據描述危險區域的不同提出了第二重現期,并將其與聯合重現期進行了比較分析,確定兩者有較強的指數相關關系;史黎翔等[8]進一步研究了上述3種重現期標準的差異,論證了第二重現期在危險區域描述方面的合理性。目前確定設計洪水值的方法[9鄄19]大致分為兩類:第一類為基于重現期等值線確定設計洪水值的方法,常用的有最可能超越法、極大似然法、雙變量同頻率法;另一類為先確定滿足條件的洪水峰量組合線,再求出該曲線與洪水重現期等值線交點(即為對應重現期的設計洪水值)的方法,常用的有條件最可能組合法、條件期望組合法、雙變量同頻率組合法。

本文基于Archimedean Copula函數建立峰量組合二維分布模型,選擇C測度、聯合、同現重現期,采用極大似然法、條件最可能組合法和雙變量同頻率法,研究不同重現期標準下雙變量設計洪水的計算方法。

1 基于Copula函數的洪水峰量聯合分布模型

洪水具有多特征變量的屬性,如洪峰流量、洪量、歷時、洪水過程線等,本文選擇洪峰流量、洪量兩個變量構建洪水聯合分布模型。

定義q、w分別為一場洪水過程的洪峰流量和洪量,設FQ(q)、FW(w)為洪峰流量、洪量的邊緣分布函數,其聯合分布函數為F(q,w),按照Sklar定理[1],洪峰流量、洪量的聯合分布函數可以用對稱Archimedean Copula函數來表示:

其中

式中:茲為 Copula函數的參數;u、v分別為洪峰流量、洪量的邊緣分布函數。

2 雙變量重現期標準計算方法

目前描述雙變量的常見重現期標準為聯合重現期和同現重現期。不同重現期標準的差異本質為定義危險事件方式的不同。

若以洪峰流量和洪量兩個特征量來表述洪水事件,對于聯合重現期T胰,其危險事件E胰表示洪峰流量或洪量其中一個變量超過設定的閾值(q0或w0),對應的危險域AE胰={(q,w)沂R:q>q0胰w>w0}。聯合重現期可用Copula函數表示為

對于同現重現期T疑,其危險事件E疑表示洪峰流量和洪量同時超過設定的閾值,對應的危險域為同現重現期可用Copula函數表示為 假設一個洪水事件的雙變量重現期為T,聯合重現期和同現重現期的危險區域如圖1所示,圖中的實線即為重現期等值線。圖1(a)中玉區域為在聯合重現期等值線上任選一點A(u1,v1)所對應的危險區域,圖1(b)中域區域為在同現重現期等值線上任選一點B(u2,v2)所對應的危險區域,通過公式u=FQ(q)、v=FW(w)可以推求出等值線上不同點所對應的洪峰流量q以及洪量w。結合圖1(a)和(b)進行分析,在重現期等值線上選擇不同的點會得到相對應的危險區域,選定的點不同,對應的危險區域也不同,對于上述兩種重現期標準,其重現期水平與危險區域并非一一對應關系。另外,聯合、同現重現期對危險區域的定義也存在不合理性,一般而言,重現期越大相應的危險區域越小,重現期大的事件的危險區域應包含于重現期小的事件中,但在聯合、同現重現期標準下存在不滿足上述規律的情況,結合圖2可看出,B點處的重現期大于A點處的重現期(TB>TA),B點處的危險區域包括事件E,即事件E在B點處被識別為危險事件,但是在A點,其危險區域并不包括危險事件E。采用上述兩種重現期標準,均出現了重現期較大的洪水組合識別的危險事件在重現期較小的洪水組合中被識別為非危險事件的矛盾的情景。

圖1 重現期的危險區域

圖2 重現期危險區域識別分析

基于上述兩種重現期標準的局限性,將危險事件定義為洪水峰量組合在重現期等值包絡線以上的區域,即圖1(c)中的芋區域。假設重現期等值線滿足C茲(u,v)=t,其對應的危險區域為AEC(t)={(u,v)沂玉2:C茲(u,v)>t},其中0t)可以用C測度函數KC(t)[13]來表示:

式中(t)為Copula函數的生成元。

在C測度函數下對應的重現期為C測度重現期TC,表達式為 C測度重現期等值線上任意一點對應的危險區域均相同,且重現期水平越大,危險區域越小,重現期水平小的危險區域必定會覆蓋重現期水平大的危險區域。可見C測度重現期對危險區域的定義更為合理。

C測度重現期是以聯合重現期等值線為基礎通過修正危險區域得到的新的重現期標準,因此其等值線線型與聯合重現期類似。C測度函數是以Archimedean Copula函數為前提條件推導得到的,在此基礎上定義的C測度重現期僅適用于以Archimedean Copula函數構建的聯合分布模型,對于其他類型的Copula函數的適用性還需論證。

3 設計洪水值計算

對于選定的洪水雙變量重現期水平,存在一系列峰量組合滿足防洪標準要求,需根據一定的準則從中選擇設計洪水值組合,本文選擇極大似然法、條件最可能組合法以及雙變量同頻率法作為推求設計洪水值的方法,并以Archimedean Copula函數中的G鄄H Copula函數為例進行分析。

3.1 極大似然法

在相應重現期等值線上,要求設計洪水值組合(qd,wd)滿足權重值GML最大,其目標函數為

式中:c茲(u,v)為 Copula函數的聯合密度函數;fQ(q)、fW(w)分別為洪峰流量、洪量對應的邊緣密度函數。

式(6)應滿足重現期等值線約束,即保證所求設計洪水值組合(qd,wd)均在對應重現期等值線上,約束條件為

3.2 條件最可能組合法

當洪峰流量q取qi時,要求洪量w對應的密度函數值fW|Q(w)最大,其目標函數為

求出洪峰流量qi對應的洪量wi,其對應的洪水峰量組合為(qi,wi),得出一系列的峰量組合點連成條件最可能洪水峰量組合線,其與對應重現期等值線的交點即為設計洪水值組合(qd,wd)。

式(8)以洪峰流量為條件進行分析計算。

3.3 雙變量同頻率法

采用雙變量同頻率法求解設計洪水值組合,即要在相應重現期T上求u=v的點為設計洪水值組合(qd,wd),應滿足以下條件:

也可先求出雙變量同頻率洪水降量組合線,該線洪峰流量與洪量的邊緣分布函數相同,其目標函數為

雙變量同頻率洪水峰量組合線與對應重現期等值線的交點即為設計洪水值組合(qd,wd)。

4 算例分析

新安江水庫建在錢塘江支流新安江銅官峽,位于建德市新安江鎮以北約4郾5 km,壩址以上河長為323 km,壩址控制流域面積10 442 km2,主要承擔發電、防洪、供水等任務[20鄄21]。羅桐埠站是新安江水庫的設計依據站,本文選取羅桐埠站1961—2008年的洪水資料,根據新安江流域的洪水特征,提取歷時為48 a的年最大洪峰流量q和年最大10 d洪量w10的峰量數據,進行峰量聯合設計洪水研究。

4.1 邊緣分布

首先分別建立年最大洪峰流量q和年最大10 d洪量w10的邊緣分布函數。GEV分布被廣泛用于洪峰流量、洪量的概率統計分析擬合中并取得了較好的擬合效果,本文選擇GEV概率分布函數作為邊緣分布函數,采用極大似然法估計 FQ(q)和FW(w)的統計參數,運用K鄄S擬合優度檢驗方法對邊緣分布進行檢驗,計算q和w10的樣本檢驗統計值和顯著水平琢=0郾05的檢驗臨界值,結果見表1。對于q和w10,樣本檢驗統計值均小于相應的檢驗臨界值,表明FQ(q)和FW(w)均服從GEV分布。

表1 洪水特征量邊緣分布統計參數

4.2 聯合分布

選擇Archimedean Copula函數中的G鄄H Copula函數、Clayton Copula函數和 Frank Copula函數對q鄄w10聯合分布進行模擬,采用極大似然估計法計算參數茲,并采用赤池信息量準則(Akaike information criterion,AIC)、貝葉斯信息量準則(Bayesian information criterion,BIC)和均方根誤差(RMSE)3種擬合優度檢驗方法選擇最優Copula函數,AIC、BIC、RMSE值越小代表Copula函數的擬合度越好。計算結果見表2,其中G鄄H Copula的AIC、BIC值均最小,RMSE值略大于Frank Copula函數,基于擬合結果選擇G鄄H Copula構建q鄄w10聯合分布。圖3為經驗聯合分布頻率與理論聯合分布(G鄄H Copula函數)頻率擬合結果,樣本數據均在45毅線附近,表明G鄄H Copula函數擬合度較好。

表2 q鄄w10聯合分布Copula函數相關參數

圖3 經驗聯合分布頻率與理論聯合分布頻率擬合結果

4.3 不同重現期標準比較分析

根據同現重現期、聯合重現期、C測度重現期計算公式,繪制不同重現期標準下新安江流域10年、20年、50年、100年、200年一遇洪水的重現期等值線如圖4所示。由圖4分析可知,對相同重現期標準,不同重現期水平的等值線走勢一致;對不同重現期標準,聯合重現期等值線與C測度重現期等值線走勢相同,聯合、C測度重現期與同現重現期等值線走勢相反。從危險事件定義的角度分析,同現重現期與聯合重現期等值線的走勢差異性體現在前者考慮的為并集危險事件,后者考慮的為交集危險事件;而C測度重現期是基于聯合重現期等值線提出的修正,因此二者走勢相似。

圖4 不同重現期標準下重現期等值線及設計洪水值組合

圖5 不同重現期水平下重現期等值線及設計洪水值組合

在給定重現期水平下,不同重現期標準之間的關系見圖5。同現、聯合和C測度重現期在同一重現期水平下其重現期等值線無交點;對于聯合重現期與同現重現期,其重現期等值線的漸近線均為單變量(q或w10)在相同重現期水平下的設計洪水值;對于C測度重現期,其重現期漸近線均小于單變量(q或w10)在相同重現期水平下的設計洪水值。對于同一洪峰流量值(洪量值),聯合重現期等值線對應的洪量值(洪峰流量值)最大,C測度重現期等值線對應的次之,同現重現期等值線對應的最小。這是由于同現重現期標準下的峰量組合值僅考慮了并集危險事件發生的概率,使得其峰量組合值均偏小;而聯合重現期標準下的峰量組合值卻相對偏大,因為其考慮的是交集危險事件發生的概率,即放大了實際的危險區域;C測度重現期危險區域的范圍在二者之間,因此其對應的峰量組合值也在二者范圍之內。峰量組合值的分布區域與危險區域的劃分方式呈對應關系,選擇不同的重現期計算標準即確定了不同的峰量組合值的分布區域,是確定設計洪水值的前提條件。

4.4 設計洪水值計算

以 10 a、20 a、50 a、100 a、200 a 為重現期,采用極大似然法、條件最可能組合法和雙變量同頻率法,分別計算同現重現期、聯合重現期、C測度重現期標準下的設計洪水值組合,結果見圖4。

4.4.1 不同重現期標準設計洪水

對不同重現期標準,比較相同方法推求設計洪水值組合的差異性。圖5中A、B、C點分別為在不同重現期標準下利用極大似然法所推求的設計洪水值,以百年一遇洪水為例,同現、聯合、C測度重現期的設計洪水峰量值組合分別為(15 606郾06 m3/s,511郾76 億 m3)、 (18741郾60 m3/s,580郾81 億 m3)和(16436郾38 m3/s,524郾24 億 m3),聯合重現期的設計洪水值組合中的峰量值最大,C測度重現期次之,同現重現期的最小。由圖5可看出,對于任意給定的重現期水平,3種重現期標準的設計洪水值組合間的大小關系均滿足上述規律。

重現期的選擇對設計洪水值組合起著決定性的作用。在同現重現期標準下得到的設計洪水值組合最小,危險范圍考慮得最小,相應的設計標準、工程造價最低,安全性最低。在聯合重現期標準下得到的設計洪水值組合最大,危險范圍考慮最大,相應的設計標準和工程造價最高,因此安全性也最高。水利工程設計過程中采用同現重現期標準,即以洪峰流量、洪量值同時大于某一峰量設計值作為危險事件,忽略了短時強降雨以及連續降雨引發的洪水的危險因素,低估了洪水的危險范圍,確定的防洪設計標準偏低。若以洪峰流量或洪量值中的其中一個變量大于某一峰量設計值作為危險事件,在此基礎上確定設計標準,以其對應的重現期為聯合重現期,在危險事件的確定上,則考慮了洪峰流量較大而洪量較小、洪峰流量較小而洪量較大這類不構成相應洪水威脅的條件,過分強調安全性,提高了設計標準。上述兩種水利工程設計標準都過于極端,僅適用于有特殊設計要求的工程,未能體現出雙變量聯合分布在水利工程設計中的優勢。從定義角度分析,C測度重現期是同現重現期和聯合重現期危險事件的折衷考慮,既考慮了洪水特征變量之間的聯系,又達到了經濟性與安全性的平衡。在雙變量頻率分析中,推薦采用C測度重現期作為雙變量聯合分布的重現期計算標準。

4.4.2 不同計算方法設計洪水

以雙變量同頻率法為基準,比較采用極大似然法、條件最可能組合法推求設計洪水值組合差異,結果見表3~5以及圖4。從整體上看,各組合的差異均維持在2%左右,最低為0郾37%,最高為2郾51%,3種設計洪水計算方法具有較強的穩定性,差異性較小。在同現重現期標準下采用3種不同的方法,由極大似然法推求的設計洪水值組合位于雙變量同頻率法和條件最可能組合法之間,且隨著重現期水平的提高,極大似然法推求的設計洪水值組合越接近雙變量同頻率法推求的設計洪水值組合。在聯合重現期以及C測度重現期標準下采用3種不同的方法,由雙變量同頻率法推求的設計洪水值組合位于極大似然法和條件最可能組合法之間,且極大似然法與雙變量同頻率法的接近程度更高。可通過選擇中間值的方法來確定設計洪水值的計算方法,即在同現重現期標準下,可采用極大似然法,在聯合重現期及C測度重現期標準下,可采用雙變量同頻率法。

表3 同現重現期標準下設計洪水值組合差異

表4 聯合重現期標準下設計洪水值組合差異

表5 C測度重現期標準下設計洪水值組合差異

圖6為采用雙變量同頻率法與條件最可能組合法的峰量組合曲線,條件最可能組合法的設計洪水值組合均大于雙變量同頻率法的對應值,且差距越來越大,但在重現期等值線上的差異維持在2%左右。

圖6 峰量組合曲線

5 結 論

a.采用C測度重現期作為雙變量頻率分析的標準,可確定合理的危險事件的范圍,有效規避以同現重現期為標準帶來的安全隱患,解決以聯合重現期為標準造成的性價比低的問題。

b.在不同的重現期標準下確定不同重現期水平下的設計洪水值組合時,采用的3種方法計算結果較為接近,可通過選擇中間值的方法來確定設計洪水值的計算方法。在同現重現期標準下,可選擇極大似然法作為雙變量設計洪水的計算方法,在聯合重現期及C測度重現期標準下,可選擇雙變量同頻率法作為雙變量設計洪水的計算方法。

[1]SKLAR M.Functions de repartition a n dimensions et leurs marges[J].Publ Inst Statist Univ Paris,1959(8):229鄄231.

[2]王鵬新,馮明悅,孫輝濤,等.基于主成分分析和Copula函數的干旱影響評估研究[J].農業機械學報,2016,47(9):334鄄340.(WANG Pengxin,FENG Mingyue,SUN Huitao,etal.Droughtimpactassessmentbased on principal component analysis and copula function[J].Transactions ofthe Chinese Society forAgricultural Machinery,2016,47(9):334鄄340.(in Chinese))

[3]關帥,林穎妍,查悉妮,等.基于Copula函數的韓江流域干支流洪水遭遇分析[J].中山大學學報(自然科學版),2015,54(5):130鄄137.(GUAN Shuai,LIN Yingyan,ZHA Xini,et al.Copula function鄄based flood coincidence probability analysis for mainstream and tributary of the Hanjiang River Basin[J].Acta Scientiarum Naturalium Universitatis Sunyatseni,2015,54(5):130鄄137.(in Chinese))

[4]高玉琴,陳釔西,趙立梅,等.秦淮河流域不同頻率降雨聯合概率分析[J].水電能源科學,2016(3):1鄄5.(GAO Yuqin,CHEN Yixi,ZHAO Limei,et al.Joint probability analysis of different frequencies rainfall of Qinhuaihe River Basin[J].Water Resources and Power,2016(3):1鄄5.(in Chinese))

[5]SALVADORI G,MICHELEC D,DURANTE F.On the return period and design in a multivariate framework[J].Hydrology and Earth System Sciences,2011,15(11):3293鄄3305.

[6]SALVADORI G,DE MICHELE C.Multivariate multiparameter extreme value models and return periods:a copula approach[J].Water Resources Research,2010,46(10):W10501.

[7]SALVADORI G.Bivariate return periods via 2鄄Copulas[J].Statistical Methodology,2004,1(1):129鄄144.

[8]史黎翔,宋松柏.基于Copula函數的兩變量洪水重現期與設計值計算研究[J].水力發電學報,2015,34(10):27鄄34.(SHI Lixiang,SONG Songbai.Calculation of return periods and design values of bivariate floods based on Copula function[J].Journal of Hydroelectric Engineering,2015,34(10):27鄄34.(in Chinese))

[9]劉和昌,梁忠民,姚軼,等.基于Copula函數的水文變量條件組合分析[J].水力發電,2014,40(5):13鄄16.(LIU Hechang,LIANG Zhongmin,YAO Yi,et al.Analysis on conditional compositions of hydrological variants based on Copula function[J].Water Power,2014,40(5):13鄄16.(in Chinese))

[10]李天元,郭生練,劉章君,等.基于峰量聯合分布推求設計洪水[J].水利學報,2014,45(3):269鄄276.(LI Tianyuan,GUO Shenglian,LIU Zhangjun,et al.Design flood estimation based on bivariate joint distribution of flood peak and volume[J].JournalofHydraulic Engineering,2014,45(3):269鄄276.(in Chinese))

[11]郭生練,劉章君,熊立華.設計洪水計算方法研究進展與評價[J].水利學報,2016,47(3):302鄄314.(GUO Shenglian,LIU Zhangjun,XIONG Lihua.The research progress and the estimate about design flood calculation method[J].Journal of Hydraulic Engineering,2016,47(3):302鄄314.(in Chinese))

[12]芮孝芳,張超.論設計洪水計算[J].水利水電科技進展,2014,34(1):20鄄26.(RUI Xiaofang,ZHANG Chao.Research on design flood calculation[J].Advances in Science and Technology of Water Resources,2014,34(1):20鄄26.(in Chinese))

[13]陳甜,董增川,賈本有,等.考慮上游梯級水庫群調度的小南海水電站壩址設計洪水分析[J].河海大學學報(自然科學版),2014,42(6):476鄄480.(CHEN Tian,DONG Zengchuan,JIA Benyou,et al.Analysis of design flood of dam site of Xiaonanhai Hydropower Station considering flood control operation of upstream cascade reservoir group[J].Journal of Hohai University(Natural Sciences),2014,42(6):476鄄480.(in Chinese))

[14]黃振平,林小麗,侯云青,等.P鄄芋型分布參數的矩估計與設計洪水的計算頻率[J].河海大學學報(自然科學版),2005,33(1):49鄄51.(HUANG Zhenping,LIN Xiaoli,HOU Yunqing,etal.Momentestimation for parameters of Pearson Type鄄芋 distribution and calculated frequency of design flood[J].Journal of Hohai University(Natural Sciences),2005,33(1):49鄄51.(in Chinese))

[15]李杰友,陸波,穆方馳.棉花灘水庫建成后對青溪水庫設計洪水影響研究[J].水利水電科技進展,2005,25(1):25鄄27.(LI Youjie,LU Bo,MU Fangchi.Effect of Mianhuatan Reservoir on the design flood of Qingxi Reservoir[J].Advances in Science and Technology of Water Resources,2005,25(1):25鄄27.(in Chinese))

[16]徐長江,郭生練,熊立華,等.英國設計洪水研究進展[J].水利水電科技進展,2004,24(2):66鄄69.(XU Changjiang,GUO Shenglian,XIONG Lihua,et al.Research advances in England design flood[J].Advances in Science and Technology of Water Resources,2004,24(2):66鄄69.(in Chinese))

[17]黃振平,王春霞,馬軍建.歷史洪水重現期的誤差對設計洪水的影響[J].河海大學學報(自然科學版),2002,30(1):79鄄82.(HUANG Zhenping,WANG Chunxia,MA Junjian.Effectofestimation errorof historical flood recurrence interval on design flood[J].Journal of Hohai University(Natural Sciences),2002,30(1):79鄄82.(in Chinese))

[18]華家鵬.小流域設計洪水的計算方法[J].河海大學學報 (自 然 科 學 版),1998,26(2):100鄄104.(HUA Jiapeng.Calculation method of design flood for small watershed[J].Journal of Hohai University(Natural Sciences),1998,26(2):100鄄104.(in Chinese))

[19]NELSENR B.An introduction to Copulas[M].New York:Springer,1999.

[20]吳書悅,趙建世,雷曉輝,等.氣候變化對新安江水庫調度影響與適應性對策[J].水力發電學報,2017,36(1):50鄄58.(WU Shuyue,ZHAO Jianshi,LEI Xiaohui,et al.Impacts of climate change on operation of Xin爺an River Rreservoirand adaption strategies[J].Journalof Hydroelectric Engineering,2017,36(1):50鄄58.(in Chinese))

[21]鄭艷妮,聞昕,方國華.新安江流域氣候變化及徑流響應研究[J].水資源與水工程學報,2015,26(1):106鄄110.(ZHENG Yanni,WEN Xin,FANG Guohua.Research on climate change and runoff response in Xin爺an River Basin[J].JournalofWaterResources & Water Engineering,2015,26(1):106鄄110.(in Chinese))