航道工程作用下福姜沙水道灘槽演變特性及對水沙過程響應的數值研究

孫孟喆，胡 鵬, 2，李有為，劉奇峰, 曹志先

(1. 浙江大學 海洋學院，浙江 舟山 316021; 2. 浙江大學 舟山海洋研究中心，浙江 舟山 316021; 3.長江航道規劃設計研究院，湖北 武漢 430040; 4. 武漢大學 水利水電學院，湖北 武漢 430072)

福姜沙水道位于長江口感潮河段，上起江陰鵝鼻嘴，下迄護槽港，全長約30 km。水動力受上游徑流和下游潮流的共同作用，泥沙輸移和河床演變復雜，是航道整治的難點河段[1-4]。福姜沙水道在平面上呈現“兩級分汊、三汊并存”的形態：主流經福姜沙分為左汊和右汊(福南水道)，左汊水流經雙澗沙分為福北水道和福中水道，形成長江下游典型的多級分汊多級匯流的復雜河勢格局。福北水道和福中水道為主航道，分別為12.5 m深水航道的上行通航分道和上下行通道。歷史上雙澗沙的沖淤變化以及灘面橫流、竄溝發展對福中、福北水道的航行條件產生了不利影響[5-6]，為了限制雙澗沙的不穩定對各汊道航行條件的不利影響，2009年至2011年期間實施了雙澗沙守護工程，分為雙澗沙頭部潛堤、北順堤、南順堤3部分。工程雖然達到“封堵竄溝、固灘護沙”的整治目的，但由于工程規模有限，灘面流未完全遏制，兩側深槽發展仍不能有效控制[7-9]，甚至對航道產生淤積等不良影響[4]。為解決這一問題，將長江口深水航道進一步向上延伸，航道部門對福姜沙水道實施了深水航道二期工程，在福姜沙左緣建設4道丁壩，在雙澗沙北側建設4道丁壩、南側建設8道丁壩，優化雙澗沙頭部潛堤結構。2017年3月二期工程完工后，由于底沙下泄和彎道環流輸沙影響，仍需對航槽邊緣及航槽內淤淺部位進行維護以達到12.5 m深水航道通航尺度[10-11]。研究福姜沙水道的沖淤特征，并揭示其在航道工程作用下的灘槽演變特性顯得尤為重要和緊迫。

眾多學者基于實測資料分析、物理模型試驗、數值模擬的方法，從河流動力學和河流地貌學的角度對福姜沙水道的灘槽演變特性進行了研究[1, 3-4, 6-14]。研究表明，水流動力軸線擺動是河床調整變化的主要動力因子[12]，邊心灘變形、航道回淤與特征流量級的大小和持續時間有關[14]，河道整體沖淤與上下游河道演變密切相關[15]。已有研究成果通常基于定性的描述，且研究背景多集中于二期工程實施前，對航道工程實施后的河道演變特性研究較為匱乏。基于實測地形、水文資料分析表明，航道工程實施前后河道的演變特征發生了改變[16]，正在向與工程相適應的河勢調整，福姜沙左緣邊灘和福姜沙沙體體積增大，福中水道深槽擴容，福北水道深槽淤淺[11]。過去十余年，水沙數學模型被大量應用于福姜沙水道的整治方案設計、航槽回淤、灘漕演變規律[12]等問題的研究[4, 12, 17-21]。基于TVD—有限體積法的水沙模型發展已有二十余年的歷史[22-26]，但其在野外尺度的應用非常少見。這主要是因為：第一，TVD—有限體積法為顯式格式，計算時間步長受限于CFL穩定性條件；第二，為了編程的方便，以往數學模型習慣上采用滿足CFL條件的全局最小時間步長。當計算網格尺度或水流狀態不均勻時，整體時間步長被最小時間步長限制，大大降低模型的計算效率[27-28]。對于清水模型，局部時間步長方法(LTS)可以使每個單元的變量采用滿足穩定性條件盡可能大的時間步長進行更新，從而降低計算耗時[29-33]。進一步采用混合局部時間步長和整體最大時間步長方法(Hybrid LTS/GMaTS)將局部時間步長技術推廣至動床模型中，已應用于長江中游水道的沖淤數值計算[34]，但尚未應用于徑潮流共同作用的水道。

