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

舒遠麗 胡婷婷 王 拓

(三峽大學水利與環境學院，湖北 宜昌 443000)

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

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

1 流域介紹

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

圖1 深圳河灣水系

梧桐山流域屬于亞熱帶海洋季風氣候。夏季受海陸氣溫差異影響，冬季則受到西伯利亞冷空氣影響，全年雨水充沛。由于地處珠江口地區，臺風數量多，登陸頻繁，持續時間長。根據1953—1999年的統計資料，登陸珠江口的臺風達到每年4.1次。臺風通常在5—11月出現，并且在7—9月更為頻繁。該流域屬于丘陵地區，整體地勢東高西低、北高南低。土壤地質狀況為：東北部為灰黃土、紅土等壤土和黏土，西南部為第四紀殘破積土、沖洪水土和海積黏性土。

2 理論模型

2.1 模型

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

2.1.1 水動力模塊

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

(1)

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

動量方程為

(2)

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

2.1.2 泥沙輸移模塊

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

(3)

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

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

(4)

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

通過分離變量并進行積分，可以建立以下關系，公式為

(5)

泥沙侵蝕的計算為

(6)

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

通過分離變量并進行積分，可以建立以下關系：

(7)

2.2 模型建立及驗證

2.2.1 模型范圍與剖分

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

圖2 河道斷面模擬圖

2.2.2 模型參數

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

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

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

d.泥沙粒徑：根據當地土壤條件，選擇的泥沙粒徑見表1。

表1 泥沙粒徑

2.2.3 模型驗證

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

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

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

圖3 水面線對比

2.2.4 含沙量觀測

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

表3 河道含沙量數據

3 模擬與分析

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

3.1 河道泥沙淤積分析

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

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

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

樁號820斷面河床底部高程值由36.26m變為了36.48m，升高了0.22m，在樁號900～820斷面河段間泥沙淤積量達到最大值22.13t，樁號900～80斷面間河段屬于整個河段淤積強度最大的河段。樁號540斷面河床底部高程值由33.05m變為了33.16m，升高了0.11m，樁號580～820斷面河段間泥沙淤積量為10.51t。樁號900～820、580～540斷面間河段處于轉彎處，彎道形成緩流區，而泥沙由于彎道環流的作用，導致泥沙淤積。

樁號80出口斷面處河床底部高程值由28.83m變為了28.95m，升高了0.12m，樁號200～80斷面間河段泥沙淤積量為4.96t。由于出口斷面下游為深圳水庫，過水斷面變寬，下游水庫頂托作用導致水流流速變緩，形成泥沙淤積。

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

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

表4 分段泥沙淤積量

3.2 年內輸沙總量估算

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

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

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

4 結 語

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

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

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

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