武桂芝 楊帆 黃明翔

摘要：針對北方干旱或半干旱區(qū)季節(jié)性河道滲漏現(xiàn)象嚴(yán)重問題，為提高流域洪水演進(jìn)模擬精度，在傳統(tǒng)水動力學(xué)洪水演進(jìn)模型基礎(chǔ)上耦合河道滲漏項，構(gòu)建季節(jié)性河流洪水演進(jìn)數(shù)學(xué)模型。采用霍頓入滲模型描述季節(jié)性河道滲漏過程，解決河道滲漏量計算及與基于圣維南方程組洪水演進(jìn)模型耦合問題；并利用Preissmann四點偏心隱格式與追趕法對模型離散和計算，開發(fā)FORTRAN計算程序；最后以大沽河流域?qū)崪y資料進(jìn)行驗證。結(jié)果表明，沙灣莊和南村斷面水位與流量過程的模擬值與實測值吻合且變化趨勢一致，洪水峰面到達(dá)斷面模擬時間和實際時間基本一致，研究區(qū)各河段來水量、滲漏量的實測值與模擬值相對誤差范圍分別為-166%～036 %、-309%～372 %，所建模型可滿足在現(xiàn)狀條件下進(jìn)行洪水模擬預(yù)報的需要。

關(guān)鍵詞：季節(jié)性河道；河道滲漏；洪水演進(jìn)；圣維南方程；耦合模型

中圖分類號：TV139.16文獻(xiàn)標(biāo)志碼：A文章編號：16721683（2018）03003305

Study on seasonal channel flood routing model coupled with leakage term

WU Guizhi，YANG Fan，HUANG Mingxiang

（School of Environmental & Municipal Engineering，Qingdao University of Technology，Qingdao 266033，China）

Abstract：In view of the serious problem of seasonal river leakage in arid or semiarid areas of northern China，in order to improve the simulation accuracy of flood routing in the basin，we established a seasonal channel flood routing model which couples the traditional hydrodynamic flood routing model with channel leakage.We used the Horton infiltration model to describe seasonal river infiltration and solved the coupling problem between river channel leakage and flood routing model based on SaintVenant equations.Then， the model was discretized and solved with the Preissmann fourpoint partialnode implicit scheme and chasing method， and the FORTRAN calculation program was developed.Finally，the model was tested with measured data of the Dagu River basin.The results indicated that the simulated values of water level and flow process in Shawan and Nancun were consistent with the measured values，and their trends of variation were also consistent.The simulated time for the flood peak to reach the crosssection was basically the same as the actual time.The relative error ranges between the measured value and the simulated value of water inflow and leakage in all the channel segments of the study area were -166%036% and -309%372% respectively.It can be concluded that this coupled model can satisfy the need of flood forecasting under the present condition.

Key words：seasonal channel；channel leakage；flood routing；SaintVenant equations；coupling model

中國北部干旱或半干旱地域河流多為季節(jié)性的河流，季節(jié)性河流汛期水量大，非汛期河道斷流、干枯或少水狀態(tài)，河床基質(zhì)滲透能力較強，行洪期河道滲漏現(xiàn)象嚴(yán)重，汛期滲漏主要影響到洪水流量與地下水增補量，洪水演進(jìn)過程和規(guī)律則發(fā)生顯著變化[13]。所以，在洪水模擬預(yù)報中必需探討河道滲漏的影響，耦合滲漏項的洪水演進(jìn)模擬成為干旱或半干旱區(qū)域洪水模擬中關(guān)鍵問題[46]。

