吳明珠，潘 昀

(1.湖南省航務(wù)勘察設(shè)計研究院, 湖南 長沙 410005; 2.浙江海洋大學(xué) 國家海洋設(shè)施養(yǎng)殖工程技術(shù)研究中心， 浙江 舟山 316000)

整治工程對概化分汊河口分流和分沙影響分析

吳明珠1，潘 昀2

(1.湖南省航務(wù)勘察設(shè)計研究院, 湖南 長沙 410005; 2.浙江海洋大學(xué) 國家海洋設(shè)施養(yǎng)殖工程技術(shù)研究中心， 浙江 舟山 316000)

河口是連接河流和海洋的紐帶，也是海上運輸?shù)膯⒊厅c和目的地。河口航道的自然水深難以滿足航運的需求，因此需要對航道進(jìn)行整治和疏浚，但整治工程對分汊河口的影響機制并不明確。以明渠均勻流為驗證基礎(chǔ)，建立分汊河口不同整治工程概化水沙數(shù)學(xué)模型，探討均勻流和正規(guī)半日潮作用下一汊開挖航道和建設(shè)丁壩群后對另一汊河道的水沙動力影響。研究結(jié)果表明，斷面泥沙通量與斷面流量變化不一致，存在最優(yōu)丁壩長度，通過回歸性分析得到優(yōu)良相對長度0.62和極限長度0.78使得斷面流量和泥沙通量受丁壩影響較小。此計算模式僅能得出航道的開挖對分汊河道的分流和分沙有較大的影響。均勻流和潮流作用下周期內(nèi)斷面的分流比存在一致性，分沙比在0.5左右時各工況潮流影響略有不同，潮流對匯流處分沙比的影響值得進(jìn)一步的深入研究。研究結(jié)果為減少航道維護(hù)和整治分汊河口提供參考，也為分汊河口物理模型試驗提供理論依據(jù)。

航道； 整治工程； 分汊河口； 水沙動力

0 前言

分汊河口的江心洲往往空間尺度比較大，由于護(hù)岸工程的存在也較穩(wěn)定，所以與分汊河道在空間上有許多不同。徑流和潮流共同作用的動力條件也較為復(fù)雜，并且河口均存在一定的放寬率。航道整治工程在影響局部水沙運動的同時長時間作用后可能會對整體系統(tǒng)帶來影響，而這部分影響往往卻被忽略。

我國對航道整治工程影響下的河道水沙運動研究較多，研究深入數(shù)值計算，水槽試驗和相似物理模型試驗。談廣鳴[1,2]對分汊河道沖淤和懸移質(zhì)分沙進(jìn)行了系統(tǒng)的研究，得出單個丁壩對水沙運動的影響機制和懸移質(zhì)的分布情況。近3年對于長江口水沙動力問題逐漸開展，竇潤青[3]將長江口南槽和北槽作為一個整體系統(tǒng)，利用實測地形資料和FVCOM數(shù)值模型計算流場，研究近10余年北槽落潮分流比和分沙比的變化原因。童朝峰[4]利用分汊水槽進(jìn)行分汊口水流結(jié)構(gòu)和底沙輸沙試驗，研究不同支汊流量、不同泥沙粒徑及分汊角度對分流分沙比的影響和分汊口泥沙輸沙機理。Amaurthy等[5,6]通過分汊水槽試驗，驗證其由動量方程得出的分流比、兩分汊Fr數(shù)和水深的關(guān)系。Barkdoll等[7]對直角分汊口的底沙運動做了初步測試。Neary等[8]人實驗測定了直角分汊口的水流結(jié)構(gòu)，得到一些分流面的變化規(guī)律。

從國內(nèi)外對分汊河口水沙運動研究來看，動力從均勻流到潮流逐漸由簡單到復(fù)雜對分汊河口進(jìn)行認(rèn)識。以往的研究多抓住航道整治工程內(nèi)部，比如丁壩附近水流結(jié)構(gòu)和懸疑質(zhì)分布或者航道灘槽環(huán)流結(jié)構(gòu)和泥沙運動等細(xì)節(jié)，對內(nèi)在機制進(jìn)行深入研究，為數(shù)學(xué)模型和物理模型尋求可靠的理論依據(jù)。近年來許多河口因一汊整治工程的存在對另一汊帶來許多不利影響，并且這方面的理論基礎(chǔ)多屬空白。于是非常有必要把分汊河口的兩汊作為一個整體系統(tǒng)進(jìn)行研究，探討整治工程一汊與另一汊的一般性相互作用機制。因此本文首先詳細(xì)闡述全沙均勻流數(shù)學(xué)模型，并以此為基礎(chǔ)率定各種參數(shù)，然后建立分汊河口不同整治工程概化水沙數(shù)學(xué)模型，探討均勻流和正規(guī)半日潮作用下一汊開挖航道和建設(shè)丁壩群后對另一汊河道的水沙動力影響。為對分汊河口水沙運動的進(jìn)一步研究提供參考。