采用基于Hybrid LTS/GMaTS方法的平面二維水—沙—床耦合數學模型，以福姜沙水道徑潮流和河床沖淤過程為例，驗證其計算效率和精度。通過8組數值算例，分析二期工程實施后，福姜沙水道在不同流量特征下的水動力特征和河床沖淤特點，嘗試從洲灘特征長度變化的角度定量揭示洲灘演變對徑流特征的響應過程，在此基礎上建立航道工程影響下航道沖淤量與洲灘演變、徑流特征的數學關系式，以加深對福姜沙水道灘槽演變規律的認識，為航道工程提供參考。

1 水—沙—床耦合數學模型

1.1 控制方程

水—沙—床耦合數學模型的控制方程[35]為：

(1)

(2)

(3)

式中：t為時間；x和y為笛卡爾坐標系下的空間坐標；h為水深；u和v分別為x和y方向深度積分平均的流速分量；g為重力加速度；ck為第k粒徑組泥沙的深度平均體積含沙量；zb為床面高程；δ為床面活動層厚度；fa,k為活動層第k粒徑組泥沙體積百分比；Ek、Dk為床面與水流之間第k粒徑組泥沙交換的上揚和沉降通量，DT=∑Dk，ET=∑Ek；p0為泥沙孔隙率；Sbx和Sby分別為x、y方向的底坡坡度；Sfx和Sfy分別為x、y方向的阻力坡度；fs,k為底床存儲層和活動層交界面處第k粒徑組泥沙體積百分比；η=zb-δ為活動層下界面的高程。

1.2 經驗公式

阻力坡度采用曼寧公式計算：

(4)

(5)

其中，n為曼寧糙率系數。

底床與水流之間的泥沙上揚和沉降通量按式(6)～(7)計算：

Dk=αkckωk(1-αkck)m

(6)

Ek=αkce,kωk(1-αkce,k)m

(7)

式中：αk表示近底泥沙含量與深度平均含沙量之間的比值，考慮到該河段引起底床沖淤變化主要是臨底層的高濃度泥沙[13]，αk取3.0；ωk為第k粒徑組單顆粒泥沙沉速，采用張瑞瑾統一沉速公式計算；ce,k為第k粒徑組泥沙的水流挾沙力。采用張瑞瑾水流挾沙力計算公式計算[36]：

(8)

(9)

式中：υ為水流運動黏滯系數，取υ=1.2×10-6m2/s;dk為第k粒徑組泥沙的特征粒徑；s=(ρs-ρw)/ρw，ρs、ρw分別為清水和泥沙密度；K0、m0為經驗參數，本文取K0=0.12，m0=0.92。

1.3 基于Hybrid LTS/GMaTS的數值算法

采用非結構三角形網格如圖1所示。圖1(a)為計算區域某單元以及三個相鄰單元；圖1(b)為某個角點以及其所屬的三角形單元；圖1(c)為某個界面及其左右兩側的單元。單元變量使用Hybrid LTS/GMaTS進行更新，對控制方程采用有限體積法進行離散，挾沙水流的質量與動量守恒方程、隨流泥沙質量守恒方程[式(10)]采用局部時間步長(LTS)求解，底床所有組分泥沙總質量守恒方程[式(11)]、床沙分粒徑組的質量守恒方程[式(12)]采用整體最大時間步長(GMaTS)求解[34]：

圖1 非結構三角形網格Fig. 1 Sketches of the unstructured triangular meshes

(10)

(11)

(12)

2 模型設置與驗證

2.1 模型設置

計算區域如圖2所示。模型采用非結構化網格，對航道工程處進行局部加密，網格尺度最大為260 m，最小為56 m，平均為135 m。選取2017年2月實測地形作為初始地形，采用實測床沙級配數據。懸沙粒徑分為6個組別：0.042 mm、0.07 mm、0.085 mm、0.11 mm、0.16 mm和0.35 mm，各個粒徑的組分根據2017年8月福姜沙水道所有斷面懸移質實測數據加權平均得到。河道糙率采用經驗公式n=0.01+0.01/H計算[38]。上游邊界取在江陰河道內，采用大尺度模型為上游邊界提供流量過程，來沙采用大通站實測含沙量，下游邊界位于徐六涇水文站斷面，采用實測水位過程。