考慮河道滲漏損失的洪水演進(jìn)模擬和預(yù)報，關(guān)鍵是要解決河道滲漏模擬、滲漏量計算及其與洪水演進(jìn)模型的結(jié)合等問題。國內(nèi)齊春英等[7]探討了明渠非恒定流在考慮入滲條件下的方程組，把水力學(xué)原理和水文學(xué)方法綜合，建立了海河南部平原有河流滲漏的洪水流量計算關(guān)系式。田守波[8]通過大沽河河道野外滲漏試驗和現(xiàn)場水文觀測，研究了河道的滲漏規(guī)律，并計算出了各個河段的穩(wěn)定滲漏通量，利用MIKE11軟件建立大沽河考慮滲漏條件下的一維河道洪水演進(jìn)模型。田月等[9]綜合玉符河特點，對非恒定流連續(xù)性方程基礎(chǔ)上引入流段內(nèi)滲漏量，構(gòu)建了考慮河床滲漏的圣維南過程，利用軟件模擬獲得研究區(qū)段內(nèi)46個典型監(jiān)測斷面的水位與流量過程。程亮等[10]構(gòu)建了基于霍頓入滲模型的河道入滲模擬辦法，且將入滲作為單位區(qū)間出流，和以馬斯京根法為基礎(chǔ)的河流洪水演進(jìn)數(shù)學(xué)模型結(jié)合，建立了強烈入滲條件的洪水演進(jìn)模型。王宗志等[11]采用質(zhì)量守恒原理和微元分析方法推導(dǎo)含下滲項的圣維南方程組，并以此作為描述河網(wǎng)洪水運動的控制方程，以有限差分方法作為控制方程的離散格式，采用水文學(xué)與水力學(xué)相結(jié)合的途徑，耦合上下邊界，以及內(nèi)邊界處理模擬模塊，建立考慮河道下滲的河網(wǎng)洪水模擬模型。在國外，Mudd等[12]以穩(wěn)定入滲率模擬季節(jié)性河道入滲，建立和數(shù)值模擬了考慮滲漏損失條件下季節(jié)性河道洪水演進(jìn)模型，用縮放比例法研究了水深、洪水峰值、水文過程、河道寬度等因素對河道滲漏以及洪水演進(jìn)的影響。Efrat Morin等[13]研究了干旱區(qū)干流入滲對洪水演進(jìn)的作用、洪水運動過程對河流滲漏的作用等。

第16卷 總第96期·南水北調(diào)與水利科技·2018年6月武桂芝等·耦合滲漏項的季節(jié)性河道洪水演進(jìn)模型研究國內(nèi)外一些學(xué)者在水文學(xué)、水力學(xué)及水文與水力學(xué)結(jié)合的方法基礎(chǔ)上建立了考慮河道入滲的洪水演進(jìn)模型，一定程度上可更好地模擬洪水運動過程，彌補傳統(tǒng)洪水演進(jìn)模型的不足。但是，針對我國北方季節(jié)性河流滲漏特點，河道入滲模擬以及耦合滲漏項的洪水演進(jìn)模型研究等方面還不夠完善。因此，前期課題組武桂芝等[1415]和馮增帥等[1617]對青島大沽河滲漏及對洪水演進(jìn)影響方面進(jìn)行研究，如開展野外入滲試驗、數(shù)值模擬計算簡單入滲因素對洪水演進(jìn)影響等，為本次建立滲漏模型研究開展奠定了一定基礎(chǔ)。本文遴選了霍頓入滲數(shù)學(xué)模型模擬季節(jié)性河流滲漏過程，基于流體力學(xué)原理建立耦合滲漏項的季節(jié)性河流洪水演進(jìn)數(shù)學(xué)模型，采取了Preissmann四點偏心隱格式與追趕法對方程組離散和計算，并且利用大沽河流域?qū)崪y資料對模型驗證。

1耦合滲漏項的洪水演進(jìn)模型建立與求解

1.1耦合滲漏項的圣維南方程組建立

洪水波在河流中的演進(jìn)屬于非恒定明渠漸變流，非恒定明渠基本方程即為圣維南（SaintVenant）方程組，包含連續(xù)性方程與運動方程[1820]。本文從不可壓縮流體運動的相關(guān)基本原理出發(fā)，利用了流體力學(xué)中質(zhì)量守恒原理的連續(xù)性微分方程和牛頓第二力學(xué)定律的運動性微分方程分析構(gòu)建了耦合河道滲漏項的洪水演進(jìn)數(shù)學(xué)方程組。

連續(xù)性方程：

BZt+Qx=-f·B（1）

運動方程：

Qt+x Q2A+gAZx=-gASf-f2BS0（2）

式中：B為河道斷面寬度；Z為水位變量；Q為流量；x，t分別為距離和時間的坐標(biāo)；f為滲漏項入滲率；A為過水?dāng)嗝婷娣e；g為重力加速度；Sf為摩阻坡度；S0為河道底坡，指河底的縱向坡度。

由式（1）、式（2）可知，入滲率是耦合滲漏項的洪水演進(jìn)模擬的關(guān)鍵。本文利用季節(jié)性河流大沽河上游馬家會村、中游馬軍寨村、下游賈疃村三處野外垂向入滲試驗數(shù)據(jù)，對比了Horton、Kostiakov和Philip 三種模型，其中Horton模型對描述入滲過程效果最好，其參數(shù)具有較強的物理意義，且便于和洪水演進(jìn)數(shù)學(xué)模型進(jìn)行耦合。

入滲率公式：

f（t）=fc+（f0-fc）exp（-β·t）（3）