1 計算方法

建立非結(jié)構(gòu)網(wǎng)格概化模型，解析水動力方程和全沙非平衡輸沙方程探討均勻流和正規(guī)半日潮作用下一汊開挖航道和建設(shè)丁壩群后對另一汊河道的影響以及分汊河口因航道整治引起的水沙動力響應(yīng)。

1.1 網(wǎng)格建立

Delaunay三角剖分算法是構(gòu)建數(shù)字高程模型(DEM)的主要算法，具有快速、高效的內(nèi)存交換特點，可以穩(wěn)定的計算用戶指定的孔和凹槽等特殊的三角形區(qū)域。三角形算法在網(wǎng)站上免費提供源代碼，提供許多包括三角剖分的算法和數(shù)據(jù)結(jié)構(gòu)實施方法，以及建立和完善網(wǎng)格的步驟，解決算法中出現(xiàn)穩(wěn)定性問題。

1.2 解析方法

模型基于非結(jié)構(gòu)網(wǎng)格和非壓縮雷諾方法調(diào)用Boussinesq假設(shè)和靜水壓力解析平均的斯托克斯方程實現(xiàn)。因此模型涉及連續(xù)方程，動量方程和密度方程，通過平均紊流方程進(jìn)行閉合求解。方程組的空間離散采用中心有限體積法進(jìn)行。淺水方程和運動方程的時間積分有兩種方法低階和高階，低階的方法是一階顯式歐拉方法，高階的方法是二階Runge-Kutta方法。

沿閉邊界(陸地邊界)的所有變量的強迫控制為零。開邊界可以指定定量的流量和水位作為水動力的邊界條件。運動方程邊界可以指定定量值或者梯度值。

1.2.1 淺水方程

二維的淺水方程的動量方程(僅給出u方向，即水平方向)和連續(xù)方程的微分形式如下：

(1)

(2)

式中上面短線代表平均值，方程(2)中Txx、Txy分別為粘滯摩擦力和紊動摩擦力,省略的v方向方程存在Tyy為對流摩擦力?？梢圆捎没谄骄魉偬荻鹊奈蓜诱硿匠逃嬎悖?/p>

(3)

1.2.2 床面剪切力

(4)

(5)

二維的ub平均水深處的流速計算，阻力系數(shù)可根據(jù)謝才系數(shù)C和曼寧系數(shù)M計算：

(6)

1.2.3 泥沙運動控制方程

泥沙運動方程根據(jù)水流的對流擴散方程可推導(dǎo)如下：

(7)

泥沙的對流擴散方程根據(jù)QUICKEST格式[9]再推導(dǎo)的三階顯式有限差分ULTIMATE格式[10]計算。

1.3 參數(shù)確定

模型可以模擬水流作用下細(xì)顆粒(<63 μm )粉沙和淤泥的沖刷、運動和淤積。為了精確地模擬泥沙的沖刷過程，床面特意添加了固結(jié)層和半固結(jié)層，沖刷和淤積均以床面剪切力進(jìn)行判別。模型中涉及的主要參數(shù)有：泥沙顆粒沉速，懸浮泥沙濃度，水平流速，垂向流速，懸沙斷面分布形式和臨界沖刷淤積剪切力。

1.3.1 淤積和沖刷

模型中添加有隨機系統(tǒng)以體現(xiàn)泥沙淤積的隨機性，所以淤積公式為：SD=wscbpd，式中ws為泥沙顆粒沉速，m/s；cb為近底含沙量，kg/m3；pd為泥沙淤積幾率。淤積幾率pd可表示為：

pd=1-τb/τcd,τb≤τcd

( 8)

式中τb為床面剪切力，N/m2；τcd為臨界淤積床面剪切力，N/m2。