圖2 計算區域和初始地形Fig. 2 Location of the study area and initial bed topography

2.2 Hybrid LTS/GMaTS方法對計算結果的影響

使用均方根誤差(εRMSE)量化Hybrid LTS/GMaTS方法與傳統最小時間步長方法之間計算結果的相對誤差，使用全局泥沙質量的相對誤差ε量化泥沙計算的守恒性：

(13)

(14)

式中：f表示水位、流速、河床高程變化等物理量；N為網格節點數；fLG+GMa和fGMi分別表示Hybrid LTS/GMaTS方法和傳統最小時間步長方法的計算結果；V(t2)和V(t1)分別表示計算區域內t2和t1兩個時刻水體中泥沙的總體積；Vin和Vout分別表示兩個時刻之間從邊界流入或流出的泥沙體積；VE和VD分別表示兩個時刻之間河床與水體泥沙交換過程中上揚和沉降的泥沙體積。

模擬2017年8月15日至18日的水動力和河床沖淤過程，考慮muser=1～7的七個工況(muser是使用者給出的參數，用于限制LTS層級，muser越大，整體最大時間步長能達到的值越大)，為方便計算誤差，將最小時間步長設置為1 s。Hybrid LTS/GMaTS方法計算效率和相對誤差如表1所示。當muser=0時得到傳統整體時間步長方法，即每個網格的時間步長取全局最小值。Hybrid LTS/GMaTS方法可以大幅度提高計算效率，計算耗時隨著muser的增大而減小。Hybrid LTS/GMaTS方法(muser=7時)相比于傳統算法提升近十倍。傳統方法計算結果與實測數據的相對誤差為0.09 m(水位)、0.12 m/s(流速)、1.64 m(床面高程變化，2017年2月至2018年2月)，而Hybrid LTS/GMaTS方法相對于傳統方法的最大誤差僅12.93×10-4m(水位)、9.67×10-4m/s(流速)、7.01×10-4m(床面高程變化)，遠小于傳統方法計算結果與實測值之間的相對誤差。該方法在大幅度減少計算成本的基礎上帶來的計算誤差可忽略不計。Hybrid LTS/GMaTS方法(muser=7)的全局泥沙質量相對誤差量級僅為10-9，說明泥沙計算具有良好的守恒性。

表1 Hybrid LTS/GMaTS方法計算效率和相對誤差Tab. 1 Calculation cost and relative error using Hybrid LTS/GMaTS method

2.3 水動力及含沙量驗證

根據2017年洪季水文測驗結果，驗證福姜沙水道至通州沙水道沿岸四個潮位站(2L、4R、6R、任港)的水位(圖3)，四個測點(FZ-C、FZ1-B、JLG-A、TSG-B)的流速、流向和含沙量(圖4、圖5)。測站和測點位置如圖2所示。對比結果顯示，模型計算的水位、流速、流向、含沙量與實測過程線趨勢基本一致；水位計算結果與實測值的均方根誤差為0.07 m，計算水位與實測值的偏差均在0.1 m內；流速計算結果的均方根誤差為0.07 m/s，考慮到實測平均流速為0.75 m/s，平均流速偏差在10%以內；流向計算結果的均方根誤差為9°，平均流向偏差在10°以內；含沙量計算結果的均方根誤差為0.02 kg/m3，考慮到實測平均含沙量為0.10 kg·m3，平均含沙量偏差在30%以內。本模型的水動力和含沙量驗證滿足《水運工程模擬試驗技術規范》(JTS/T 231—2021)的規定。

圖3 計算水位與實測值對比Fig. 3 Comparison of calculated and measured water levels in different stations

圖4 計算垂向平均流速流向與實測值對比Fig. 4 Comparison of calculated and measured vertical average current speeds and directions

2.4 河床沖淤驗證

圖6分別為福姜沙水道2017年2月至2018年2月實測和計算沖淤分布。對比結果表明，模型成功復演了本河段2017年2月至2018年2月的灘漕沖淤過程的主要特征：雙澗沙整體淤漲，兩側深槽沖刷，靖江邊灘沖刷后退，計算沖淤部位和沖淤趨勢與實測結果相似。該河段實測沖淤量分別為4 953萬m3和4 391萬m3，計算沖淤量分別為4 296萬m3和4 152萬m3，沖淤量的偏差分別為13%和5%，沖淤量誤差均控制在20%內，滿足《海岸及河口潮流泥沙模擬技術規程》(JTS/T 231-2—2010)，說明本模型的泥沙模塊滿足精度要求。

