汪昊燃,王 容,黃鵬年,何 健,姚 成

(1.河海大學(xué)水文與水資源學(xué)院,江蘇 南京 210098; 2.水利部水文水資源監(jiān)測預(yù)報中心,北京 100053;3.南京信息工程大學(xué)水文與水資源工程學(xué)院,江蘇 南京 210044; 4.江蘇省水文水資源勘測局,江蘇 南京 210029)

中小河流流域的洪水預(yù)報一直是防汛工作的重點[1-3]。分布式水文模型如柵格新安江模型等,對山地丘陵地形時空模擬精細(xì)化程度較高,在匯流階段一般采用馬斯京根法或擴(kuò)散波計算匯流,計算結(jié)果受河道地形、回水頂托和閘壩建筑等影響,導(dǎo)致模擬精度降低[4-11]。水力學(xué)模型對感潮河段潮位計算具有良好的模擬精度,但上游山地丘陵地區(qū)多為天然河道,地形復(fù)雜,斷面形狀不規(guī)則,造成水力學(xué)相關(guān)參數(shù)計算困難[12],因此有必要將水文學(xué)模型與水力學(xué)模型相結(jié)合,以此來建立河流與流域之間的水力聯(lián)系。殷駿[13]在南四湖區(qū)對水文和一維水力耦合預(yù)報模型的研究明確了該模型的可行性和穩(wěn)定性;在水文學(xué)和水力學(xué)模型結(jié)合的基礎(chǔ)上,李致家等[14]對比了一維水力學(xué)模型和二維水力學(xué)模型的適用情況,發(fā)現(xiàn)一維水力學(xué)模型計算精度在低水位情況下受各水體之間連通性影響較大;孫映宏等[15]根據(jù)山區(qū)、平原區(qū)、分洪區(qū)以及河道等不同的產(chǎn)匯流條件將水文學(xué)與水力學(xué)模型耦合后應(yīng)用在山前小流域,取得了良好的模擬效果。一維水力學(xué)模型在低水位計算時,由于各水體之間水流交換減少,導(dǎo)致模擬結(jié)果精度受到影響,需采用實時校正方法對誤差結(jié)果進(jìn)行修正[16]。實時校正方法一般分為過程誤差校正方法和終端誤差校正方法,在一維水力學(xué)模型的實時校正中,可采用卡爾曼濾波法對區(qū)間入流、狀態(tài)變量等過程進(jìn)行校正[17-19],也可采用K最近鄰(K-nearest neighbor,KNN)法等對計算結(jié)果進(jìn)行校正。KNN法[20]不需要建立自回歸方程,不需要計算相關(guān)參數(shù),在實時校正中有較好的應(yīng)用前景[21]。劉開磊等[22]采用KNN法對一維水力學(xué)模型計算結(jié)果進(jìn)行校正,發(fā)現(xiàn)該方法能較大程度、穩(wěn)定地提高模擬精度;反饋模擬法對于較短的預(yù)見期有良好的校正效果,且與預(yù)報模型無關(guān),同樣是終端誤差實時校正的常用方法之一[23]。

秦淮河流域河網(wǎng)密集,上游以山區(qū)為主,坡度較大;下游河段出口與長江相連,其水位受長江潮汐作用影響顯著,屬于感潮河段且同時受水閘調(diào)度的影響。本文采用柵格新安江模型與一維河網(wǎng)水力學(xué)模型結(jié)合的模型模擬水位。在流域上游山地丘陵地區(qū)采用柵格新安江模型計算出流,流域下游平原感潮河段采用一維水力學(xué)模型計算水位,并對計算結(jié)果實時校正,以期提高洪水水位模擬的精度。

1 研究區(qū)域概況

秦淮河流域位于長江下游,地形四周高、中間低,坡度較大,總面積2684km2;上游屬山地丘陵區(qū),下游地區(qū)河道窄、坡度小的地形條件導(dǎo)致流域上游匯流快,下游泄洪緩慢[24]。

東山站位于南京市江寧區(qū)東山大橋上游200m處,是秦淮河流域下游出口的水位控制站。秦淮河在東山站下游300m處分流入老秦淮河及秦淮新河。秦淮新河全長16.8km,秦淮新河閘位于秦淮新河河口處;武定門閘位于老秦淮河中段,與東山站相距12km,見圖1。武定門閘上水位站、秦淮新河閘上水位站主要為調(diào)度作用,東山站為主要負(fù)責(zé)報汛的水位控制站,因此選取東山站作為水位模擬計算的研究水位站。東山站水位受到上游來水、武定門閘、秦淮新河閘以及長江潮汐的影響,汛期水位上漲迅速,因此水位呈現(xiàn)繩套或非單一性復(fù)雜關(guān)系。