泥沙方程通過有水動力計算出的床面剪切力來描述。床面剪切力根據(jù)對數(shù)流速分布以及床面阻力高度可表示為：

τb=ρfcV2/2

(9)

τb為床面剪切力,N/m2；ρ為流體密度，kg/m3；V為平均流速,m/s；fc為水流阻力系數(shù)，可根據(jù)經(jīng)驗公式計算：

(10)

h為水深,m；k為床面阻力高度,m。

床面沖刷分為固結(jié)層的沖刷和半固結(jié)或者不固結(jié)層的沖刷。

對于固結(jié)沖刷(SE)：

(11)

式中E為單位面積時間沖刷系數(shù),kg/(m2·s-1)；τb為床面剪切力，N/m2；τce為臨界沖刷床面剪切力,N/m2；n沖刷指數(shù)系數(shù),一般取1。

對于半固結(jié)或者不固結(jié)層的沖刷：

(12)

1.3.2 泥沙顆粒沉速

細(xì)顆粒泥沙的沉速取決于泥沙顆粒，絮凝，溫度，懸浮物濃度和有機物含量。如今一般存在2組沉速機理：沉速隨著泥沙濃度增加而增加和沉速隨著泥沙濃度增加而減小。后者主要與泥沙的沉降的阻滯作用有關(guān)，大多數(shù)河口符合前者。模型中只定義定常值沉速0.0001m/s。

1.3.3 濃度剖面形式

(13)

2 模型計算

數(shù)值模型平面布置、范圍及床面高程如圖1所示，河道上游分為南支和北支兩汊，下游相匯。數(shù)值水槽總長240 km，中間江心洲105 km，進(jìn)出口河道寬15 km，分汊河道寬7.5 km。底床床面為上游-5 m到下游-7.8 m的矩形斷面斜坡，河床等深線如圖1所示。北支設(shè)有長105 km和寬400 m的航道，以及間隔為10 km的對丁壩群N1-N5和S1-S5，初始單個丁壩的尺寸為長1800 m，寬200 m。最大網(wǎng)格長度800 m，最小網(wǎng)格長度200 m，如圖2所示。最小和最大時間步長為0.01 s和30 s，計算時間足夠長。

初始丁壩群的計算網(wǎng)格如圖2所示，先以數(shù)值水槽均勻流率定參數(shù)，然后添加泥沙和修改航道整治建筑物尺寸以探討均勻流作用下河道不沖不淤時的水沙動力變化。采用M2調(diào)和正規(guī)半日潮潮位作為下游邊界條件，計算由均勻流率定出的參數(shù)探討水沙變化規(guī)律。

2.1 均勻流驗證

上游流量14600 m3/s和下游水位1 m為初始條件計算水槽水動力，其中曼寧系數(shù)為常值32 N/m2。計算出的穩(wěn)定流速場和流速矢量如圖3(A)，除進(jìn)流口和匯流口流速有不同外，其他區(qū)域均為均勻流，大小為0.45 m/s。

圖1 模型地形和監(jiān)測斷面

圖2 模型計算網(wǎng)格

根據(jù)式(4)計算床面剪切力為0.27 N/m2。為了保證驗證模型中底床不被沖刷所以取沖刷臨界剪切力為0.27 N/m2，同時為了保證上游的懸疑質(zhì)一旦淤積不會被再懸浮所以選擇淤積臨界切應(yīng)力也為0.27 N/m2。計算中，取泥沙沉速為常值0.0001 m/s，假設(shè)只有一層底床且密度為300 kg/m3和糙率高度為0.001 m。泥沙初始邊界條件上游給定泥沙濃度為0.1 kg/m3，下游泥沙濃度0.01 kg/m3。計算出的穩(wěn)定懸疑質(zhì)泥沙濃度分布如圖3(B)，泥沙濃度從上游往下游分布依次減小。南支和北支懸疑質(zhì)泥沙濃度基本分布對稱，因網(wǎng)格疏密影響使得北支0.08 kg/m3等值線稍長。從圖3計算結(jié)果來看，上述給出的參數(shù)可以為水槽模型進(jìn)一步的計算提供合理依據(jù)。

圖3 流場和懸沙濃度分布圖

2.2 概化方案設(shè)計

