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

吳明珠，潘 昀

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

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

吳明珠1，潘 昀2

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

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

航道； 整治工程； 分汊河口； 水沙動(dòng)力

0 前言

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

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

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

1 計(jì)算方法

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

1.1 網(wǎng)格建立

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

1.2 解析方法

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

沿閉邊界(陸地邊界)的所有變量的強(qiáng)迫控制為零。開(kāi)邊界可以指定定量的流量和水位作為水動(dòng)力的邊界條件。運(yùn)動(dòng)方程邊界可以指定定量值或者梯度值。

1.2.1 淺水方程

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

(1)

(2)

式中上面短線代表平均值，方程(2)中Txx、Txy分別為粘滯摩擦力和紊動(dòng)摩擦力,省略的v方向方程存在Tyy為對(duì)流摩擦力?？梢圆捎没谄骄魉偬荻鹊奈蓜?dòng)粘滯方程計(jì)算：

(3)

1.2.2 床面剪切力

(4)

(5)

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

(6)

1.2.3 泥沙運(yùn)動(dòng)控制方程

泥沙運(yùn)動(dòng)方程根據(jù)水流的對(duì)流擴(kuò)散方程可推導(dǎo)如下：

(7)

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

1.3 參數(shù)確定

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

1.3.1 淤積和沖刷

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

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

( 8)

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

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

τb=ρfcV2/2

(9)

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

(10)

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

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

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

(11)

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

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

(12)

1.3.2 泥沙顆粒沉速

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

1.3.3 濃度剖面形式

(13)

2 模型計(jì)算

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

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

2.1 均勻流驗(yàn)證

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

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

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

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

圖3 流場(chǎng)和懸沙濃度分布圖

2.2 概化方案設(shè)計(jì)

航道整治工程有-10 m深航道和丁壩群，其中丁壩長(zhǎng)度從初始1.8 km依次每次增加0.2 km至3.6 km。為了探究丁壩壩田區(qū)淤積后再減小該汊水量和增大阻力的同時(shí)對(duì)進(jìn)流口和匯流口的水沙的影響變化，把壩田區(qū)地形增加3 m作為一種計(jì)算方案，如圖1所示。動(dòng)力分為上述已經(jīng)驗(yàn)證的均勻流(流速0.45 m/s)和M2調(diào)和的潮位0 m至2 m的正規(guī)半日潮。具體計(jì)算方案依次為：①丁壩群(丁壩長(zhǎng)度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)測(cè)和統(tǒng)計(jì)流量和泥沙通量的斷面如圖1，上游相互對(duì)稱位置的A斷面和B斷面模型中稱為分汊口，A斷面恰好在航道的左頂端處，前斷面模型中稱為分流口和A斷面相距2.5 km。下游斷面與上游相似，記為C斷面、D斷面和后斷面，模型中稱為出汊口和匯流口。

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

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

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

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

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

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

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

圖4 分流比

圖5 分沙比

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

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

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

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

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

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

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

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

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

根據(jù)分流比和分沙比的定義可知兩式的和應(yīng)為1并且各項(xiàng)系數(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)

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

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

4 結(jié)論

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

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

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

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

[4] 童朝鋒.分汊口水沙運(yùn)動(dòng)特征及三維水流數(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

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

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

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

U 617

A