基于HEC-RAS的梧桐山河流域泥沙輸移特性分析

舒遠(yuǎn)麗 胡婷婷 王 拓

(三峽大學(xué)水利與環(huán)境學(xué)院，湖北 宜昌 443000)

天然河流在徑流形成和運(yùn)行過程中，常挾帶不同程度的泥沙顆粒，造成河床淤積、河道改變，給防洪、灌溉及水利建設(shè)帶來不同程度的影響。河流泥沙不僅反映水土流失狀況，還是流域內(nèi)降水、徑流、土壤、地形地貌植被等流域特性的體現(xiàn)[1～4]，徑流泥沙是衡量土壤侵蝕的重要參數(shù)之一，可為土壤侵蝕動力過程的模擬與研究、土壤侵蝕預(yù)報模型的建立等提供基礎(chǔ)資料，為水土流失的監(jiān)測、防治等提供科學(xué)依據(jù)。數(shù)學(xué)模型用于對泥沙運(yùn)動理論問題深入研究，河流泥沙數(shù)學(xué)模型能對水流泥沙的運(yùn)動過程進(jìn)行模擬，模擬泥沙組成的各種變化及沖淤量與河床形態(tài)[5]。與物理模型比較，迅速發(fā)展的河流泥沙數(shù)值模型，明顯看出其周期短、成本低、應(yīng)用靈活的特點(diǎn)。從國內(nèi)外泥沙數(shù)值模型研究成果來看，在工程上河流泥沙數(shù)值模型已經(jīng)被廣泛地應(yīng)用[6]。近年，國內(nèi)外出現(xiàn)的較為成熟的數(shù)學(xué)模型，多用于預(yù)測解決實(shí)際工程問題。如黃永健等[7]按汛期與非汛期引水引沙條件，建立泥沙數(shù)學(xué)模型，并利用實(shí)測資料對模型進(jìn)行了驗(yàn)證；于鵬[8]應(yīng)用二維水流數(shù)學(xué)模型進(jìn)行水面線推求，應(yīng)用于東汶河工程治理當(dāng)中；劉英等[9]結(jié)合渠道比降、流量等不同的因素對U形渠道水力性能的影響進(jìn)行研究；王景元等[10]利用支渠特性，增強(qiáng)挾沙能力，減少干渠泥沙淤積量；季飛等[11]建立二維數(shù)值模型，考慮邊界條件，計算分析洪奇瀝航道整治工程后工程區(qū)的沖淤變化。MAREN.V et al.[12]考慮鹽度、波浪、潮汐等要素建立泥沙數(shù)值模型，研究渠道疏浚措施對懸沙濃度及淤積的影響。

以上研究雖從不同方面分析了河渠工程的水力特性及泥沙淤積問題，但較少從全流域角度考慮河道泥沙輸移特性及流域出口斷面輸移總量，本文基于HEC-RAS軟件研究暴雨條件下梧桐山流域河道泥沙輸移特性及模擬估算流域出口斷面年內(nèi)泥沙輸移總量。

1 流域介紹

深圳市位于廣東省東南部珠江口的東岸，北連惠州市、東莞市，南隔深圳河與香港九龍新界相鄰，東依大鵬灣、大亞灣，西瀕伶仃洋與珠海市相望。深圳依山臨海，境內(nèi)流域面積大于1km2的河流共有310條，流域面積大于100km2的河流有深圳河、茅洲河、龍崗河、觀瀾河和坪山河，分屬9大流域。梧桐山河發(fā)源于梧桐山北麓，流經(jīng)橫瀝口、坑背等7個村注入深圳水庫，沿線4條支流匯入梧桐山河，全流域?qū)贃|深供水水系，為飲用水源保護(hù)區(qū)(見圖1)。河道上游建有橫瀝口水庫，屬小(2)型。河道長度為3.87km，流域面積12.53km2，管理范圍面積14.34萬m2，河道防洪標(biāo)準(zhǔn)為50年一遇。

圖1 深圳河灣水系

梧桐山流域?qū)儆趤啛釒ШＱ蠹撅L(fēng)氣候。夏季受海陸氣溫差異影響，冬季則受到西伯利亞冷空氣影響，全年雨水充沛。由于地處珠江口地區(qū)，臺風(fēng)數(shù)量多，登陸頻繁，持續(xù)時間長。根據(jù)1953—1999年的統(tǒng)計資料，登陸珠江口的臺風(fēng)達(dá)到每年4.1次。臺風(fēng)通常在5—11月出現(xiàn)，并且在7—9月更為頻繁。該流域?qū)儆谇鹆甑貐^(qū)，整體地勢東高西低、北高南低。土壤地質(zhì)狀況為：東北部為灰黃土、紅土等壤土和黏土，西南部為第四紀(jì)殘破積土、沖洪水土和海積黏性土。