航道整治工程有-10 m深航道和丁壩群，其中丁壩長度從初始1.8 km依次每次增加0.2 km至3.6 km。為了探究丁壩壩田區(qū)淤積后再減小該汊水量和增大阻力的同時對進(jìn)流口和匯流口的水沙的影響變化，把壩田區(qū)地形增加3 m作為一種計算方案，如圖1所示。動力分為上述已經(jīng)驗證的均勻流(流速0.45 m/s)和M2調(diào)和的潮位0 m至2 m的正規(guī)半日潮。具體計算方案依次為：①丁壩群(丁壩長度0、1.8、2.0、2.2、2.4、2.6、2.8、3.0、3.2、3.4、3.6 km)；②丁壩群(0～3.6 km)&航道(航深10 m)；③丁壩群(0～3.6 km)&壩田淤積(淤積3 m)；④丁壩群(0～3.6 km)&壩田淤積(淤積高度3 m)&航道(航深10 m)。

監(jiān)測和統(tǒng)計流量和泥沙通量的斷面如圖1，上游相互對稱位置的A斷面和B斷面模型中稱為分汊口，A斷面恰好在航道的左頂端處，前斷面模型中稱為分流口和A斷面相距2.5 km。下游斷面與上游相似，記為C斷面、D斷面和后斷面，模型中稱為出汊口和匯流口。

3 水沙動力響應(yīng)

分析均勻流作用下航道整治工程對模型水沙動力的影響及分流比和分沙比的變化。通過分析均勻流作用下和潮流作用下分流比及分沙比的相關(guān)性，探尋兩者的異同。

3.1 均勻流作用下水沙動力響應(yīng)

圖4所示，均勻流作用下隨著丁壩長度的增加A斷面的泥沙通量逐漸減少，B斷面的泥沙通量逐漸增大。C斷面和D斷面的泥沙通量相比較，D斷面增大的幅度大數(shù)倍，基本與后斷面泥沙通量保持一致。

圖4和圖5可以得出工況丁壩群和工況壩田淤積的流量和泥沙通量曲線基本重合，工況航道和工況壩田淤積&航道的流量和泥沙通量曲線基本重合。但前2組和后2組工況卻有明顯的差別，說明航道的開挖對分汊河道的分流和分沙有較大的影響。

3.2 均勻流作用下分流比和分沙比變化

如圖4和圖5橫坐標(biāo)為丁壩長度與一汊河道寬度的比值，縱坐標(biāo)為單汊河道流量或泥沙通量與兩汊流量或泥沙通量和的比值，定義為分流比和分沙比。

圖4 分流比

圖5 分沙比

受相對丁壩長度的影響，有工程一汊分流比或分沙比減小，另一汊增大。前后斷面的分流比對稱性較好，分沙比不然，在兩汊道內(nèi)發(fā)生較多泥沙的沖刷或淤積。存在航道的2組工況對分汊河道的分流比和分沙比有較大的影響。

3.3 均勻流作用和潮流作用相關(guān)性分析

把均勻流作用下A斷面至D斷面分流比和分沙比作為x軸和潮流作用下的全部所有分流比和分沙比作為y軸標(biāo)記如圖6。對于分流比，僅有C斷面和D斷面2個點稍微偏離45°線，說明無論是均勻流作用下還是潮流作用下周期內(nèi)所有工況斷面的分流比一致性非常好。對于分沙比，C斷面和D斷面分沙比在0.5左右時各工況潮流影響略有不同。整體可總結(jié)，正規(guī)潮流周期內(nèi)對分汊河口的影響作用與均勻流作用相當(dāng)。潮流導(dǎo)致C斷面和D斷面分沙比的影響值得進(jìn)一步的深入研究。

圖6 潮流和均勻流作用下分流比和分沙比相關(guān)性分析

3.4 分流比和分沙比回歸分析

由于A斷面和B斷面的分流比與C斷面和D斷面的相等，所以僅對前一組斷面分流比進(jìn)行曲線擬合，而分沙比前后兩組均進(jìn)行了曲線擬合。實際河口中由于同時存在多種航道整治工程和多種動力因素，為了提高所得曲線的實用性和指導(dǎo)價值，所以在曲線擬合時選擇了所有工況的數(shù)據(jù)點。具體曲線方程和相關(guān)系數(shù)如圖7、圖8和圖9所示。

圖7 所有工況A斷面和B斷面分流比擬合曲線

圖8 所有工況A斷面和B斷面分沙比擬合曲線

圖9 所有工況C斷面和D斷面分沙比擬合曲線