2 預(yù)報模型構(gòu)建

2.1 柵格新安江模型

柵格新安江模型產(chǎn)流部分采用蓄滿產(chǎn)流模型,計算單元柵格冠層截留量、河道降水量以及蒸散發(fā)量以推求地表徑流、壤中流和地下徑流,再將3種水源按各柵格間的演算次序演算至流域出口;模型匯流部分采用主網(wǎng)格馬斯京根法,將3種產(chǎn)流水源依次演算至站點出口[25]。整理2010—2018年東山站逐時水位資料,2014—2018年秦淮河流域的逐日降雨、徑流資料,以秦淮河水系上游前垾村(秦)站為水文模型及水力學(xué)模型的分界點,流域上游采用柵格新安江模型進(jìn)行模擬。DEM網(wǎng)格大小為1km×1km,流域內(nèi)共有9座水庫,出流均按主網(wǎng)格馬斯京根法演算至前垾村(秦)站;前垾村上游包含11個雨量站,均在流域范圍內(nèi)。

2.2 水力學(xué)模型

秦淮河干流采用一維河網(wǎng)水力學(xué)模型進(jìn)行計算。模型的上邊界采用由柵格新安江模型計算的流量過程,以流域下游出口處的長江潮位作為河網(wǎng)水力學(xué)模型的下邊界。

前垾村(秦)站下游包括秦淮新河閘上、武定門閘上、東山站3個水位站,云臺山河、牛首山河和外港河3條主要支流;秦淮河干流在東山站分為人工河道秦淮新河和老秦淮河兩條河流匯入長江(圖2)。前垾村(秦)站至下游出口共設(shè)有86個斷面,其中秦淮新河糙率取0.022,老秦淮河糙率取0.024,其余支流糙率取0.021。秦淮新河閘與武定門閘為調(diào)度水閘,調(diào)度規(guī)則均為當(dāng)水位低于6.5m時引水,當(dāng)水位高于8.5m時泄洪。前垾村(秦)站以下流域地勢低洼,存在多個圩區(qū),將圩區(qū)分為水面、水稻田、旱荒地和城市道路進(jìn)行計算;產(chǎn)流計算將土壤分為兩層,采用蓄滿產(chǎn)流方法計算;圩區(qū)內(nèi)采用單位線法計算匯流,匯流至圩區(qū)河道后,以圩區(qū)排澇模數(shù)進(jìn)行量化,作為旁側(cè)入流進(jìn)入秦淮河河道內(nèi)。

圖2 河網(wǎng)概化圖

采用的河網(wǎng)水力學(xué)方程為圣維南方程組:

(1)

(2)

式中:A為過水?dāng)嗝婷娣e;Q為流量;u為流速;Z為斷面水位;g為重力加速度;S0為水流沿程損失;Sf為摩阻比降。將式(1)用四點線性Preissmann隱式格式進(jìn)行差分,可以得到任一河段差分方程:

(3)

式中Cj、Dj、Ej、Fj、Gj、φj均為系數(shù),由初值計算。

河網(wǎng)水力學(xué)模型將河網(wǎng)分為樹狀河網(wǎng)和環(huán)狀河網(wǎng)。前垾村(秦)站至東山站之間的河道及匯入河道的支流呈樹狀河網(wǎng)水系,采用樹狀河網(wǎng)追趕法進(jìn)行河道水位計算,即對于已知的流量邊界條件,上斷面L1和下斷面L2可假設(shè)如下追趕關(guān)系:

(4)

式中:j為斷面編號(j=L1,L1+1,L1+2,…,L2-1);K、T、P、V為追趕系數(shù)。由式(2)和式(3)聯(lián)立可得遞推關(guān)系,依次求得Kj+1、Tj+1、Pj+1、Vj+1,最后得到:

QL2=PL2-VL2ZL2

(5)

式(5)與下邊界條件QL2=f(ZL2)聯(lián)解可得ZL2,依次回代可求得Zj、Qj。

環(huán)狀河網(wǎng)包括內(nèi)河道和外河道,外河道計算與樹狀河網(wǎng)計算類似,均由河道一端邊界條件已知的外節(jié)點推求另一端水力要素未知的內(nèi)節(jié)點。東山站以下的老秦淮河、秦淮新河構(gòu)成環(huán)狀河網(wǎng)水系,采用環(huán)狀河網(wǎng)雙追趕系數(shù)法進(jìn)行水位計算。內(nèi)河道的計算沒有端點邊界條件可用,河段差分方程依然可表示為式(2),首、末斷面水位為基本未知量,因此需要利用雙追趕方程求解。