圖6 2017年2月至2018年2月地形沖淤變化Fig. 6 Topographic changes from February 2017 to February 2018

3 灘槽演變特征對水沙過程的響應

3.1 工況設置

長江徑流流量存在洪枯季，且每時都發生變化，河道內灘槽演變過程是對不斷變化的徑流流量的綜合響應。徑流過程對于灘槽演變的綜合作用大致相當于某一特征流量級[39]。如表2所示，選取大通站2003年后多年流量累計頻率為5%(59 800 m3/s)、25%(40 500 m3/s)、50%(25 000 m3/s)、75%(15 000 m3/s)作為洪、中、枯水期的特征流量，選取大通站2003年后多年洪、中、枯水期平均含沙量(分別為0.24 kg/m3、0.14 kg/m3、0.09 kg/m3)作為來流特征含沙量。下游邊界條件考慮徑流和潮汐的影響，將水位過程分解為日平均水位和潮波水位[40]：特征流量與基準流量之差作為徑流增量對日平均水位進行預報[41]，徐六涇站實測潮位序列做調和分析對潮波水位進行預報。

表2 工況設置Tab. 2 Simulation cases

選取福姜沙水道2017年2月實測地形作為模型計算初始地形。一方面，此時二期工程已接近竣工，考慮航道工程對水道灘槽演變的影響較為合理。另一方面，水流切割下來的灘體淤積在六助港附近，福北水道航道條件有惡化趨勢[42]，歷史上靖江邊灘多次發生切割過程，切割下來的沙體呈現形成、沖刷下移、并入福北水道及雙澗沙、再形成的周期性過程[15]，該地形下河道演變過程對航道條件的影響具有代表性。每組工況中設置兩個算例，在第一個算例中考慮航道工程措施，包括整治建筑物和基建性疏浚工程，如圖7(a)所示。疏浚區域見圖7(a)，設計挖槽深度為理論最低潮面以下12.5 m，超深取0.5 m；第二個不考慮航道工程(即自然狀態)作為對比算例。8組工況的模擬時長均為90 d。

3.2 不同流量級條件下水動力及沖淤特征

流速大小反映了水流動力作用的強弱，河流斷面或沿程流速分布不均，導致水流挾沙力差異，使河床產生沖淤變化，從而影響河流的灘槽演變特性。水流經江陰水道流入福姜沙水道時，干流受到左岸炮臺圩節點和右岸鵝鼻嘴節點共同控制[42]。該河段內落潮流是塑造主槽河床的主要動力[4]，圖8是不同流量特征下各斷面大潮落潮流流速分布，斷面位置見圖7(a)。枯水期(流量為15 000 m3/s)大潮落潮流主流位于左汊深槽處，左汊深槽處水流流速大于左岸淺灘，靖江邊灘位于緩流區而淤積。隨著流量的升高，水流逐漸脫離炮臺圩節點擺動至右岸，對鵝鼻嘴節點的頂沖作用增大，節點挑流作用增強，水動力軸線逐漸向左岸邊灘擺動。當流量大于40 500 m3/s時，左岸淺灘處水流流速大于左汊深槽，大潮落潮流主流偏向靖江邊灘處，直接作用在靖江邊灘上，使靖江邊灘沖淤加劇，如圖7(c)、(d)所示。因此，該河段進口干流具有“大水趨直，小流坐彎”[4]的水流特性：枯水時主流在邊灘束水作用下走深槽，中洪水時主流線偏靠邊灘而取直，特別是在大流量(59 800 m3/s)下，邊灘呈現頭部沖刷，中部淤積，尾部沖刷的沖淤特點，說明邊灘在主流直接作用下出現沖刷后退的現象，有進一步被切割成心灘的風險。以上水流動力軸線在不同流量級條件下的擺動造成靖江邊灘洪沖枯淤的演變特性。而各流量級下航道工程對于福姜沙入口處水動力影響較小，在福姜沙左緣邊灘丁壩群作用下，深槽及左岸淺灘處流速略有增加，最大增加約0.1 m/s，右岸淺灘處流速有所減小，減小幅度約為0.05 m/s。