根據(jù)分流比和分沙比的定義可知兩式的和應(yīng)為1并且各項系數(shù)應(yīng)相等。式(14)、式(15)分別為所有工況A斷面和B斷面分流比與分沙比數(shù)據(jù)擬合曲線。式(16)為所有工況C斷面和D斷面分沙比數(shù)據(jù)擬合曲線。

yAd=0.51e-0.17x-0.00005e9.1x

yBd=0.49 e-0.17x+0.00005 e9.1x

(14)

yAs=0.53e0.04 x-0.0088 e4.2x

yBs=0.47 e0.04 x+0.0088 e4.2x

(15)

yCs= 0.18+0.36cos(2.86x)+

0.3sin(2.86x)

yDs=0.82-0.36cos(2.86x)-

0.3sin(2.86x)

(16)

對式(14)和式(15)曲線進(jìn)行求導(dǎo)求出相對丁壩長度對應(yīng)的導(dǎo)數(shù)即曲線在該點處的斜率，如圖10。分流比隨著丁壩長度的增加在小于0.7時曲線斜率緩慢減小，之后指數(shù)性減小。分沙比隨著丁壩長度的增加接近于線性減小，比分流比的變化快。取斜率-0.5和-1作為判別丁壩長度對兩汊河道允許影響的標(biāo)準(zhǔn)得到優(yōu)良相對長度0.62和極限長度0.78。

圖10 A斷面擬合曲線的斜率

4 結(jié)論

以明渠均勻流為驗證基礎(chǔ)，建立分汊河口不同整治工程概化水沙數(shù)學(xué)模型，探討均勻流和正規(guī)半日潮作用下一汊開挖航道和建設(shè)丁壩群后對另一汊河道的水沙動力影響。斷面泥沙通量與斷面流量變化不一致，存在最優(yōu)丁壩長度使得斷面流量和泥沙通量受丁壩影響較小。斷面流量受潮流的影響，距離潮流邊界越近受潮流影響越強。均勻流和潮流作用下周期內(nèi)斷面的分流比存在一致性，分沙比在0.5左右時各工況潮流影響略有不同。潮流對匯流處分沙比的影響值得進(jìn)一步深入研究。對所有工況的分流比和分沙比進(jìn)行回歸性分析，取文中變化幅值參數(shù)-0.5和-1作為判別丁壩長度對兩汊河道允許影響的指標(biāo)得到優(yōu)良相對長度0.62和極限長度0.78。

[1] 張悅, 談廣鳴, 余新明, 等. 丁壩對分汊河道分流比的影響[J]. 中國農(nóng)村水利水電, 2007(4):91-93.

[2] 余新明, 談廣鳴, 張悅, 等. 分汊河道水沙輸移特征試驗[J]. 武漢大學(xué)學(xué)報(工學(xué)版), 2007,40(4):9-13.

[3] 竇潤青, 郭文云, 葛建忠, 等. 長江口北槽落潮分流比變化原因分析[J]. 華東師范大學(xué)學(xué)報(自然科學(xué)版), 2014(3):93-104.

[4] 童朝鋒.分汊口水沙運動特征及三維水流數(shù)學(xué)模型應(yīng)用研究[D].南京: 河海大學(xué), 2005.

[5] Amurthur S. Ramamurth, Mysore.G. Satish. Division of Flow in Short Open Channel Branches[J].Journal of hydraulic engineering, 1988, 114(44):428-438.

[6] Amurthur S. Ramamurth, Due Minh Tran, Luis B.Carballada.Division Flow in Open Channels[J].Journal of hydraulic engineering. 1990, 116(3):449-455.

[7] Brian D. Barkdoll, Robert Ertema, A. Jacob Odgaard. Sediment Control at Lateral Divisions: Limit Sand Enhancements to Vane Use[J].Journal of Hydraulic Engineering, 1999, 125(8):362-870.

[8] V. S. Neary, A. Jacob Odgaard. Three Dimensional Flow Structure at Open Channel Divisions[J].Journal of Hydraulic Engineering, 1993, 119(11):1223-1230.

[10] Leonard B P. A stable and accurate convective modelling procedure based on quadratic upstream interpolation[J].Computer methods in applied mechanics and engineering, 1979, 19(1): 59-98.

2017-01-13

國家自然科學(xué)基金( 51479122)

吳明珠( 1979-) ，男，工程師，主要從事港口與航道工程研究。

1008-844X(2017)01-0166-06

U 617

A