式中：f（t）為入滲率；t為入滲時間；fc為穩(wěn)定入滲率；f0為初始入滲率；β為入滲參數(shù)。式（1）至式（3）組成了耦合河道入滲項的洪水演進(jìn)方程組。

1.2耦合滲漏項的方程組求解

四點隱式差分格式計算過程雖較復(fù)雜，但解法是無條件穩(wěn)定的，穩(wěn)定性與精度較好，可設(shè)置較大的步長，計算速度較快，是目前普遍采納的方法[2123]。采用Preissmann四點偏心隱格式對含有滲漏項的連續(xù)性與運動性方程離散，得到以增量描述的非線性數(shù)學(xué)關(guān)系式，忽略二階微量，簡化成線性代數(shù)關(guān)系式，再利用追趕法計算出結(jié)果，見圖1。

圖1四點隱式格式差分網(wǎng)格圖

Fig.1Fourpoint implicit difference grids

簡化線性隱式格式對關(guān)系式推導(dǎo)如下：

f（x，t）=12（fnj+1+fnj）（4）

fx=θΔfj+1-ΔfjΔx+fnj+1-fnjΔx（5）

ft=Δfj+1-Δfj2Δt（6）

式中：上標(biāo)n，n+1為時間標(biāo)識；下標(biāo)j，j+1為空間標(biāo)識；θ為權(quán)重系數(shù)，且0≤θ≤1。

將式（4）至式（6）的差分關(guān)系式依次代入式（1）和式（2）中，并對其線性化處置后，最終得到水流連續(xù)性方程和運動性方程式的差分關(guān)系式：

ΔQj+1-ΔQj+CjΔZj+1+CjΔZj=Dj（7）

EjΔQj+1+GjΔQj+FjΔZj+1-FjΔZj=Φj（8）

式中：Cj、Dj、Ej、Gj、Fj、Φj均為離散系數(shù)。

式（7）與式（8）構(gòu)成的方程關(guān)系式是基于第j個長方形網(wǎng)格構(gòu)建的方程組，其未知數(shù)有ΔQj、ΔQj+1、ΔZj、ΔZj+1，僅在一個網(wǎng)格時方程式是不封閉的，基于整個研究河道時，將其分為L2-L1個單元河段，則存在L2-L1+1個斷面（圖2），可寫出2（L2-L1）個差分方程，分別給定上下邊界關(guān)系條件，從而構(gòu)成封閉的數(shù)學(xué)方程組，然后利用追趕法迭代計算。

1.3FORTRAN計算程序開發(fā)

運用FORTRAN95語言設(shè)計將追趕法迭代求解程序化，實現(xiàn)河道計算斷面水位與流量批量計算，模擬研究河道滲漏對洪水演進(jìn)影響效應(yīng)。程序框圖見圖3。

2模型實例驗證與分析

2.1大沽河流域概況

大沽河起源地為煙臺市轄招遠(yuǎn)阜山鎮(zhèn)，橫穿青

島市轄萊西、平度、即墨、膠州等區(qū)市，青島市境內(nèi)流域距離為1799 km，河流總面積達(dá)到6 1313 km2，主要有小沽河、豬河、落藥河、流浩河、城子河、桃源河、云溪河等多個支流匯入[24]。大沽河流域古峴以北由花崗巖與變質(zhì)巖組成的構(gòu)造剝蝕低山丘陵地形，古峴以南由碎屑巖和火山巖組成的構(gòu)造剝蝕平原，沿大沽河中下游河床兩側(cè)為河谷沖積平原，自北向南呈不規(guī)則的帶狀分布。古峴至店埠以南河谷開闊，地形平坦，地形標(biāo)高4～40 m，微向南傾斜，坡降2‰～07‰，下游河谷寬度一般6 km。古河谷被掩埋于地下，現(xiàn)代河床上疊在古河谷堆積物之上，河床深度一般2 m。河谷兩側(cè)為剝蝕——堆積準(zhǔn)平原，基底巖石主要為白堊系王氏組砂頁巖，上覆較薄的殘坡積層[25]。大沽河河道入滲過程大致分為非穩(wěn)定與穩(wěn)定入滲階段。非穩(wěn)定入滲階段，入滲率開始很大，然后快速下降，再緩慢下降；穩(wěn)定入滲階段，入滲率逐漸穩(wěn)定。河道包氣帶砂樣越粗糙，透水性越好，非穩(wěn)定入滲階段持續(xù)時間越短，河床包氣帶砂層為細(xì)顆粒時，滲透系數(shù)較小，入滲率較?。活w粒粗糙、松散時，滲透系數(shù)較大，入滲率較大[14]。