圖7 各工況下地形沖淤計算結果Fig. 7 Calculated results of topographic erosion and deposition

圖8 各斷面大潮落潮流流速分布Fig. 8 Cross-section ebb velocity distribution of spring tide

主流經雙澗沙分為福北、福中水道兩汊，主流進入福中水道，分流比約為70%[17]，如圖7(b)～(d)所示，在自然狀態下，流量為15 000 m3/s時，福中水道主流位于靠近福姜沙左緣邊灘的深槽內，福北水道主流于入口處位于偏左岸的深槽內，而后偏向雙澗沙一側。隨著流量增大，入口處主流擺至靠近雙澗沙沙頭，當流量超過40 500 m3/s時，主流進一步擺動至雙澗沙上部，而福中水道水流動力軸線擺向福姜沙北緣淺灘、福北水道水流動力軸線則擺動至雙澗沙左緣。這種主流擺動特征與兩汊河道呈現洪季“凸岸沖刷、凹岸淤積”、枯季“凹岸沖刷、凸岸淤積”[42]的演變特征一致，進而導致雙澗沙沙頭“洪沖枯淤”、雙澗沙兩側淺灘在不同流量級下沖淤交替的演變特征，帶來雙澗沙灘體左右擺動的不利影響。在福姜沙左緣丁壩和雙澗沙丁壩群的共同作用下，即使流量超過40 500 m3/s，福中水道主流也控制在主槽內，僅在丁壩挑流作用下向右岸淺灘一側輕微擺動。各流量級下福北水道入口處水動力條件均得到顯著改善，航槽刷深顯著。流量為15 000 m3/s時，福中水道中段總體呈現淤積趨勢，特別是福中、福南水道交匯處受支汊入流頂托效應影響，淤淺尤為顯著。然而隨著上游來流流量增大，深槽處水動力增強，福中水道總體上轉淤為沖，呈現出“洪沖枯淤”的特點。福北水道入口處流速相比于自然狀態下反而有所降低，受到雙澗沙北緣丁壩群挑流作用，福北水道主流穩定在深槽內，隨著流量增大，水流動力軸線偏向左岸，是靠近雙澗沙一側深槽轉沖為淤[圖7(a)～(d)]的主要原因。上游來流流量增大、入口處水動力作用顯著增強的情況下，淤積反而增加，說明水動力作用不是影響福北水道沖淤變化的主要因素。作用于雙澗沙灘體上的水動力條件減弱，說明灘面流在航道工程作用下得到控制，且位于丁壩壩田區的灘體在各流量級下均呈現淤積趨勢，特別是洪季雙澗沙灘頭沖刷趨勢得到抑制，雙澗沙灘體形態趨于穩定，對于兩汊的發展產生有利影響。

3.3 不同流量級條件下灘槽演變特征

3.3.1 灘槽演變特征

將12.5 m等深線(基于理論最低潮面)作為劃分灘槽的依據，繪制各流量級下灘槽分布圖，如圖9所示。流量為15 000 m3/s時，在自然狀態下躍進港附近邊灘沖刷萎縮，而蟛蜞港附近邊灘淤積展寬，說明在水流作用下，靖江邊灘整體上呈現沖刷后退的趨勢，并隨著流量的增大，沖刷趨勢增強。當流量達到59 800 m3/s時，靖江邊灘尾部深入主流的低矮灘體被水流切割脫離邊灘，進入航槽形成心灘，符合靖江邊灘在高于60 000 m3/s流量下完成切割過程的歷史演變特征[13]。中枯水時期，雙澗沙沙舌形態幾乎不發生變化，但當流量大于40 500 m3/s，雙澗沙沙舌受水流沖擊后退。福姜沙左緣邊灘和雙澗沙右緣淺灘在各流量級下灘面形態幾乎沒有發生大的變化，僅在福姜沙尾部有向福南水道左岸淺灘淤積展寬的趨勢，因此福中水道航槽條件無明顯影響。初始地形下福北水道夏仕港附近12.5 m深槽有中斷，隨著流量增大，中斷長度有所減小。但福北水道深槽在洪季(流量大于40 500 m3/s)發生較大程度的變形，左岸淺灘向內側擠壓航槽，而雙澗沙左緣淺灘則受沖刷后退，福北水道航槽整體有向南擺動的趨勢，對于航槽穩定性及左岸港口水深產生了不利影響。