3 實時校正方法

3.1 KNN法

KNN法校對計算主要包括4個步驟:①計算模擬時段前70%的誤差數(shù)據(jù)存放為全部歷史誤差數(shù)據(jù),用于以后時段的實時校正;②根據(jù)誤差自回歸理論,考慮前面多個時段(NT),設(shè)定預(yù)見期為ext,則可以從NT+ext時刻開始至m-1+NT+ext建立m個歷史預(yù)報誤差向量集;③識別最近鄰樣本,以向量間的歐式距離為評價相似相關(guān)性的指標(biāo),找出與當(dāng)前預(yù)報誤差相關(guān)向量最相似的X個預(yù)報誤差向量;④對該X個預(yù)報誤差向量對應(yīng)的預(yù)報誤差加權(quán),使用反距離權(quán)重法對當(dāng)前誤差進(jìn)行估計。

3.2 反饋法

反饋法基于洪水預(yù)報誤差的相似性,對洪水預(yù)報誤差進(jìn)行實時校正,在洪水預(yù)報實時校正方面已取得較好的結(jié)果[26]。反饋法通過計算實測流量與預(yù)報流量之間的反饋因子(F),將流量過程分為漲洪段和落洪段,進(jìn)而對預(yù)報流量進(jìn)行修正。

先從第2個時刻起計算各時刻與前一時刻的流量差:

ΔQobs,i=Qobs,i-Qobs,i-1

(6)

ΔQcal,i=Qcal,i-Qcal,i-1

(7)

式中:Qobs,i為本時刻實測流量;Qcal,i為本時刻預(yù)報流量。

通過各時刻與前一時刻的預(yù)報流量差和實測流量差,計算F:

(8)

可以得到漲水段反饋模擬實時校正計算公式為

Qobs,i=Qobs,i-1+ΔQcal,iF

(9)

退水段反饋模擬實時校正計算公式為

(10)

4 計算結(jié)果對比

選取秦淮河流域2014—2018年9場典型洪水,采用柵格新安江模型計算前垾村(秦)站上游至該站的產(chǎn)匯流過程,選定前8場進(jìn)行次洪模型參數(shù)率定,最后1場用來檢驗。部分分布參數(shù)由輸入的土壤類型、植被覆蓋類型等數(shù)據(jù)直接計算,如自由水蓄水容量SM、張力水蓄水容量WM等。部分敏感參數(shù)率定結(jié)果為:蒸散發(fā)折算系數(shù)1.00;上層張力水容量占比0.17;下層張力水容量占比0.50;地下水消退系數(shù)0.993;壤中流消退系數(shù)0.90;河網(wǎng)水流消退系數(shù)0.83。洪水模擬結(jié)果如表1所示。

表1 前垾村(秦)站洪水模擬成果統(tǒng)計

根據(jù)GB/T 22482—2008《水文情報預(yù)報規(guī)范》:模擬計算結(jié)果中,前垾村(秦)站以上流域8場洪水,徑流深合格率與洪峰合格率均為87.5%;預(yù)報檢驗的1場洪水徑流深和洪峰都合格。所選9場洪水的徑流深誤差與洪峰誤差都達(dá)到乙等精度要求。由柵格新安江模型計算的前垾村(秦)站洪水流量過程可用于一維水力學(xué)模型上邊界流量輸入過程。

選取次洪模擬中的5場洪水進(jìn)行水位模擬計算。將前垾村(秦)站的模擬流量過程作為一維水力學(xué)模型上邊界流量輸入條件,下游出口處長江潮位作為一維水力學(xué)模型的下邊界,并采用KNN法和反饋法進(jìn)行實時校正。對模擬水位的評價標(biāo)準(zhǔn)主要考慮其與實測水位的絕對誤差,因此對東山站水位計算結(jié)果的評價采用確定性系數(shù)(CNSE)和洪峰水位誤差(δPE)2種指標(biāo)。對未校正的柵格新安江-河網(wǎng)水力學(xué)模型計算成果、反饋模擬法校正的計算成果和KNN法校正的計算成果分別計算兩種指標(biāo),計算結(jié)果如圖3所示。部分模擬水位過程如圖4所示。

圖3 校正方法指標(biāo)對比

圖4 模擬水位過程