2 理論模型

2.1 模型

HEC-RAS (River Analysis System)模型主要用于天然或人造河網(wǎng)的一維水力學(xué)計算。HEC-RAS主要由恒定流水面線計算、非恒定流模擬、水質(zhì)分析、可動邊界泥沙輸移計算4部分組成。此次計算只涉及水動力及泥沙模塊。運(yùn)行采用水動力模塊與泥沙輸移計算串聯(lián)模式進(jìn)行，即泥沙模型的計算結(jié)果可以通過河床變形、邊界糙率修改等途徑反作用于水動力模塊。該模式計算涉及河床變形、糙率變化等，根據(jù)選用的相應(yīng)計算公式，模擬動態(tài)的河床沖淤過程。

2.1.1 水動力模塊

HEC-RAS 的非恒定流模擬基于連續(xù)方程和動量方程，其中連續(xù)方程為

(1)

式中：ρ為流體密度,kg/m3；u為流速,m/s；t為時間,s。

動量方程為

(2)

式中：f為質(zhì)量力,N；P為壓力,N；ν為流體運(yùn)動黏滯系數(shù)；xi，xj為離下一斷面的距離,m。

2.1.2 泥沙輸移模塊

泥沙輸送方程：Meyer-Peter公式在其試驗(yàn)資料范圍內(nèi)有很高的精度，Meyer-Peter公式為

(3)

式中：Qb=BRbU，Q=BhU；Kb和K′b是阻力系數(shù)；B為河道寬度，m；h為水深,m；J為比降；U為平均流速，m/s；γs為泥沙比重；γ為水的比重；g為重力加速度,9.8m/s2；D為泥沙粒徑,mm；gb為推移質(zhì)單寬輸沙率,kg/(s·m)；Rb為與河床床面有關(guān)的水力半徑，m。

泥沙沉積的計算是由Krone(1962)公式來計算沉積速率量化，具體公式為

(4)

式中：C為沉積物濃度,kg/m；t為時間,s；τb為河床切應(yīng)力,N/m2；τc為沉積臨界切應(yīng)力,N/m2，Vs為下降速度,m/s；y為水深,m。

通過分離變量并進(jìn)行積分，可以建立以下關(guān)系，公式為

(5)

泥沙侵蝕的計算為

(6)

式中：m為沉積物的質(zhì)量,kg；t為時間，s；τb為河床切應(yīng)力，N/m2；τc為侵蝕臨界切應(yīng)力,N/m2；M為沖刷經(jīng)驗(yàn)侵蝕率,kg/s。

通過分離變量并進(jìn)行積分，可以建立以下關(guān)系：

(7)

2.2 模型建立及驗(yàn)證

2.2.1 模型范圍與剖分

考慮梧桐山河河段的地形特點(diǎn)、河流洪水特點(diǎn)等，選取模型河段范圍為:上游邊界進(jìn)口斷面(樁號3800)離橫瀝口水庫約80m，上游邊界距下游邊界出口斷面(樁號80)3.8km。為了較好反映河道地形，滿足流場計算精度要求，根據(jù)研究的工程問題與河道平面的特點(diǎn)，實(shí)測地形，橫斷面采用實(shí)際測量數(shù)據(jù)，通過斷面編輯器逐一輸入每個橫斷面的坐標(biāo)、距下游相鄰斷面的距離、主槽及邊灘的糙率等參數(shù)，即可完成河道橫斷面的輸入。梧桐山河段共布設(shè)57個斷面，其中斷面間最大間距為150m，最小間距為50m。所繪制河段橫斷面及河系顯示在幾何數(shù)據(jù)編輯器的繪圖區(qū)中，用于模擬河段。建好的河系見圖2。

圖2 河道斷面模擬圖

2.2.2 模型參數(shù)

a.時間步長：模型計算時間步長1h，輸出時間步長1h。