圖9 各工況下灘槽格局變化Fig. 9 Change in morphology of bars and channels

在航道工程作用下，靖江邊灘尾部和福北水道兩側淺灘的灘槽演變趨勢產生了根本性變化。在一次性疏浚條件下，福北水道12.5 m航槽被打通。夏仕港至焦港航槽在各流量級下均保持穩定，且在雙澗沙左緣丁壩的控制下，12.5 m深槽有向左岸展寬的趨勢。焦港至如皋港段水流脫離雙澗沙航道工程控制，雙澗沙一側淺灘向內擠壓深槽，隨著流量增大，深槽內回淤程度加劇。雙澗沙沙舌向福北水道航槽內淤漲，當流量達到40 500 m3/s時，雙澗沙沙舌與淤積在六助港附近的淺灘相連，依靠疏浚打通的12.5 m深水航槽被阻斷。當流量達到59 800 m3/s時靖江邊灘尾部深入主流的低矮灘體沒有被水流切割，而是繼續向航槽內發育，在福北水道入口處形成淺灘，進一步影響福北水道的航行條件。說明航道工程建成后，靖江邊灘的沖刷下移和雙澗沙沙舌的橫向發育是影響福北水道入口處航槽穩定的主要因素。

3.3.2 洲灘演變及特征流量對航道條件的影響

將靖江邊灘垂直流向最寬處至炮臺圩節點的直線距離作為靖江邊灘特征尺度，用于表征靖江邊灘沖刷后退的程度，繪制靖江邊灘特征尺度關于特征流量的散點圖，如圖10(a)所示。將雙澗沙沙舌垂直流向最寬處作為雙澗沙沙舌的特征尺度，用于表征雙澗沙沙舌橫向發育的強弱，繪制雙澗沙沙舌特征長度關于特征流量的散點圖，如圖10(b)所示。靖江邊灘特征長度和流量級呈現正相關關系，且流量大于40 500 m3/s時增長速率顯著提高，靖江邊灘會出現較大規模的下移，航道工程下靖江邊灘下移速率大于自然狀態。雙澗沙沙舌在自然狀態下與流量級相關較弱，航道工程下特征尺度隨著流量增大而快速增長。通過以上分析，一方面說明年內流量過程中洪季流量對靖江邊灘和雙澗沙沙舌變形起到主導作用，另一方面，航道工程實施后，加快了靖江邊灘沖刷后退、雙澗沙沙舌橫向發育的速率，福北水道入口處水動力條件減弱，不利于沖刷下移的泥沙通過福北水道向下游輸移，造成局部區域形成礙航淺灘。

圖10 靖江邊灘及雙澗沙沙舌特征尺度隨流量級變化關系Fig. 10 Relationship between discharge and characteristic scale of point /channel bar deformation

靖江邊灘—雙澗沙沙舌—福北水道是主要的輸沙通道，基于2014年7月至2016年11月的實測數據，王建軍等[13]建立了福北水道航槽沖淤量V福北水道航槽與靖江邊灘沖淤量V靖江邊灘、雙澗沙沙舌沖淤量V雙澗沙沙舌及徑流量V徑流量的關系，見式(15)。該擬合式的實際意義是：靖江邊灘和雙澗沙沙舌沖刷下來的泥沙進入福北水道航槽并淤積，且隨著徑流量的增加航槽淤積量增多。將計算結果代入式(10)進行驗證，計算結果總體上符合擬合值的變化趨勢(相關性系數R=0.932 1)；而航道工程影響下的計算值與擬合結果偏差較大，且與擬合值呈現負相關(相關性系數R=-0.804 3)，反映了航道工程完工后，靖江邊灘、雙澗沙沙舌、福北水道的灘槽演變關系和對徑流特征響應過程發生根本性變化。式(15)不適用于描述航道工程影響下福北水道航槽沖淤量與靖江邊灘和雙澗沙沙舌、徑流特征的關系。