由圖3中可以看出:未校正的水位模擬計算結(jié)果存在較大的不確定性,經(jīng)過反饋法和KNN法校正后的模擬水位過程確定性系數(shù)均有較大提升,且KNN法校正效果相比反饋法更好。對于洪峰水位絕對誤差,未校正的模擬洪峰水位計算結(jié)果波動較大,經(jīng)過校正后精度得到提高,且KNN法校正后的模擬洪峰水位誤差更加穩(wěn)定,誤差值在0.04～0.55m之間,模擬精度較高。2015062508號洪水未校正的模擬水位確定性系數(shù)為0.388,該場洪水的上邊界前垾村(秦)站的入流模擬結(jié)果較好,但經(jīng)過河網(wǎng)水力學(xué)模型計算后,東山站的模擬水位較實測水位相差較大,造成確定性系數(shù)偏小。經(jīng)過調(diào)查,本場洪水在南京市區(qū)造成嚴(yán)重內(nèi)澇,因此在三汊河口處人工開鑿臨時泄洪道,對下邊界條件產(chǎn)生影響,導(dǎo)致高水位時模擬結(jié)果與實測結(jié)果誤差較大。同時根據(jù)調(diào)度規(guī)則,下游兩處水位站會在水位高于8.5m時泄洪,因此模擬水位在8.5m左右時會出現(xiàn)模擬調(diào)度實際未調(diào)度或模擬未調(diào)度實際調(diào)度的情況,造成在低水位時模擬結(jié)果不佳。經(jīng)過實時校正后,2015062508號洪水洪峰/水位模擬結(jié)果和確定性系數(shù)均大幅度提高,可以滿足模擬精度要求。

2017060920號洪水未校正的模擬水位誤差只有0.13m,但確定性系數(shù)僅為0.322。同時根據(jù)水位過程線可知,由于模擬的水位過程漲水階段相較于實測水位過程滯后約12 h,造成確定性系數(shù)偏小。分析前垾村(秦)站入流情況可知,在洪水發(fā)生初始階段該站流量為負(fù)值,即出現(xiàn)逆流情況,而柵格新安江模型模擬流量過程不會出現(xiàn)逆流情況,因此上邊界條件出現(xiàn)較大誤差,造成模擬水位過程與實際水位過程在時間尺度出現(xiàn)偏移。經(jīng)過KNN法校正后的水位模擬過程在最高水位時出現(xiàn)了一些精度損失,而峰現(xiàn)時差以及確定性系數(shù)有了較大的提高。

在實時校正前,柵格新安江模型與河網(wǎng)水力學(xué)模型結(jié)合的預(yù)報方法對水位的模擬計算結(jié)果受平原圩區(qū)地形影響較大,且因模型的局限性導(dǎo)致在低水位時模擬結(jié)果不佳,洪峰水位模擬精度較高但峰現(xiàn)時間有較大誤差。而對模擬結(jié)果采用KNN法和反饋法進(jìn)行實時校正后,均能有效提高東山站水位模擬精度。對比2種實時校正方法,反饋法計算簡單,KNN法可以根據(jù)實際預(yù)報需要靈活調(diào)整參數(shù),且校正后的水位模擬結(jié)果相較于反饋法精度更高。

5 結(jié) 語

秦淮河流域地勢四周高中間低,形成了許多圩區(qū),導(dǎo)致上游匯流快,下游泄洪慢,且受水閘泵站和長江潮汐的影響,洪水預(yù)報工作具有不確定性和復(fù)雜性。本文構(gòu)建以柵格新安江模型計算的前垾村(秦)站流量過程為一維河網(wǎng)水力學(xué)模型上邊界、以下游出口處長江潮位為模型下邊界的水位模擬計算模型,選取秦淮河流域2014—2017年的5場次洪水進(jìn)行水位模擬計算。結(jié)果表明,河網(wǎng)水力學(xué)模型對平原圩區(qū)河道匯流計算具有較高的精度,結(jié)合柵格新安江模型的預(yù)報徑流過程作為上邊界輸入條件后,實現(xiàn)了對受水工建筑物影響的感潮河段水位的模擬計算。在經(jīng)過實時校正之后,解決了低水位情況下水位預(yù)報精度不高的問題,其中KNN法校正后確定性系數(shù)為0.718～0.975,達(dá)到了乙等預(yù)報精度要求,可為秦淮河流域防汛提供參考。由于本文所選洪水場次較少,還需要在收集更多資料后對該模型進(jìn)行進(jìn)一步計算驗證和校正,逐步形成預(yù)報可靠性較高的秦淮河流域柵格新安江與河網(wǎng)水力學(xué)結(jié)合的預(yù)報模型。