2.2耦合滲漏項的大沽河流域洪水演進(jìn)模型

（1）模型條件。

①基本參數(shù)。本文選擇研究區(qū)域為大沽河干流河段上起產(chǎn)芝水庫和北墅水庫，下至南莊橡膠壩，距離為1027 km。斷面劃分為103份，邊界糙率為0039，主槽糙率為0041，漫灘糙率為005，權(quán)重系數(shù)為08，時間步長為1 min。

②邊界條件。起點產(chǎn)芝水庫水文站流量過程確定為上邊界關(guān)系條件；終點南莊橡膠壩的水位過程確定為下邊界條件；大沽河流域各個支流當(dāng)作點源流量邊界條件進(jìn)行處理。

③初始條件。對大沽河研究河段進(jìn)行恒定非均勻流水面線的計算。

④入滲率。以季節(jié)性河流大沽河上、中、下游三處野外垂向入滲試驗數(shù)據(jù)為基礎(chǔ)，建立入滲率公式，構(gòu)建耦合滲漏項的大沽河洪水演進(jìn)模型。

（2）模型驗證與分析。

利用2003年4月大沽河大規(guī)模調(diào)水工程收集的水文資料數(shù)據(jù)對考慮滲漏項的洪水演進(jìn)模型進(jìn)行驗證。選取沙灣莊和南村作為代表性斷面，分別對水位過程、流量過程、過水時間、滲漏量進(jìn)行了模擬驗證。模擬結(jié)果見圖4-圖7與表1。

圖4、圖5分別為沙灣莊監(jiān)測斷面水位和流量過程的模擬值和實測值對比，圖6、圖7分別為南村監(jiān)測斷面水位和流量過程的模擬值與實測值對比。可以得出，沙灣莊與南村兩斷面的模擬值與實測值變化趨勢基本一致，洪水峰面達(dá)到斷面模擬時間與實際時間基本一致。

由表1可知，產(chǎn)芝水庫至賈疃閘區(qū)間各監(jiān)測斷面來水量實測值與模擬值相對誤差范圍為-166%～036 %，各監(jiān)測斷面滲漏量的實測值和模擬值相對誤差范圍為-309%～372 %，表明了模擬值和實測值基本吻合。

3結(jié)論

季節(jié)性河道汛前期常年處于干枯和斷流狀態(tài)，河床基質(zhì)滲透能力較強，汛期洪水呈現(xiàn)水位高、流量大、暴漲暴落等特征，且沿程河道水量損失比較大，大部分為河道滲漏損失。使用傳統(tǒng)的不考慮滲漏洪水演進(jìn)數(shù)學(xué)模型進(jìn)行洪水模擬預(yù)報，結(jié)果精度受到較大影響。本文建立了耦合滲漏項的季節(jié)性河道洪水演進(jìn)模型，解決了河道滲漏損失演算、實際過程入滲率參數(shù)確立、滲漏與洪水演進(jìn)模型耦合等問題。

利用青島市大沽河流域?qū)崪y水文資料對所建模型進(jìn)行應(yīng)用與驗證，沙灣莊與南村斷面模擬斷面水位和流量過程的模擬值與實測值基本吻合且曲線趨勢一致，洪水峰面到達(dá)斷面模擬時間和實際時間基本一致；研究區(qū)各個河段來水量、滲漏量的實測值和模擬值相對誤差范圍分別為-166%～036 %、-309%～372 %，表明所建模型精度較高，可滿足在季節(jié)性河道現(xiàn)狀條件下進(jìn)行洪水演進(jìn)模擬和預(yù)報的需要。

參考文獻(xiàn)（References）：

[1]范朝暉.論北方城市河道多樣形態(tài)斷面可行性[J].新農(nóng)村：黑龍江，2014，14（8）：288-289.（FAN Z H.On the feasibility of diversified river channel forms in northern cities[J].New Countryside：Heilongjiang，2014，14 （8）：288289.（in Chinese）） DOI：10.3969/j.issn.16748409.2014.08.272.

[2]郝鵬.季節(jié)性河流河道滲漏試驗與模擬研究[D].青島：中國海洋大學(xué)，2008.（ HAO P.Experimental and simulation study on channel seepage in seasonal rivers[D].Qingdao：Ocean University of China，2008.（in Chinese）） DOI：10.7666/d.y1337881.

