(長江水利委員會(huì) 長江上游水文水資源勘測局,重慶 400020)
試驗(yàn)與研究
連續(xù)API產(chǎn)匯流模型的改進(jìn)及應(yīng)用
曹磊陳川建
(長江水利委員會(huì) 長江上游水文水資源勘測局,重慶 400020)
結(jié)合長江流域上游水文特性和實(shí)際預(yù)報(bào)工作,采用連續(xù)的日消退系數(shù)計(jì)算方法、連續(xù)的地下水分割方式、自動(dòng)計(jì)算降雨徑流相關(guān)曲線(P~Pa~R曲線)和瞬時(shí)單位線計(jì)算方法,并結(jié)合降雨徑流相關(guān)計(jì)算流程進(jìn)行產(chǎn)匯流計(jì)算和評(píng)定,得到一套較完整的連續(xù)API降雨徑流率定參數(shù)。通過該方法將API模型歸納成11個(gè)數(shù)字參數(shù),方便利用智能優(yōu)化算法進(jìn)行優(yōu)化率定計(jì)算,提高了API洪水預(yù)報(bào)模型的率定精度。提出的連續(xù)API模型的計(jì)算方法,為降雨徑流相關(guān)模型方案的建立、率定和實(shí)際應(yīng)用提供了一套新思路。
產(chǎn)匯流計(jì)算;連續(xù)API模型;蓄滿產(chǎn)流;模型參數(shù)率定
API模型區(qū)間產(chǎn)匯流模型作為一種洪水預(yù)報(bào)模型在大多數(shù)流域得到廣泛使用。該模型的參數(shù)率定涉及到場次洪水的選峰、地下水分割、凈流量計(jì)算、時(shí)段單位線率定等多方面,率定過程較復(fù)雜,主要靠人工結(jié)合計(jì)算機(jī)交互技術(shù),采用單場次洪水進(jìn)行交互分析和數(shù)據(jù)預(yù)處理才能完成,特別是利用人工退水曲線的方式來分割場次凈雨具有一定的隨意性。隨著地下水的分割[1]和還原[2-4]方法的使用,以及利用瞬時(shí)單位線取代時(shí)段單位線,逐漸與計(jì)算機(jī)技術(shù)結(jié)合后,提出了連續(xù)API模型并得到改進(jìn)[5-6],已在區(qū)間流計(jì)算方面得到推廣和應(yīng)用[7],取得了一定的效果。然而,目前還是通過選擇場次洪水人工交互分析得到連續(xù)API模型中退水曲線及場次洪水的分割,以及核心參數(shù)(P~Pa~R曲線)。同時(shí),日消退系數(shù)K值,仍然靠人為來確定,這些都約束了連續(xù)API模型的率定和自動(dòng)計(jì)算的應(yīng)用,沒有真正達(dá)到連續(xù)率定和連續(xù)計(jì)算的目的。本文以蓄滿產(chǎn)流原理為基礎(chǔ),結(jié)合地下水分割和分布水文模型中的植被指數(shù)等,結(jié)合計(jì)算機(jī)技術(shù)提出了一套真正意義上連續(xù)率定、連續(xù)計(jì)算的連續(xù)API模型,包括自動(dòng)進(jìn)行地下水還原計(jì)算、P~Pa~R曲線率定、計(jì)算瞬時(shí)單位線和每月K值的計(jì)算過程。
API模型是前期降雨指數(shù)(Antecedent Precipitation Index)模型的縮寫,其使用的表達(dá)式是傳統(tǒng)的降雨徑流相關(guān)圖,故又稱為降雨徑流經(jīng)驗(yàn)相關(guān)法。在20世紀(jì)40年代,R.K.林斯雷等提出了以季節(jié)、歷時(shí)、前期降雨指數(shù)API為參數(shù)的五變數(shù)降雨徑流相關(guān)圖,配合時(shí)段單位線技術(shù),成為一種有效的降雨洪水預(yù)報(bào)模型。1969年,西納(Sittner)等提出了模型地下徑流的建議,結(jié)合謝爾曼時(shí)段單位線即構(gòu)成了可模擬流域降雨凈雨過程的API模型。
API模型屬于多輸入、單輸出靜態(tài)的系統(tǒng)數(shù)學(xué)模型。模型建立之初,主要用于一次洪水徑流量預(yù)報(bào)。在我國,模型經(jīng)簡化,形成了四變數(shù)關(guān)系線、三變數(shù)關(guān)系圖(P~Pa~R)或兩變數(shù)關(guān)系圖(P+Pa~R)。簡化后的API模型便于方案的率定和使用。其降雨產(chǎn)流函數(shù)關(guān)系可表示為:R=f(P,Pa,T,季節(jié)) 、R=f(P,Pa,T) 、R=f(P,Pa)或R=f(P,P+Pa)。其中,北方流域使用較多P+Pa~R產(chǎn)流方案,而南方、特別是長江流域多采用P~Pa~R三軸曲線。API模型的主要參數(shù)有:最大初損(IM)、前期降雨指數(shù)(Pa)、日消退系數(shù)(K)、降雨凈流相關(guān)曲線(P~Pa~R)和流域時(shí)段單位線(UH),計(jì)算時(shí)采用流域平均降雨來推算產(chǎn)匯流,參數(shù)的求解和含義可參考其他相應(yīng)文獻(xiàn)。
在模型使用過程中,水文預(yù)報(bào)工作者對(duì)該模型進(jìn)行了不斷改進(jìn)[6],包括Pa計(jì)算方式的改進(jìn),并將分塊計(jì)算方式引入到了API模型中[6],考慮到降雨分布不均勻情況下的產(chǎn)流計(jì)算;將API的計(jì)算方式與計(jì)算機(jī)結(jié)合提出API模型的連續(xù)計(jì)算方法[7];將納西(Nash)瞬時(shí)單位線(IUH)引入模型計(jì)算中,利用調(diào)蓄滯時(shí)(Ku)和串聯(lián)水庫個(gè)數(shù)(N)兩個(gè)參數(shù)來確定流域瞬時(shí)單位線,再根據(jù)流域面積和計(jì)算時(shí)段長也可將瞬時(shí)單位線轉(zhuǎn)化成時(shí)段單位線,使得單位線的率定更方便。
不管是傳統(tǒng)的API模型還是改進(jìn)后的模型,在進(jìn)行模型率定時(shí),主要是靠人工結(jié)合計(jì)算機(jī)交互功能,經(jīng)過場次樣本選取、基流分割、徑流計(jì)算、單位線推求、樣本評(píng)定5個(gè)主要的復(fù)雜步驟來完成。特別是地下水分割時(shí)具有一定的隨意性和不確定性(如圖1),使得模型的率定和實(shí)際使用過程中出現(xiàn)不一致的情況,這一定程度上制約了API模型的率定與使用,影響了模型計(jì)算精度。
圖1 地下水不同分割示意
針對(duì)API模型的不足做了改進(jìn),提出如下連續(xù)API模型的率定和計(jì)算方法。
2.1.1 常規(guī)分割方法的局限性
在日常的API模型參數(shù)率定過程中,通常采用圖1中列出的方法來進(jìn)行基流分割,分割方式方法有一定的隨意性,計(jì)算的凈雨量和實(shí)際凈雨量有時(shí)差異較大。在API模型方案參數(shù)的實(shí)際使用過程中,降雨未結(jié)束洪峰未出現(xiàn)時(shí)(預(yù)報(bào)根據(jù)時(shí)間后的過程),無法確定直接徑流停止點(diǎn)(B點(diǎn))(平割法除外),同時(shí)也很難確定計(jì)算的直接徑流和地下水過程的總流量過程。不少使用者在直接徑流的基礎(chǔ)上加一個(gè)基流定值作為流域出口總出流,相當(dāng)于不管采用哪種分割方式確定的方案都當(dāng)成平割法進(jìn)行處理,這種處理在一定程度上降低了預(yù)報(bào)的精度,計(jì)算流程也不是太合理。
2.1.2 數(shù)字濾波
濾波起源于通信,它是將信號(hào)中特定波段頻率濾除得到有用信號(hào)的一種技術(shù)。內(nèi)森(Nathan)和麥克馬洪(McMahon)在1990年首次將數(shù)字濾波技術(shù)應(yīng)用到基流分割中,該方法已成為近年應(yīng)用最廣泛的基流分割方法[8]。分割方程如下:
(1)
Qb2=Q2-Qd2
(2)
式中,Qd2和Qd1分別為濾波得到的第2時(shí)刻和第1時(shí)刻的地表徑流,m3/s;Q2和Q1分別為第2時(shí)刻和第1時(shí)刻的徑流,m3/s;Qb2為第2時(shí)刻的基流,m3/s;α為濾波系數(shù)。
這是一個(gè)遞推公式,基于較長的系列資料進(jìn)行計(jì)算時(shí),特別是枯季開始時(shí),Qd1比較小,可設(shè)置一個(gè)初值,利用公式(1)計(jì)算得到Qd2,再利用公式(2)計(jì)算得到第2時(shí)刻的基流,以此類推計(jì)算得到整個(gè)洪水過程的基流。這種方法在水文中得到不斷的應(yīng)用和推廣,取得了較好的效果[9-11]。
2.1.3 地下水還原計(jì)算[3]
1969年美國天氣局提出API模型中地下水出流采用馬斯京根法線性水庫(x=0)演算,模型中假定入滲到地下水庫中的下滲量f與直接徑流Q呈線性關(guān)系,即:
f=Z(Q-Qg)
(3)
地下水匯流用馬法演算,可以得到如下公式:
Qg2=C0(f1+f2)+C2Qg1
(4)
將式(3)帶入式(4)得到如下地下水推算公式:
(5)
式中,Qg1為上時(shí)段初地下水流量,m3/s;Q1、Q2為時(shí)段始末的直接徑流量,m3/s;C0、C2為馬法演算參數(shù);Qg2為要計(jì)算的時(shí)段末地下水流量,m3/s;Z為地下水入流比例系數(shù)。
該方法基于直接徑流進(jìn)行地下水還原計(jì)算,避開了地下水的分割難題,在實(shí)際預(yù)報(bào)過程中采用公式(5)即可推導(dǎo)得出。
根據(jù)新安江三水源模型的蓄水容量曲線和降雨徑流關(guān)系的詮釋,蓄水容量曲線與降雨徑流具有對(duì)應(yīng)關(guān)系,故流域蓄水容量曲線的線形經(jīng)確定后降雨徑流關(guān)系曲線便隨之確定[3]。蓄水容量曲線的線性在濕潤地區(qū)廣泛采用拋物線,因此,給定流域最大平均蓄水容量Wm(API模型里為IM)和拋物線的指數(shù)b,根據(jù)不同的前期蓄水量(在蓄滿產(chǎn)流中為W0,在API模型中為Pa),可對(duì)應(yīng)生成一簇P~Pa~R曲線或一條P+Pa~R曲線,把復(fù)雜的關(guān)系曲線族歸結(jié)于IM、b兩個(gè)參數(shù)。
根據(jù)蓄滿產(chǎn)流的原理,只要有降雨,根據(jù)蓄水容量分布曲線總有部分流域產(chǎn)流。但在實(shí)際工作中發(fā)現(xiàn),在降雨較小的情況下,流域基本無產(chǎn)流或產(chǎn)流較小,連續(xù)API模型中設(shè)置了初始損失降雨P(guān)s,使得在Pa較小時(shí),曲線直接交Y軸與最大初損,Pa越大,初損逐步遞減;有的流域土壤厚度較大,很難蓄滿達(dá)到全流域產(chǎn)流的情況,在連續(xù)API模型中還設(shè)置了將整個(gè)曲線族按照曲線與Y軸交線進(jìn)行旋轉(zhuǎn)的參數(shù)β,得到最后相關(guān)線族。歸納起來,P~Pa~R曲線簇由流域最大缺水(IM)、蓄水容量曲線指數(shù)(b)、流域最大損失(Ps)和旋轉(zhuǎn)角度(β)4個(gè)參數(shù)可推導(dǎo)出來。
早期的日消退系數(shù)K采用一個(gè)常值(通常取值0.85左右[4]),為了反映不同季節(jié)不同蒸散發(fā),后來多采用每月一個(gè)值,在一定程度上更趨于合理,但相當(dāng)于把原來單獨(dú)的一個(gè)參數(shù)無形中變成了12個(gè)參數(shù),增加了模型率定的難度。從K值反應(yīng)的蒸散發(fā)的含義出發(fā),力求算法簡便。K.本旺(Keith Benven)[12-13]歸納了前人的研究成果,整理得到估算潛在蒸散發(fā)的季正弦曲線公式,公式如下:
(6)
另外,程根偉[14]等提出物候指數(shù)公式:
(7)
式中,Day為時(shí)間日序列的日數(shù),即從年初開始的天數(shù),d;A、B分別為與植物有關(guān)系的系數(shù);91為水文年與自然年的日差值。
根據(jù)以上理論和API模型的日消退系數(shù)的含義,將該模型的日消退系數(shù)公式定義為:
(8)
圖2 日消退系數(shù)變化
最初的API模型中Pa值按照場次降雨進(jìn)行計(jì)算。在連續(xù)API模型中,利用計(jì)算機(jī)連續(xù)計(jì)算時(shí),對(duì)于連續(xù)分散的小雨的情況,很難確定降雨開始時(shí)間及Pa值,這也會(huì)給預(yù)報(bào)計(jì)算帶來一定的誤差。連續(xù)API模型沿用原API模型的日模型計(jì)算結(jié)構(gòu),進(jìn)行每日滾動(dòng)計(jì)算日Pa值。每日切換一個(gè)Pa值對(duì)應(yīng)的一條P~Pa~R關(guān)系曲線,每時(shí)段采用累積降雨進(jìn)行查線的方法計(jì)算。
圖3 連續(xù)API模型結(jié)構(gòu)示意
該流程利用多個(gè)數(shù)據(jù)化參數(shù),概括了連續(xù)API模型結(jié)構(gòu)和率定、計(jì)算的順序,可利用遺傳算法等AI算法進(jìn)行優(yōu)化計(jì)算,得到較為精確的模型參數(shù)。
表1 連續(xù)API模型優(yōu)化計(jì)算參數(shù)
注:采用優(yōu)化計(jì)算方法自動(dòng)率定不包括濾波系數(shù)α。
連續(xù)API模型可以應(yīng)用人工逐步率定,也可以利用優(yōu)化算法進(jìn)行自動(dòng)優(yōu)化計(jì)算。
人工交互率定可以分為兩個(gè)步驟完成。
(1) 利用人工進(jìn)行交互率定時(shí),可先通過調(diào)整濾波系數(shù)(α),將整個(gè)過程分成地下徑流和直接徑流,再調(diào)整地下水馬法參數(shù)(Kmas)和地下水入流比例系數(shù)(Z),使得模擬的地下徑流跟分割的地下徑流過程吻合且具有一定的合理性。
連續(xù)API模型11個(gè)參數(shù)(自動(dòng)率定為10個(gè)參數(shù))都已數(shù)字化,可以直接利用如遺傳算法等優(yōu)化計(jì)算方法來進(jìn)行優(yōu)化率定。優(yōu)化計(jì)算的適應(yīng)度函數(shù)除了利用預(yù)報(bào)過程的確定性系數(shù)外,還可以結(jié)合計(jì)算過程的洪峰高度和洪峰時(shí)間誤差等系數(shù)統(tǒng)一進(jìn)行約束,達(dá)到洪峰高度和長系列過程都能滿足要求的最優(yōu)參數(shù)。也可以先通過人工交互的方式調(diào)整好幾個(gè)不敏感的參數(shù)并進(jìn)行固定,其他剩下的參數(shù)再利用優(yōu)化計(jì)算法進(jìn)行優(yōu)化計(jì)算,最終得到合理的優(yōu)化模型參數(shù)。
以長江上游某流域?yàn)槔浼娣e1 441 km2,采用雨洪對(duì)應(yīng)較好的1990年洪水摘錄資料,計(jì)算時(shí)段長Δt=1 h,利用遺傳算法,設(shè)置50個(gè)基因個(gè)體,對(duì)該連續(xù)API模型進(jìn)行優(yōu)化計(jì)算。采用完全確定性系數(shù)[15]控制和多洪峰(峰時(shí)和峰量)結(jié)合確定性系數(shù)控制(三者權(quán)值分別為20%,20%和60%)兩種方法分別進(jìn)行優(yōu)化計(jì)算,每種方式各計(jì)算10次,最后對(duì)計(jì)算結(jié)果進(jìn)行分析。
分析得出,按照確定性系數(shù)控制條件每次計(jì)算得到的最低指標(biāo)都在乙級(jí)以上,按照多目標(biāo)權(quán)值控制計(jì)算模式計(jì)算結(jié)果適應(yīng)度函數(shù)平均值有所提高,但整個(gè)長系列的確定性系數(shù)稍有減小,對(duì)于洪峰的高度控制較好。計(jì)算結(jié)果見圖4,模型參數(shù)如表1所示。
圖4 連續(xù)API模型優(yōu)化計(jì)算結(jié)果
利用該表中優(yōu)化后的參數(shù),逐年對(duì)1990~2012年23 a長系列資料進(jìn)行檢驗(yàn),方案合格率達(dá)到乙級(jí)的年份所占百分比在75%以上,效果良好。
連續(xù)API模型采用11個(gè)數(shù)字參數(shù)很好地模擬了API模型的率定和實(shí)際使用計(jì)算,也為采用如遺傳算法等優(yōu)化計(jì)算方法進(jìn)行API模型參數(shù)的優(yōu)化計(jì)算提供了可能,一定程度簡化了API模型率定的過程。通過控制優(yōu)化計(jì)算的控制條件,采用傳統(tǒng)的“確定性系數(shù)”結(jié)合洪峰高度和洪峰時(shí)間的適應(yīng)度函數(shù),提高并優(yōu)化了計(jì)算的精度。
但由于連續(xù)API模型對(duì)一些參數(shù)進(jìn)行了概化,也存在以下一些不足之處。
(1) 瞬時(shí)單位線由調(diào)蓄滯時(shí)(Ku)和串聯(lián)水庫個(gè)數(shù)(N)決定,在整個(gè)汛期過程采用的都是一條單位線,不能對(duì)單位線進(jìn)行分類和分組。
(2) 公式(3)中,下滲量f和直接徑流成比列關(guān)系只是一個(gè)假設(shè),實(shí)際上f的量值既與土壤下滲特性有關(guān),還與包氣帶調(diào)蓄作用有關(guān)[3]。
(3) 該方法只能提供三參數(shù)和兩參數(shù)的API模型計(jì)算,未能考慮四變數(shù)和五變數(shù)的模型。
(4) 連續(xù)計(jì)算后,計(jì)算結(jié)果為某個(gè)流域的整個(gè)汛期產(chǎn)匯流過程,不是單場次洪水樣本,只能使用確定性系數(shù)法來評(píng)定。若使用該連續(xù)API模型還需要進(jìn)行場次評(píng)定,也可在整個(gè)長系列計(jì)算完成后,選擇每場降雨和對(duì)應(yīng)的洪水過程作為樣本,設(shè)定樣本的開始時(shí)間和結(jié)束時(shí)間對(duì)每場洪水的預(yù)見期、變幅、凈雨、洪峰和峰現(xiàn)時(shí)間等進(jìn)行場次精度評(píng)定和統(tǒng)計(jì)。
(5) 由于受到優(yōu)化計(jì)算算法的限制,每次得到的參數(shù)不一定為最優(yōu),計(jì)算完成后最好進(jìn)行合理性分析,并進(jìn)行多次計(jì)算,以選擇精度高并且又合理的參數(shù)為最后結(jié)果,在多次計(jì)算結(jié)果差別不大的情況下也可取平均值,需要具體問題具體分析。
隨著計(jì)算機(jī)技術(shù)在洪水預(yù)報(bào)中的不斷應(yīng)用,如何在原API模型中結(jié)合智能優(yōu)化算法實(shí)現(xiàn)自動(dòng)率定和參數(shù)自動(dòng)優(yōu)選,本文提出改進(jìn)后的“連續(xù)API模型”不失為一種有效、實(shí)用的方法。
[1] 袁作新.流域水文模型[M].北京:中國水利電力出版社,1990.
[2] 林三益.水文預(yù)報(bào)[M].北京:中國水利水電出版社(高等學(xué)校教材第二版),2001.
[3] 張泉生.水文預(yù)報(bào)[M].南京:河海大學(xué)出版社, 2000.
[4] 包為民.水文預(yù)報(bào)(第三版)[M].北京:中國水利水電出版社,2006.
[5] 劉國富.連續(xù)API水文模型的研究及應(yīng)用[J].大壩與安全,2008(2):4-8.
[6] 張恭肅,王成明.對(duì)API模型的改進(jìn)[J].水文.1996(4):20-25.
[7] 鐘小燕,文磊,余鐘波.連續(xù)API模型子沂河臨沂站徑流預(yù)報(bào)中的應(yīng)用[J].人民長江,2017,48(13):26-30.
[8] 熊立華,郭生練.采用非線性水庫假設(shè)的基流分割方法及應(yīng)用[J].武漢大學(xué)學(xué)報(bào)(工學(xué)版),2005,38(1):27-29.
[9] 黃國如.流量過程線的自動(dòng)分割方法探討[J].灌溉排水學(xué)報(bào),2007,26(1):73-76.
[10]徐磊磊.水文過程的基流分割方法研究進(jìn)展[J].應(yīng)用生態(tài)學(xué)報(bào),2011,22(11):3073-3080.
[11]陳文艷.數(shù)字濾波法分割基流的論證[J].水力發(fā)電,2014,40(2):37-40.
[12]Keith Beven.Rainfall-Runoff Modelling Lancaster University[M]. UK:John Wiley & Sons, Ltd, 2012.
[13]張行南,劉建芬,夏達(dá)忠.降雨-徑流模擬[M].北京:中國水利水電出版社,2016.
[14]程根偉,舒棟材.水文預(yù)報(bào)的理論與數(shù)學(xué)模型[M].北京:中國水利水電出版社,2006.
[15]GBT 22482-2008水文情報(bào)預(yù)報(bào)規(guī)范.
2017-10-10
曹 磊,男,長江水利委員會(huì)長江上游水文水資源勘測局,工程師.
1006-0081(2017)12-0031-05
P456.8
A
(編輯:朱曉紅)