b.糙率是影響水面線計算極為敏感的因素之一，在有歷史洪水調(diào)查資料時，一般通過歷史洪水水面線進(jìn)行率定。缺乏資料時，則根據(jù)河床、河灘情況，對照天然河道糙率表確定。梧桐山河無歷史洪水資料，根據(jù)現(xiàn)場勘探及資料確定河床及河灘糙率，并利用實(shí)測流量水位進(jìn)行驗(yàn)證。

c.HEC-RAS非恒定流模擬提供了1種上游邊界條件和3種下游邊界條件，分別為上游流量過程曲線、正常水深、階段序列、下游水位-流量關(guān)系曲線。本文梧桐山河河段上游邊界采用實(shí)測流量曲線及斷面地形，下游邊界則利用正常水深。

d.泥沙粒徑：根據(jù)當(dāng)?shù)赝寥罈l件，選擇的泥沙粒徑見表1。

表1 泥沙粒徑

2.2.3 模型驗(yàn)證

為了使模型能夠正確模擬計算區(qū)域的河道洪水演進(jìn)及反映河道泥沙輸移特征，根據(jù)2019年梧桐山河實(shí)測來水過程及水位觀測資料對模型進(jìn)行參數(shù)率定。模型在率定過程中，梧桐山河段分別采用不同糙率模擬流場阻力。采用RTK采集的梧桐山河道9個斷面水位數(shù)據(jù)整理得到河道沿程水面線進(jìn)行參數(shù)率定，沿程水位高程見表2。

表2 梧桐山河沿程水面線 單位：m

經(jīng)模型率定，模型的曼寧系數(shù)一般取0.025～0.060。由圖3對比河段各斷面的實(shí)測水位與模擬的計算水位可知：計算值與實(shí)測值擬合較好，表明數(shù)學(xué)模型對河道阻力的模擬是適宜的。

圖3 水面線對比

2.2.4 含沙量觀測

水體含沙量的測量是水文觀測中的重要內(nèi)容，輸沙量大小隨著河道含沙量的變化而變化。泥沙含量是反映水中固體顆粒物泥沙多少的物理量。本次研究在梧桐山河道開展了多次取樣工作，獲得了含沙量數(shù)據(jù)。含沙量的測量主要通過烘干法又稱稱重法。使用取樣瓶采取水樣帶回梧桐山項(xiàng)目部實(shí)驗(yàn)室，采取水樣后，稱量水樣原始的質(zhì)量與烘干后的質(zhì)量，由此確定水樣的含沙量。首先分3點(diǎn)測量分流桶內(nèi)水深，然后將分流桶內(nèi)泥沙攪勻。攪勻后，邊攪動邊用取樣瓶進(jìn)行渾水取樣 (約3000mL)，取回的渾水樣在實(shí)驗(yàn)室內(nèi)用過濾法測定含沙量。即取500mL的渾水樣用濾紙過濾，然后將附有泥沙的濾紙置于烘箱內(nèi)在105℃恒溫條件下烘24h后,測量濾紙和泥沙的重量，減去濾紙的烘干重量即為泥沙干重，除以水樣體積(500mL)，則得到水樣含沙量。再做一個重復(fù)渾水樣測量。取二者含沙量的平均值為該泥沙重量。表3為2019年7月8—18日期間河道含沙量數(shù)據(jù)，數(shù)據(jù)表明河道含沙量在暴雨期和無雨期差異很大。

表3 河道含沙量數(shù)據(jù)

3 模擬與分析

利用2019年8月1—2日的20190801場次暴雨資料，采用HEC-RAS對暴雨情況下梧桐山河流域河道泥沙輸移進(jìn)行模擬。

3.1 河道泥沙淤積分析

選擇河道樁號2050、樁號820、樁號540、樁號80共4個典型斷面進(jìn)行分析。利用HEC-RAS模擬計算場次暴雨后河道沿程累計淤積量及出口斷面泥沙輸移量。暴雨后沿程河床典型斷面淤積情況對比見圖4。

圖4 河床典型斷面淤積情況對比

河道上游進(jìn)口斷面至下游出口斷面河段總的泥沙淤積量為68.83t。其中，樁號2050斷面河床底部高程值由48.62m變?yōu)榱?8.76m，升高了0.14m，中游彎道樁號2100～2050斷面間河段泥沙淤積量達(dá)到13.59t，樁號2100～2050斷面間河段處于彎道處，流速降低，造成泥沙堆積。

