馬海波 董曉華 劉 冀
(三峽大學 水利與環(huán)境學院,湖北 宜昌 443002)
瞬時單位線由納希1957年提出,指在流域上分布均勻的一個單位瞬時脈沖入流,在流域出口斷面所形成的地面徑流過程線,也就是流域降雨徑流系統的單位脈沖響應函數[1-6].與經驗單位線相比,瞬時單位線理論性強,數學推導嚴密,在一定程度上反映了流域的匯流特性,且該單位線只有兩個參數,并且適用于資料缺乏地區(qū),從而利于進行流域降雨徑流關系的理論研究,且模型的地區(qū)綜合實用性非常強,目前在國內外許多流域都得到了廣泛的應用[5,7-8].
目前納希瞬時單位線的推求方法主要是采用統計數學中的矩法,由于矩法本身受計算時段長度的影響較大,導致利用該方法計算出來的瞬時單位線經過S曲線轉換得到的時段單位線,進行還原洪水計算時,計算得到的洪水過程與實測的洪水過程差別較大,這在很大程度上限制了瞬時單位線的進一步應用[8-9].為了解決這一問題,目前有學者提出了應用遺傳算法求解瞬時單位線的參數,計算結果表明,利用這種方法得到的瞬時單位線對地面徑流過程的還原效果要優(yōu)于應用矩法推得的瞬時單位線還原效果.但是較之矩法,遺傳算法的計算繁瑣,且計算時間長,計算效率不高.為了解決目前單位線優(yōu)化方法精度、效率不高的問題,本研究擬建立一套對瞬時單位線參數優(yōu)化的方法:即先由矩法分析出參數的大致取值,然后再由梯度搜索法在參數的大致取值范圍內對參數進行尋優(yōu),最后通過實例分析驗證本方法的可靠性.
目前,流域洪水預報中常用的納西瞬時單位線模型為[1,10-11]
式中,u(t)為瞬時單位線與時間t相對應的縱坐標;t為時間;n,k為模型的兩個參數,n的物理意義為瞬時單位線定義設定的線性調解水庫的個數;K的物理意義為瞬時單位線定義的線性水庫的傳播時間;Γ(n)為伽馬函數,其定義為Γ(n)=dt;e為自然對數的底.
實際應用中,一般是先把瞬時單位線u(0,t)轉化為S(t)曲線,然后由S(t)曲線推求無因次時段單位線u(Δt,t),再由無因次時段單位線u(Δt,t)推求時段單位線q(Δt,t),具體推算公式為
本次研究中,對瞬時單位線積分在Matlab平臺上實現,Matlab函數庫提供了3種用于數值積分的函數[12],分別為梯形法積分函數(trapz)、Simpson遞歸法積分函數(quad)以及Lobatto積分法函數(quadl),其中Lobatto積分法較之前兩種積分方法精度更高.故本次研究采用Lobatto積分法進行計算,具體算法實現為
1)建立瞬時單位線函數
function y=iuh(u)
n=2.07;
y=1/gamma(n)*u.^(n-1).*exp(-u)
2)對瞬時單位線進行積分
st(j,1)= quadl('iuh',a,b)
式中,j為單位線時段數,a,b為積分上下限.當S(j,1)≥0.999 9時,可停止計算,此時,由于S曲線的最大值為1[6],因此當積分計算停止時,若S曲線的最大值為1,則此積分結果即為S曲線;若S曲線的最大值不為1,則可在S曲線的尾部再加上一個數值1即得到最終的瞬時單位線的S曲線.
式中,S(t-Δt)是由S(t)向后錯Δt時間(凈雨時段長).
式中,q(Δt,t)為時段長為Δt的時段單位線(m3/s);F為流域面積(km2).
設h(t)為地面凈雨過程,則由(4)式求得的時段單位線,可還原地面徑流過程Q(t)為
式中,hj為第j個時段的凈雨,k1、k2為對流域出口斷面流量有影響的凈雨時段的累積上下限,其值取決于地面徑流時序t與凈雨時段數m以及單位線時段數n之間的大小關系,即
由于實際洪水預報中,一般對洪峰點附近的預報精度要求較高,因此,本研究以流域出口斷面實測流量與用瞬時單位線還原流量之間的加權殘差絕對值最小為原則,建立最優(yōu)化問題來估計瞬時單位線的兩個參數,該最優(yōu)化問題可以表述為
式中,Qc(t)為用瞬時單位線還原計算得到的流域出口斷面流量過程,Qs(t)為實測的流域出口斷面流量過程,w(t)為權函數.
該優(yōu)化問題用常規(guī)方法處理很困難,首先借助矩法對參數的取值范圍進行確定,然后再由梯度搜索法求解該最優(yōu)化問題[13-14].
表1 某站的凈雨過程、實測和還原徑流過程
圖1 瞬時單位線還原徑流過程圖
從表1及圖1可以看出,梯度搜索法優(yōu)化得到的瞬時單位線以及由矩法得到的瞬時單位線還原徑流過程與實測徑流過程絕對誤差加權之和分別為32和43,且在洪峰附近,梯度搜索法較之矩法還原的精度更高,矩法在洪峰處的還原的相對誤差為10.45,而梯度搜索法的相對誤差僅為3.85%,這對于在實際的洪水預報中準確預報洪峰值具有重要意義.
鑒于應用矩法得到的瞬時單位線還原洪水過程時,還原的洪水過程與同次實測洪水過程擬合效果不佳,尤其是洪峰處的擬合效果較差的現象,本文提出了參數優(yōu)選的最優(yōu)化模型,并提出了在應用矩法得出基本參數的基礎上,用梯度搜索法來求解最優(yōu)化問題的一套優(yōu)化瞬時單位線的方法.實例研究結果表明,該方法的計算精度高于矩法,能夠得到基本反映流域匯流特性的瞬時單位線.在Matlab平臺上,通過應用Lobatto積分法對瞬時單位線進行積分來求解S曲線,進而推求時段單位線.這種方法避免了用近似數值公式以及查表帶來的誤差和不便,提高了計算效率.
[1] 李麗琴,周惠成,顧妍平.流域時變瞬時單位線非線性匯流模型及應用[J].哈爾濱工業(yè)大學學報,2009(2):180-182.
[2] 王桂林,伊學農,劉遂慶.遺傳算法推求瞬時單位線參數并計算流量過程線[J].環(huán)境污染與防治,2003(6):367-369.
[3] 李志龍.新安江模型在資料缺乏的寒區(qū)流域的應用研究[D].南京:河海大學,2006.
[4] 雷 璐,宋星原,羅 鵬,等.變動單位線與納西瞬時單位線在洪水預報中的應用[J].水電能源科學,2011(3):51-53.
[5] 耿鴻江,齊松茹.瞬時單位線的Excel快速算法[J].人民長江,2003(2):6-7.
[6] 包為民.水文預報[M].北京:中國水利水電出版社,2009.
[7] 姚 成,章玉霞,李致家,等.無資料地區(qū)水文模擬及相似性分析[J].河海大學學報:自然科學版,2013(2):108-113.
[8] 金菊良,丁 晶,魏一鳴.瞬時單位線的優(yōu)化估計[J].水力發(fā)電學報,2003(1):70-75.
[9] 張 寧.基于遺傳算法的試錯法推求單位線[J].河南水利與南水北調,2013(4):25-27.
[10]詹道江,葉守澤.工程水文學[M].北京:中國水利水電出版社,2000.
[11]董四輝,周惠成.遺傳算法在估計瞬時單位線參數中的應用[J].大連鐵道學院學報,2006(4):73-77.
[12]宋葉志.MATLAB數值分析與應用[M].北京:機械工業(yè)出版社,2009.
[13]唐發(fā)明.基于統計學習理論的支持向量機算法研究[D].武漢:華中科技大學,2005.
[14]張國云.支持向量機算法及其應用研究[D].長沙:湖南大學,2006.