束慧連,包為民,2,楊?yuàn)檴?劉可新,趙麗平,楊小強(qiáng)
(1.河海大學(xué)水文水資源學(xué)院,江蘇南京210098;2.河海大學(xué)水文水資源與水利工程科學(xué)國(guó)家重點(diǎn)實(shí)驗(yàn)室,江蘇南京210098)
?
產(chǎn)流誤差動(dòng)態(tài)系統(tǒng)響應(yīng)曲線修正方法的應(yīng)用
束慧連1,包為民1,2,楊?yuàn)檴?,劉可新1,趙麗平1,楊小強(qiáng)1
(1.河海大學(xué)水文水資源學(xué)院,江蘇南京210098;2.河海大學(xué)水文水資源與水利工程科學(xué)國(guó)家重點(diǎn)實(shí)驗(yàn)室,江蘇南京210098)
摘要:動(dòng)態(tài)系統(tǒng)響應(yīng)曲線修正方法(Dynamic System Response Curve,DSRC)是一種基于微分響應(yīng)的有明確物理意義并且有效的實(shí)時(shí)修正新方法,屬于向誤差源頭進(jìn)行修正的方法。通過計(jì)算產(chǎn)流的系統(tǒng)響應(yīng)曲線對(duì)產(chǎn)流系列進(jìn)行修正,實(shí)現(xiàn)流量過程的修正。將該方法與新安江模型相結(jié)合應(yīng)用于閩江上游的水吉流域,并將修正效果與遞推最小二乘法進(jìn)行了比較。利用水吉流域29場(chǎng)洪水資料進(jìn)行模擬和修正,結(jié)果表明該方法可以有效提高洪水預(yù)報(bào)的精度,且其效果優(yōu)于遞推最小二乘法。
關(guān)鍵詞:產(chǎn)流誤差;動(dòng)態(tài)系統(tǒng)響應(yīng)曲線;水吉流域
0引言
實(shí)時(shí)洪水預(yù)報(bào)作為防洪減災(zāi)的非工程措施,有其重要的經(jīng)濟(jì)和民生意義。然而,水文模型預(yù)報(bào)洪水習(xí)慣上是采用觀測(cè)的歷史水文資料來(lái)確定模型參數(shù),然后用于未來(lái)的洪水預(yù)報(bào)中。但這樣的預(yù)報(bào)方案往往在實(shí)時(shí)洪水預(yù)報(bào)中得不到滿意的結(jié)果。因此,為了提高實(shí)時(shí)洪水預(yù)報(bào)的精度,有必要進(jìn)行實(shí)時(shí)誤差修正。所謂實(shí)時(shí)洪水預(yù)報(bào)誤差修正就是對(duì)在水文模型中沒有考慮或是考慮不得當(dāng)卻對(duì)預(yù)報(bào)結(jié)果有影響的因素,利用實(shí)測(cè)信息和其他可利用信息進(jìn)行實(shí)時(shí)校正,達(dá)到彌補(bǔ)水文模型不足的目的[1]。按修正
內(nèi)容劃分,實(shí)時(shí)修正方法可以分為:模型誤差修正,如自回歸模型[2];模型參數(shù)修正,如神經(jīng)網(wǎng)絡(luò)算法[3];模型輸入修正,如抗差估計(jì)理論[4];模型狀態(tài)修正,如卡爾曼濾波方法[5- 6];綜合修正,就是前四種方法的結(jié)合[7]。
現(xiàn)有的誤差修正方法種類很多[8-12],源于海洋與氣象學(xué)的數(shù)據(jù)同化方法也被應(yīng)用于水文預(yù)報(bào)的誤差修正中[13]。然而,常見的誤差修正方法存在各種各樣的問題,如沒有明確物理意義、預(yù)見期的損失、資料不能合理利用等,再加上很多修正方法存在假設(shè)條件,當(dāng)實(shí)際情況不能滿足假設(shè)條件的時(shí)候修正效果顯然不會(huì)理想。例如:卡爾曼濾波方法雖然在修正過程中明確地考慮了模型和資料的不確定性,但是假設(shè)條件不易得到滿足,應(yīng)用于分布式水文模型、多參數(shù)模型等時(shí)計(jì)算量大,結(jié)構(gòu)復(fù)雜,實(shí)時(shí)洪水預(yù)報(bào)難以提供誤差修正所需的信息量;自回歸模型假設(shè)誤差序列按時(shí)間排序是相關(guān)的,但是當(dāng)預(yù)報(bào)量發(fā)生突變時(shí)破壞了相關(guān)關(guān)系,因此不再適用,修正能力較弱,而且隨著預(yù)見期的增加修正效果迅速下降;神經(jīng)網(wǎng)絡(luò)算法雖然能模擬變量間的非線性關(guān)系,但其結(jié)構(gòu)不唯一,且收斂極為緩慢。
動(dòng)態(tài)響應(yīng)曲線修正方法(Dynamic System Response Curve,DSRC)是一種追溯至誤差源頭的修正方法,沒有采用傳統(tǒng)的誤差修正方法的修正模式——在預(yù)報(bào)結(jié)果中進(jìn)行誤差修正;而是通過修正用于進(jìn)行洪水預(yù)報(bào)的資料達(dá)到誤差修正的目的。新安江模型分為四個(gè)部分,產(chǎn)流誤差動(dòng)態(tài)響應(yīng)曲線修正方法將產(chǎn)流以下的部分作為響應(yīng)系統(tǒng),并用差分近似代替微分,得到時(shí)段產(chǎn)流量所對(duì)應(yīng)的系統(tǒng)響應(yīng)曲線,將實(shí)測(cè)流量和計(jì)算流量的差值作為修正信息,用最小二乘估計(jì)修正產(chǎn)流過程,達(dá)到修正計(jì)算流量過程的目的[14]。本文在閩江上游的水吉流域?qū)⒃撜`差修正方法與新安江模型相結(jié)合,比較分析了修正前后的結(jié)果,驗(yàn)證了該方法的可行性和有效性。
1修正方法介紹
圖1 系統(tǒng)示意
新安江模型將模型結(jié)構(gòu)分散為蒸散發(fā)計(jì)算、產(chǎn)流計(jì)算、分水源計(jì)算和匯流計(jì)算四個(gè)層次,匯流計(jì)算又包含坡地匯流和河道匯流兩部分,坡地匯流計(jì)算采用線性水庫(kù),河道匯流采用馬斯京根分段連續(xù)演算。文獻(xiàn)[14]提出將產(chǎn)流誤差修正的動(dòng)態(tài)系統(tǒng)響應(yīng)曲線修正方法是將新安江模型的分水源計(jì)算和匯流計(jì)算這兩個(gè)層次作為一個(gè)系統(tǒng),如圖1所示,系統(tǒng)響應(yīng)曲線由線性差分近似代替微分得到,修正信息為實(shí)測(cè)流量系列,采用最小二乘估計(jì)原理,對(duì)產(chǎn)流量R進(jìn)行修正,最終實(shí)現(xiàn)計(jì)算流量QC的修正。
1.1動(dòng)態(tài)系統(tǒng)響應(yīng)曲線
由文獻(xiàn)[14]可知,所謂動(dòng)態(tài)系統(tǒng)響應(yīng)曲線,就是給某時(shí)段的產(chǎn)流量ri增加一個(gè)單位的值,用變化后的產(chǎn)流系列帶入模型計(jì)算得到新的流量過程,這個(gè)新的流量過程與改變之前得到的模型計(jì)算流量過程Qc相減得到的曲線就是對(duì)應(yīng)ri的系統(tǒng)響應(yīng)曲線。由于該曲線受到其他時(shí)段產(chǎn)流量的影響實(shí)時(shí)變化,因此稱作動(dòng)態(tài)系統(tǒng)響應(yīng)曲線。文獻(xiàn)[14]對(duì)其計(jì)算過程進(jìn)行了推導(dǎo)。
方法的關(guān)鍵是U矩陣的求解,即
(1)
式中,RC為模型初始計(jì)算徑流系列[rc1,…,rcn];t為時(shí)間;θ為模型參數(shù);L為樣本長(zhǎng)度;Q(RC,θ,t)為模型初始計(jì)算流量過程QC。
U矩陣的每一項(xiàng)都可以通過線性差分來(lái)代替微分近似求得
(t
(2)
得到U矩陣的具體計(jì)算步驟為:
(1)首先直接通過新安江模型進(jìn)行計(jì)算,得到初始的計(jì)算徑流過程RC和流量過程QC。
(2)在以計(jì)算得到的徑流系列的基礎(chǔ)上給i時(shí)段的徑流量增加1個(gè)單位值,將變化后的徑流系列帶入新安江模型的分水源計(jì)算和匯流計(jì)算中,求得對(duì)應(yīng)的流量過程。
(3)用步驟(2)中得到的流量過程與初始計(jì)算流量過程相減,得到對(duì)應(yīng)時(shí)段徑流的響應(yīng)曲線,即U矩陣中的第i列。
(4)完成所有列的計(jì)算則得到U矩陣。
計(jì)算得到U矩陣后,產(chǎn)流量的修正量ΔR通過最小二乘法得到,即
ΔR=(UTU)-1UT(Q(R,θ,L)-Q(RC,θ,t))
(3)
將修正后的產(chǎn)流量R′=RC+ΔR帶入新安江模型,計(jì)算得到修正后的計(jì)算流量過程,可以實(shí)現(xiàn)誤差修正。
1.2遞推最小二乘法
描述時(shí)間系列的線性回歸模型[7],其觀測(cè)系列為(x11,x12,…,x1p;y1),(x21,x22,…,x2p;y2),…,(xt1,xt2,…,xtp;yt)。其中,x為自變量;y為因變量;t為時(shí)刻;p為模型回歸階數(shù)。若已知t時(shí)刻之前的自變量和因變量的觀測(cè)值,要估計(jì)t+1時(shí)刻的因變量,可根據(jù)這些觀測(cè)信息先用最小二乘法估計(jì)出參數(shù),然后再估計(jì)y值。遞推最小二乘法可以解決當(dāng)中計(jì)算量的問題和隨著信息量的增大計(jì)算機(jī)容量的問題。
將每個(gè)觀測(cè)值帶入回歸方程
(4)
寫成矩陣形式為
Yt=XtCt+Et
(5)
表1水吉流域使用新安江模型的參數(shù)
WMUMLMKCBIMCSMEXKGKICGCICSKEXE16320801.350.7900.1613.751.50.450.40.9940.920.7451.99770.2175
注:表中WM為流域平均張力水容量,mm;UM為上土層張力水容量,mm;LM為下土層張力水容量,mm;KC為蒸散發(fā)折算系數(shù);B為蓄水容量面積分布曲線指數(shù);IM為不透水面積占全流域面積比例;C為深層蒸散發(fā)擴(kuò)散系數(shù);SM為自由水蓄水容量,mm;EX為自由水蓄水容量面積分布曲線指數(shù);KG和KI分別為自由水蓄水庫(kù)對(duì)地下水和壤中流的日出流系數(shù);CG、CI、CS分別為地下水、壤中流、地面徑流;KE和XE為馬斯京根法參數(shù),其中KE單位為h。
在t-1時(shí)刻和t時(shí)刻可得參數(shù)的最小二乘估計(jì)
(6)
(7)
經(jīng)過一系列的簡(jiǎn)化和計(jì)算可以得到參數(shù)的遞推最小二乘估計(jì)
(8)
2應(yīng)用與檢驗(yàn)
2.1流域簡(jiǎn)介
閩江流域處于亞熱帶季風(fēng)氣候區(qū),雨量豐沛,氣候溫和。流域多年平均降水量1 724 mm,年內(nèi)分布不均,年際變化大,其中60%左右的年降水量主要集中在4月~6月,豐水年降雨量是枯水年的2~3倍。受氣候的影響,流域內(nèi)多年平均蒸發(fā)量為915.0 mm[15],夏季蒸發(fā)量較大,冬季蒸發(fā)量較小。
水吉流域地處閩江上游地區(qū),南平以上閩江支流南浦溪流經(jīng)該流域,流域面積3 767 km2,多年平均氣溫15~20℃。本文選用水吉流域1988年~1999年之間的資料進(jìn)行分析,流域年平均降水量為1 820 mm,年平均徑流深972.0 mm[16],徑流系數(shù)大于0.4,屬于典型的濕潤(rùn)地區(qū),使用三水源新安江模型進(jìn)行洪水模擬預(yù)報(bào)是較為合理的。
水吉流域境內(nèi)有11個(gè)雨量站,分別為:游楓、東坑、浦城、仙陽(yáng)、高坊、山下、洋源、舊館、管坦、水吉、外屯,有一個(gè)蒸發(fā)站,以水吉站作為預(yù)報(bào)站。
2.2應(yīng)用檢驗(yàn)
本文選用水吉流域1988年~1999年間29場(chǎng)洪水資料,使用新安江模型進(jìn)行模擬,模型使用的參數(shù)采用日模參數(shù)率定獲得的結(jié)果,具體見表1。
驗(yàn)證產(chǎn)流誤差的動(dòng)態(tài)響應(yīng)曲線修正方法的修正效果,具體通過以下3項(xiàng)指標(biāo)來(lái)衡量修正效果:
(1)確定性系數(shù)DC
(9)
(2)洪峰相對(duì)誤差ΔQP
ΔQP=((QCP-QOP)/QOP)×100%
(10)
(3)洪量相對(duì)誤差ΔQT
ΔQT=((RC-RO)/RO)×100%
(11)
29場(chǎng)洪水修正前和修正后的結(jié)果比較及修正效果見表2、表3。43880620和43960317場(chǎng)次洪水的修正結(jié)果如圖2和圖3所示。
表2水吉流域29場(chǎng)洪水實(shí)時(shí)修正前后計(jì)算結(jié)果
洪號(hào)實(shí)測(cè)值ROP/mmQOP/m3·s-1新安江模型ΔQP/%ΔRT/%DC動(dòng)態(tài)系統(tǒng)響應(yīng)曲線法ΔQP/%ΔRT/%DC遞推最小二乘法ΔQP/%ΔRT/%DC4388022876.61810-15.91-4.230.914-12.451.50.974-5.841.440.9404388050934.98313.981.180.7231.411.430.8285.02-0.320.7904388062094.82680-32.57-4.270.6787.112.080.9145.65-4.70.8264389061846.52020-26.5-20.710.8382.652.740.993-9.55-6.980.94543890629119.614109.8319.940.761-0.2-1.660.9402.763.480.8894390060123.4636-32.360.450.663-0.10.310.9770.83-1.550.83543900608114.12085-19.63-14.940.864-17.83-1.270.924-18.53-0.260.8694391032651.97442.0829.860.225-1.04-3.310.897-0.951.90.8264391042634.5540.0618.8311.80.5037.12-4.370.65221.8614.890.50943920321130.3101012.8819.70.6717.659.030.75411.99.390.8304392050125.458030.374.610.424-4.18-0.430.9920.79-3.60.92043920616542010-0.32.530.953-0.46-1.80.9916.152.930.93343920704158.434702.168.80.9540.2-1.040.982-6.281.780.9874393050280.2112026.21-1.820.80824.921.660.87618.333.20.90643930615223.923603.873.950.9463.663.070.9524.812.70.9574393063094.42380-27.66-13.420.8222.690.440.99312.7-0.770.9834394042562.21410-31.52-15.430.748-0.53-0.660.94117.320.480.94943940614127.73250-1.89-15.970.865-0.011.750.9814.051.630.9444395060346.31896.675.760.640.9140.620.250.99512.311.080.9804395061485.4175010.79-11.430.897-0.67-2.520.9889.28-0.090.9784395081312.622524.7924.16-0.1311.320.740.88230.121.130.93043960317791190-8.1915.070.7030.515.540.94514.582.630.9234397060554172023.0919.760.7614.79-0.750.9606.111.520.94543980301114.61196.019.314.860.9507.46-1.380.96411.12-0.330.9744398050967.7201012.4919.440.9317.22-0.050.9627.413.310.9664399041554.41840-15.444.310.930-14.143.430.9452.261.120.94043990515661177.14-0.2719.990.849-8.7-1.160.966-1.671.030.9494399071549.7940-12.63-19.890.558-2.27-0.740.90412.26-1.210.9134399082591.51390-0.0820.150.820-5.970.010.9071.115.770.904
表3水吉流域29場(chǎng)洪水實(shí)時(shí)修正平均修正效果
評(píng)價(jià)模型評(píng)價(jià)指標(biāo)ΔQp/%ΔR/%DC新安江模型14.5312.180.743動(dòng)態(tài)系統(tǒng)響應(yīng)曲線法5.101.900.930遞推最小二乘法9.022.800.905
圖2 43880620場(chǎng)次洪水的誤差修正結(jié)果
圖3 43960317場(chǎng)次洪水的誤差修正結(jié)果
從圖2、3和表2、3可以看出:采用基于動(dòng)態(tài)系統(tǒng)響應(yīng)曲線的產(chǎn)流誤差修正方法之后,洪峰相對(duì)誤差、洪量相對(duì)誤差均有所減小,確定性系數(shù)均有所提高,平均修正效果也好于遞推最小二乘法。計(jì)算過程中發(fā)現(xiàn),采用遞推最小二乘法修正時(shí)會(huì)出現(xiàn)時(shí)刻流量值突漲的現(xiàn)象。如43960317場(chǎng)次的洪水,修正后雖然各項(xiàng)指標(biāo)有所提高但洪峰偏大且峰現(xiàn)時(shí)間提前,實(shí)際上修正后并不符合要求。而動(dòng)態(tài)系統(tǒng)響應(yīng)曲線方法則沒有這一缺點(diǎn)。從后者模擬效果得到了明顯的提升可以認(rèn)為,產(chǎn)流誤差動(dòng)態(tài)響應(yīng)曲線修正方法可以有效應(yīng)用于水吉流域。
3結(jié)論
產(chǎn)流誤差動(dòng)態(tài)系統(tǒng)響應(yīng)曲線修正方法是一種新的實(shí)時(shí)校正方法,其物理意義明確、結(jié)構(gòu)簡(jiǎn)明,既不損失預(yù)見期也不增加模型的參數(shù)個(gè)數(shù),且能有效保證預(yù)報(bào)精度。在水吉流域的應(yīng)用進(jìn)一步驗(yàn)證了該方法的合理性和實(shí)用性,將29場(chǎng)洪水的洪峰平均相對(duì)誤差從14.53%降低到5.10%,洪量平均相對(duì)誤差從12.18%降低到1.90%,平均確定性系數(shù)從0.743提高到0.930。與遞推最小二乘法比較也體現(xiàn)了其優(yōu)越性,整體修正效果相對(duì)較好。雖然該方法也存在問題,如U矩陣的穩(wěn)定性問題等,但值得在更多流域進(jìn)行應(yīng)用驗(yàn)證推廣和改進(jìn),以推進(jìn)實(shí)時(shí)洪水預(yù)報(bào)誤差修正方法的發(fā)展。
參考文獻(xiàn):
[1]包為民. 水文預(yù)報(bào)[M]. 北京: 中國(guó)水利水電出版社, 2006.
[2]周軼, 李致家. 改進(jìn)最小二乘遞推算法的洪水預(yù)報(bào)應(yīng)用研究[J]. 水力發(fā)電, 2006, 32(8): 14- 16.
[3]覃光華, 丁晶. 基于人工神經(jīng)網(wǎng)絡(luò)的卡爾曼濾波實(shí)時(shí)校正技術(shù)[J]. 水力發(fā)電, 2002(11): 9- 12.
[4]包為民, 嵇海祥, 胡其美, 等. 抗差理論及在水文學(xué)中的應(yīng)用[J]. 水科學(xué)進(jìn)展, 2003, 14(4): 528- 532.
[5]WU Xiaoling, WANG Chuanhai, CHEN Xi, et al. Kalman Filtering Correction in real-time forecasting with hydrodynamic model[J]. Journal of Hydrodynamics, Ser. B, 2008, 20(3): 391- 397.
[6]KOMMA J, BL?SCHL G, RESZLER C. Soil moisture updating by Ensemble Kalman Filtering in real-time flood forecasting[J]. Journal of Hydrology, 2008, 357(3): 228- 242.
[7]瞿思敏, 包為民. 實(shí)時(shí)洪水預(yù)報(bào)綜合修正方法初探[J]. 水科學(xué)進(jìn)展, 2003, 14(2): 167- 171.
[8]HINO M. Runoff forecasts by linear predictive filter[J]. Journal of the Hydraulics Division, 1970, 96(3): 681- 702.
[9]BIONDI D, De LUCA D L. A Bayesian approach for real-time flood forecasting[J]. Physics and Chemistry of the Earth, Parts A/B/C, 2012, 42/44: 91- 97.
[10]HSU M H, FU J C, LIU W C. Flood routing with real-time stage correction method for flash flood forecasting in the Tanshui River, Taiwan[J]. Journal of Hydrology, 2003, 283(1): 267- 280.
[11]MADSEN H, SKOTNER C. Adaptive state updating in real-time river flow forecasting—A combined filtering and error forecasting procedure[J]. Journal of Hydrology, 2005, 308(1): 302- 312.
[12]PAGANO T C, WANG Q J, HAPUARACHCHI P, et al. A dual-pass error-correction technique for forecasting streamflow[J]. Journal of Hydrology, 2011, 405(3): 367- 381.
[13]YU Z, LIU D, Lü H, et al. A multi-layer soil moisture data assimilation using support vector machines and ensemble particle filter[J]. Journal of Hydrology, 2012, 475: 53- 64.
[14]司偉, 包為民, 瞿思敏. 洪水預(yù)報(bào)產(chǎn)流誤差的動(dòng)態(tài)系統(tǒng)響應(yīng)曲線修正方法[J]. 水科學(xué)進(jìn)展, 2013, 24(004): 497- 503.
[15]林虹. 閩江流域?qū)崟r(shí)洪水預(yù)報(bào)系統(tǒng)模型研究[D]. 南京: 河海大學(xué), 2005.
[16]張章新. 閩江流域水文特性分析[J]. 水文, 2000, 20(6): 55- 58.
(責(zé)任編輯陳萍)
收稿日期:2015- 06- 23
基金項(xiàng)目:國(guó)家自然科學(xué)基金面上項(xiàng)目(51279057/41371048/51479062);國(guó)家自然科學(xué)基金重大項(xiàng)目(51190090);國(guó)家自然科學(xué)基金重大項(xiàng)目子項(xiàng)目(51190091);中央高?;究蒲袠I(yè)務(wù)費(fèi)專項(xiàng)資金資助(2014B35314)
作者簡(jiǎn)介:束慧連(1989—),女,江蘇鹽城人,碩士研究生,主要從事水文物理規(guī)律模擬、水文預(yù)報(bào)方面研究.
中圖分類號(hào):P338(257)
文獻(xiàn)標(biāo)識(shí)碼:A
文章編號(hào):0559- 9342(2016)04- 0022- 05
Application of Runoff Error Correction Based on Dynamic System Response Curve
SHU Huilian1, BAO Weimin1,2, YANG Shanshan1, LIU Kexin1, ZHAO Liping1, YANG Xiaoqiang1
(1. College of Hydrology and Water Resources, Hohai University, Nanjing 210098, Jiangsu, China; 2. State Key Laboratory of Hydrology-Water Resources and Hydraulic Engineering, Hohai University, Nanjing 210098, Jiangsu, China)
Abstract:The correction method of Dynamic System Response Curve (DSRC) is a new real-time correction method that is based on differentiation and is with clear physical meaning. This method can trace back to the source of error. The dynamic system response curve of runoff is calculated and then it is used to correct the runoff, and eventually, the correction of flood process is realized. To verify the correction effectiveness of DSRC method, it is combined with the Xin’anjiang model and is tested in Shuiji Basin of Minjiang River, and the method is also compared to the recursive least square method. Total 29 floods in Shuiji Basin are calculated and corrected, and the results of correction show that the correction method of dynamic system response curve has a remarkable effect in improving the accuracy of flood forecasting and it is more effective than the recursive least square method.
Key Words:runoff error; dynamic system response curve; Shuiji Basin