樁號820斷面河床底部高程值由36.26m變?yōu)榱?6.48m，升高了0.22m，在樁號900～820斷面河段間泥沙淤積量達(dá)到最大值22.13t，樁號900～80斷面間河段屬于整個河段淤積強(qiáng)度最大的河段。樁號540斷面河床底部高程值由33.05m變?yōu)榱?3.16m，升高了0.11m，樁號580～820斷面河段間泥沙淤積量為10.51t。樁號900～820、580～540斷面間河段處于轉(zhuǎn)彎處，彎道形成緩流區(qū)，而泥沙由于彎道環(huán)流的作用，導(dǎo)致泥沙淤積。

樁號80出口斷面處河床底部高程值由28.83m變?yōu)榱?8.95m，升高了0.12m，樁號200～80斷面間河段泥沙淤積量為4.96t。由于出口斷面下游為深圳水庫，過水?dāng)嗝孀儗挘掠嗡畮祉斖凶饔脤?dǎo)致水流流速變緩，形成泥沙淤積。

樁號80斷面為梧桐山河河道出口斷面，下連深圳水庫，模型模擬河道出口斷面在20190801場次暴雨的泥沙輸出量為12.70t。

表4為梧桐山河河段沖淤總量實(shí)測值與計算值統(tǒng)計，實(shí)測沖淤總量與計算值有一定的誤差，實(shí)測值總體較計算值偏大。模型誤差最大發(fā)生在樁號2100～2050斷面間河段，實(shí)測值約為模型計算值的90.95%。在樁號2100～2050、樁號580～540、樁號200～80斷面之間，其模擬結(jié)果分別為實(shí)測值的93.20%、95.03%、94.49%。

表4 分段泥沙淤積量

3.2 年內(nèi)輸沙總量估算

選取平日及暴雨前后含沙量代表值，由表3可知含沙量范圍在4.00～1100.00g/m3，無雨時含沙量在3.92～4.95g/m3，與暴雨開始到結(jié)束泥沙含沙量96.20～1095.00g/m3相差24～280倍，假定無雨?duì)顟B(tài)下，梧桐山河出口斷面泥沙輸移總量近似為零。通過模型模擬計算所得場次典型暴雨后的泥沙輸移量用于估計年內(nèi)泥沙輸移總量。由模型模擬結(jié)果知20190801場次暴雨流域出口斷面(樁號80)向深圳水庫輸沙量為12.70t。由全年雨量監(jiān)測站資料可知，梧桐山河流域2019年有15場暴雨，據(jù)此估算流域出口斷面(樁號80)年內(nèi)輸移總量約190.50t。

參考《廣東水資源》(1986年8月)成果以及梧桐山河流域上游集雨面積、地形和地貌相近的測站資料統(tǒng)計，確定流域內(nèi)多年平均含沙量為0.09kg/m3，多年平均年徑流量220.75萬m3，由此推得梧桐山河下游出口斷面多年輸沙量為198.68t。

由模型模擬出流域出口斷面年內(nèi)輸沙量190.50t，與深圳市該地區(qū)所查多年平均含沙量所計算出的年內(nèi)輸沙量198.68t結(jié)果相近。

4 結(jié) 語

本文以梧桐山河河段為例，利用HEC-RAS泥沙計算模塊對研究河段進(jìn)行了數(shù)值模擬，初步探索了暴雨條件下河道內(nèi)泥沙輸移規(guī)律特性，并利用模型模擬計算結(jié)果估算流域出口斷面年內(nèi)泥沙輸移量，得到了以下結(jié)論:

a.模型模擬結(jié)果表明，河道的淤積部位主要位于河道中下游及河道彎道處，淤積最大值出現(xiàn)在臨近中下游彎道河段。

b.由模型模擬結(jié)果可知單場暴雨河道總的淤積量為68.83t，流域出口斷面泥沙淤積量為12.70t，假定無雨期梧桐山流域泥沙輸移量近似為0，則梧桐山河流域全年15場暴雨年內(nèi)輸移總量為190.50t。

模擬結(jié)果與實(shí)際對比表明HEC-RAS模型能夠較好地模擬梧桐山河流域暴雨后泥沙沖淤過程，這也表示基于HEC-RAS模型的泥沙計算能適應(yīng)梧桐山河流域泥沙模擬。水沙輸運(yùn)的累積效應(yīng)也會使梧桐山河的地形地貌發(fā)生改變，對于梧桐山河流域泥沙輸移規(guī)律特性的研究可服務(wù)于梧桐山河流域河道及河口整治，為流域水土流失的監(jiān)測、防治等提供一定科學(xué)依據(jù)。