根據算例A1、B1、C1、D1(考慮航道工程)的計算結果，可以得到航道工程影響下的福北水道航槽沖淤量與靖江邊灘和雙澗沙沙舌沖淤量、徑流量的擬合關系式(16)(R=0.849 8)。在靖江邊灘、雙澗沙沙舌沖淤及特征流量的綜合作用下，靖江邊灘和雙澗沙沙舌每沖刷約100萬 m3泥沙，約11.9萬 m3會淤積在福北水道航槽內；徑流量每減少1 000億 m3，福北水道航槽淤積量會增加約17萬 m3。如圖11所示，航道工程作用下福北水道的淤積量隨著流量的增長而減小，特別是洪水期轉淤為沖，而維護性疏浚量反而增大，說明對航道整體而言淤積程度放緩，洪季航槽大規模回淤的現象得到顯著改善，但受到靖江邊灘沖刷下移、雙澗沙沙舌橫向淤寬和局部區域水動力減弱的不利影響，隨著流量增長，福北水道入口處區域的航道條件發生惡化，應作為深水航道維護的重點關注區域。

圖11 福北水道整體沖淤量和維護性疏浚量Fig. 11 Total erosion and deposition volume and maintenance dredging volume of Fubei channel

V福北水道航槽=-0.348(V靖江邊灘+V雙澗沙沙舌)+0.033V徑流量+93.82

(15)

V福北水道航槽=-0.119(V靖江邊灘+V雙澗沙沙舌)-0.017V徑流量+29.383

(16)

4 結 語

采用基于Hybrid LTS/GMaTS方法的水沙—床耦合數學模型對工程河道沖淤過程進行數值模擬研究。模型使用非結構三角形網格，對航道工程區域進行局部加密，通過Hybrid LTS/GMaTS方法提高了計算效率。Hybrid LTS/GMaTS方法和傳統整體最小時間步長方法的水位、流速、含沙量、床面高程變化的計算結果相對誤差較小，但該方法大幅度提高了計算效率，在本文計算工況下提高了近10倍。模型復演了2017年7月17日至18日的徑潮流水動力過程和2017年2月至2018年2月的河床沖淤過程，說明該模型可以應用于長江感潮河段的水沙數值模擬。基于8組數值算例，得到以下結論：

1)水動力作用是福姜沙水道洲灘變形的主要動力，體現在洲灘整體沖淤和局部形態改變兩方面。一方面，靖江邊灘受主流擺動影響沖淤交替，枯水期處于緩流區而淤積、洪水期處于主流區而沖刷；雙澗沙在航道工程作用下水動力條件減弱，整體上呈現淤漲趨勢。另一方面，靖江邊灘受水動力作用向下游移動，雙澗沙沙舌受水力輸沙作用橫向淤寬，洪季流量對靖江邊灘和雙澗沙沙舌變形起到主導作用。

2)航道工程作用下，福中水道入口處水動力條件明顯改善，各特征流量級下均呈現刷深趨勢；航道中段呈現“洪沖枯淤”的沖淤特征；航道工程作用下福中水道的主流被歸順在深槽內，僅在丁壩挑流作用下輕微向福姜沙左緣邊灘擺動，因此隨著流量增長，航槽向南擺動。

3)邊心灘演變和水動力作用共同影響福北水道航道沖淤變化。航道工程建成后，靖江邊灘、雙澗沙沙舌、福北水道的灘槽演變關系和對徑流特征響應過程發生根本性變化。建立了航道工程影響下的福北水道航槽沖淤量與靖江邊灘和雙澗沙沙舌、徑流特征的擬合關系式。在靖江邊灘、雙澗沙沙舌沖淤及特征流量的綜合作用下，靖江邊灘和雙澗沙沙舌每沖刷約100萬 m3泥沙，約11.9萬 m3泥沙會淤積在福北水道航槽內；徑流量每減少1 000億 m3，福北水道航槽淤積量會增加約17萬m3。航道工程實施后，對福北水道整體而言，其淤積程度放緩，特別是洪季航槽大規模回淤得到顯著改善。然而受到靖江邊灘沖刷下移、雙澗沙沙舌橫向淤寬和局部區域水動力減弱的不利影響，隨著流量增長，福北水道入口處航道條件將發生惡化，是重點疏浚維護區域。