[3]劉江俠，王佰偉.考慮地下水入滲的洪水演進(jìn)模型構(gòu)建及應(yīng)用[J].水資源與水工程學(xué)報，2017，28（4）：168172.（LIU J X，WANG B W.Study on flood routing model and application considering groundwater leakage[J].Journal of Water Resources and Water Engineering，2017，28（4）：168172.（in Chinese）） DOI：10.11705/j.issn.1672643X.2017.04.28.

[4]DAHAN O，TATARSKY B，ENZEL Y，et al.Dynamics of flood water infiltration and ground water recharge in hyperarid desert[J].Ground Water，2008，46（3）：450461.

[5]SHENTSIS I，ROSENTHAL E.Recharge of aquifers by flood events in an arid region[J].Hydrological Processes，2003，17（4）：695712.

[6]田守波，鄭西來，林國慶，等.大沽河河道滲漏對洪水演進(jìn)的影響研究[J].水資源保護(hù)，2009（S1）：1619.（TIAN S B，ZHENG X L.LIN G Q.Study on the effect of channel leakage on flood routing of Dagu River[J].Water Resources Protection，2009（S1）：1619.（in Chinese））

[7]齊春英，劉克巖.沿程滲漏河道的洪水流量演算模型[J].水文，1997（6）：2730.（QI C Y，LIU K Y.A flood flow calculation model for seepage channel along channel[J].Hydrology，1997 （6）：2730.（in Chinese））

[8]田守波.大沽河干流滲漏對河道洪水演進(jìn)的影響研究[D].青島：中國海洋大學(xué)，2009.（TIAN S B.Study on the effect of channel leakage on flood routing of Dagu River[D].Qingdao： Ocean University of China，2010.（in Chinese））

[9]田月，田淳，郝瑞霞.河道水力學(xué)模型及最小生態(tài)需水量的估算[J].太原理工大學(xué)學(xué)報，2016，47（4）：536540.（TIAN Y，TIAN C，HAO R X.Estimation of river hydraulics model and minimum ecological water demand[J].Journal of Taiyuan University of Technology，2016，47 （4）：536540.（in Chinese）） DOI：10.16355/j.cnki.issn10079432tyut.2016.04.020.

[10]程亮，王宗志，胡四一，等.強烈下滲條件下天然河道洪水演進(jìn)模擬方法[J].中國科學(xué)：地球科學(xué)，2015（2）：207215.（CHENG L，WANG Z Z，HU S Y，et al.Simulation method of natural river flood routing under severe infiltration[J].China Science：Earth Science，2015 （2）：207215.（in Chinese）） DOI：10.1007/s114300145018x.

[11]王宗志，程亮，王銀堂，等.高強度人類活動作用下考慮河道下滲的河網(wǎng)洪水模擬[J].水利學(xué)報，2015，46（4）：414424.（WANG Z Z，CHENG L，WANG Y T，et al.Flood simulation of river network under the influence of high intensity human activity.[J].Journal of Hydraulic Engineering，2015，46 （4）：414424.（in Chinese）） DOI：10.13243/j.cnki.slxb.2015.04.005.

[12]MUDD S M.Investigation of the hydrodynamics of flash floods in ephemeral channels：scaling analysis and simulation using a shockcapturing flow model incorporating the effects of transmission losses[J].Journal of Hydrology，2006，324（1）：6579.

[13]MORIN E，GRODEK T，DAHAN O，et al.Flood routing and alluvial aquifer recharge along the ephemeral arid Kuiseb river，Namibia[J].Journal of Hydrology，2009，368（1）：262275.

[14]武桂芝.大沽河河道滲漏及其對洪水演進(jìn)影響研究[D].青島：中國海洋大學(xué)，2010.（WU G Z.River seepage in Dagu river and its influence on flood routing[D].Qingdao： Ocean University of China，2010.DOI：10.7666/d.y1830553.（in Chinese））

[15]武桂芝，鄭西來，林國慶，等.季節(jié)性河道滲漏野外試驗研究[J].干旱區(qū)研究，2011，28（1）：7478.（ WU G Z，ZHENG X L，LIN G Q，et al.Field experimental study on seasonal river seepage.[J].Arid Zone Research，2011，28 （1）：7478.（in Chinese）） DOI：10.13866/j.azr.2011.01.020.

[16]馮增帥.季節(jié)性河道非飽和帶現(xiàn)場入滲試驗及室內(nèi)試驗數(shù)值模擬研究[D].青島：青島理工大學(xué)，2016.（FENG Z H.Numerical simulation of field infiltration test and laboratory test on unsaturated channel of seasonal channel[D].Qingdao：Qingdao University of Technology，2016.（in Chinese））