姜治兵,陸 虹,楊青遠(yuǎn)

（長(zhǎng)江科學(xué)院水力學(xué)研究所,武漢 430010)

?

海水經(jīng)船閘入侵淡水運(yùn)河的三維數(shù)值模擬

姜治兵,陸 虹,楊青遠(yuǎn)

（長(zhǎng)江科學(xué)院水力學(xué)研究所,武漢 430010)

摘 要：為研究船閘運(yùn)行過(guò)程中鹽分的輸運(yùn)規(guī)律及其對(duì)淡水水域的鹽度影響,建立了三維k-ε兩相混合流數(shù)值模型,模型控制方程組采用有限體積法進(jìn)行離散,流速與壓力耦合采用SIMPLEC算法,時(shí)間項(xiàng)采用一階隱式格式,流項(xiàng)采用二階迎風(fēng)格式,計(jì)算區(qū)域采用六面體網(wǎng)格劃分。采用某船閘的海水入侵原型試驗(yàn)成果對(duì)數(shù)值模型進(jìn)行了驗(yàn)證,閘室內(nèi)鹽度的模擬值與實(shí)測(cè)值吻合較好。海水入侵淡水水域模擬結(jié)果表明:上游航道鹽度分布可分為異重流段、過(guò)渡段和擴(kuò)散段,各段鹽分輸運(yùn)速度依次減?。桓鲾嗝娴柠}度呈周期性變化,在船閘運(yùn)行一段時(shí)間后逐步趨于動(dòng)態(tài)平衡；一個(gè)循環(huán)內(nèi)的鹽分入侵量在船閘運(yùn)行初期較大,隨著運(yùn)行時(shí)間的增加逐漸減小,并趨于恒定值。研究結(jié)果可為船閘的設(shè)計(jì)與運(yùn)行調(diào)度提供科學(xué)依據(jù)。

關(guān)鍵詞：船閘；海水入侵；數(shù)值模型；鹽度分布；鹽分入侵量

2016,33（02):52-56

1 研究背景

船舶過(guò)閘過(guò)程中,海水中的鹽分在異重流作用下向淡水水域入侵,咸化淡水。如果水體鹽濃度超標(biāo),將會(huì)影響淡水水域的水質(zhì)及淡水生物的生存環(huán)境。因此,研究分析船閘運(yùn)行過(guò)程中鹽分的輸運(yùn)規(guī)律及其對(duì)淡水水域的鹽度影響,對(duì)鹽水入侵緩解措施的制定及船閘的設(shè)計(jì)與運(yùn)行調(diào)度有著重要的意義。

擬建的尼加拉瓜運(yùn)河連接太平洋、尼加拉瓜湖與大西洋,通航后將大大縮短太平洋至大西洋的航運(yùn)距離。運(yùn)河上擬建4個(gè)船閘,其中1?！?＃船閘依次位于連接尼加拉瓜湖和大西洋的河段；4＃船閘位于連接太平洋與尼加拉瓜湖的河段,上距尼加拉瓜湖22 km,下距太平洋約2 km。4＃船閘分為三級(jí),為節(jié)約淡水資源,每級(jí)閘室配備一個(gè)省水池。

由于上下游航道、閘室、船閘充水管路三者尺度存在數(shù)量級(jí)上的差異,且船閘運(yùn)行程序復(fù)雜,建立完全的三維數(shù)值模型需要數(shù)量龐大的計(jì)算網(wǎng)格,在此基礎(chǔ)上進(jìn)行長(zhǎng)時(shí)間的非恒定流模擬,計(jì)算量巨大,不具備可行性。在充分分析鹽分輸運(yùn)特點(diǎn)與前人研究成果［1-5］的基礎(chǔ)上,建立了耦合的鹽分交換分析模型和三維k-ε兩相混合流數(shù)值模型,研究分析船閘運(yùn)行過(guò)程中鹽分的輸運(yùn)規(guī)律。鹽分交換分析模型的計(jì)算范圍從下游引航道至第一閘室,包括各級(jí)閘室的省水池；三維k-ε多相流數(shù)值模型模擬區(qū)域包括第一閘室與22 km的上游航道,以及部分尼加拉瓜湖湖區(qū)。本文著重介紹三維k-ε兩相混合流數(shù)值模型的研究成果。

2 控制方程與數(shù)值方法

2.1 控制方程

模型采用多相流混合相的控制方程:

連續(xù)方程

動(dòng)量方程

2.2 數(shù)值方法

將控制方程寫(xiě)為通用格式

對(duì)通用變量在控制體上取平均,則方程（6)變?yōu)?/p>

式中:m為單元控制體的單元面總數(shù)；Aj為單元面j的面積；為單元控制體的源項(xiàng)平均值；Fj（Φ)Aj為單元面的法向通量,包括對(duì)流通量與擴(kuò)散通量。

方程組采用有限體積法進(jìn)行離散,流速與壓力耦合采用SIMPLEC算法,時(shí)間項(xiàng)采用一階隱式格式,對(duì)流項(xiàng)采用二階迎風(fēng)格式。

3 模型驗(yàn)證

海水入侵閘室或閘室咸水入侵淡水水域?yàn)楫愔亓黩?qū)動(dòng)下的非恒定過(guò)程,入侵鹽量與時(shí)間密切相關(guān),因此數(shù)值模型需滿足時(shí)間與入侵鹽量的精度要求。

參考文獻(xiàn)［5］在一座兩端分別與淡水湖和海洋相接的單級(jí)船閘上進(jìn)行了海水入侵淡水的試驗(yàn),船閘閘室長(zhǎng)196 m、寬30 m、水深10 m。試驗(yàn)分2種情況進(jìn)行,一種為海水入侵淡水閘室,另一種為閘室內(nèi)海水入侵淡水湖。

3.1 海水入侵淡水閘室

此試驗(yàn)的現(xiàn)場(chǎng)布置見(jiàn)圖1（a),閘門(mén)開(kāi)啟前閘室內(nèi)為淡水,閘門(mén)分別開(kāi)放5,10,15,20,30 min后關(guān)閉,待閘室內(nèi)淡水和咸水混合均勻,再測(cè)量閘室內(nèi)鹽度。閘門(mén)開(kāi)放時(shí)間與閘室內(nèi)平均鹽度的變化曲線見(jiàn)圖1（b)中的虛線。

圖1 海水入侵淡水閘室Fig.1 Seawater intrusion into freshwater lock chamber

三維k-ε數(shù)值模型模擬結(jié)果見(jiàn)圖1（b)中的實(shí)線。從圖1（b)可以看出:閘室內(nèi)的平均鹽度隨著閘門(mén)開(kāi)放時(shí)間的增加而增大,閘門(mén)開(kāi)放30 min條件下,閘室內(nèi)混合均勻后的水體鹽度已超過(guò)30‰；模型計(jì)算值與實(shí)測(cè)值吻合較好。

3.2 閘室內(nèi)海水入侵淡水湖

此試驗(yàn)的現(xiàn)場(chǎng)布置見(jiàn)圖2（a),閘門(mén)開(kāi)啟前閘室內(nèi)為鹽度35‰的海水,閘門(mén)分別開(kāi)放5,10,15,20, 30 min后關(guān)閉,待閘室內(nèi)淡水和咸水混合均勻,再測(cè)量閘室內(nèi)鹽度,閘門(mén)開(kāi)放時(shí)間與閘室內(nèi)平均鹽度的變化曲線見(jiàn)圖2（b)中的虛線。

圖2 閘室內(nèi)海水入侵淡水湖Fig.2 Seawater intrusion into freshwater lake

三維k-ε數(shù)值模型模擬結(jié)果見(jiàn)圖2（b)中的實(shí)線。從圖2（b)可以看出,閘室內(nèi)的平均鹽度隨著閘門(mén)開(kāi)放時(shí)間的增加而減小,閘門(mén)開(kāi)放30min條件下,閘室內(nèi)混合均勻后的水體鹽度約為5‰；模型計(jì)算值與實(shí)測(cè)值吻合較好。

以上驗(yàn)證結(jié)果表明,本文三維k-ε數(shù)值模型模擬結(jié)果精度較高,可用于后述研究。

4 模擬區(qū)域與計(jì)算網(wǎng)格及邊界條件

4.1 模擬區(qū)域與計(jì)算網(wǎng)格

模擬區(qū)域包括第一閘室、第一閘首、上游航道（22 km)、第一閘室與引航道間的充水管道以及尼加拉瓜湖的部分區(qū)域,采用六面體網(wǎng)格劃分,網(wǎng)格總數(shù)約60萬(wàn)個(gè),見(jiàn)圖3。

圖3 三維模型布置Fig.3 Layout of the three-dimensional model

船閘運(yùn)行過(guò)程中閘室、省水池和上游引航道內(nèi)的鹽度互相影響且不斷變化,分析模型和三維模型在模擬計(jì)算中需進(jìn)行有關(guān)參數(shù)實(shí)時(shí)傳遞,兩模型間的參數(shù)傳遞發(fā)生在第一閘室,如圖4所示。

圖4 模型計(jì)算區(qū)域Fig.4 Schematic diagram of area for model computation

4.2 邊界條件與初始條件

（1)鹽度邊界條件。第一閘室及尼加拉瓜湖為循環(huán)邊界。每次第一閘首處閘門(mén)開(kāi)啟前,根據(jù)鹽分交換分析模型計(jì)算結(jié)果重置第一閘室鹽度；由于模擬湖水區(qū)域有限,入湖鹽分在較小區(qū)域的聚集會(huì)造成航道出口處鹽度虛高,影響航道鹽度分布計(jì)算結(jié)果,考慮到尼加拉瓜湖可視為無(wú)限水體,入侵鹽分入湖后在風(fēng)生流等的驅(qū)動(dòng)下與湖水完全混合,船閘運(yùn)行初期,鹽分的入侵量有限,湖水鹽度相對(duì)于航道內(nèi)鹽度可以忽略,因此在每次循環(huán)結(jié)束后對(duì)湖水鹽度置0。

（2)水動(dòng)力邊界條件。尼加拉瓜湖上游入口、充水管道出口設(shè)置為流速邊界,第一閘室充水時(shí)段,兩處邊界的流速按充水流量給定,其它時(shí)段流速為零；自由水面按對(duì)稱邊界處理。

（3)動(dòng)邊界。閘門(mén)啟閉按動(dòng)邊界處理,采用動(dòng)網(wǎng)格技術(shù)進(jìn)行模擬。

（4)初始條件。各閘室、省水池和航道內(nèi)鹽度、流速均為0。

5 模擬結(jié)果

由耦合模型按設(shè)計(jì)提供的船閘運(yùn)行調(diào)度方式模擬計(jì)算各閘室、省水池以及航道內(nèi)的鹽度變化及輸運(yùn)過(guò)程。

5.1 閘首附近流態(tài)

第一閘室充水時(shí)段,第一閘室閘門(mén)開(kāi)啟時(shí)段與鹽水自由上侵時(shí)段,閘首附近呈現(xiàn)不同流態(tài)。

第一閘室充水時(shí)段的鹽度分布與流態(tài)見(jiàn)圖5。一次充水持續(xù)時(shí)間為8 min,此時(shí)段充水管進(jìn)流處呈典型的“匯”流流態(tài),縱剖面上形成匯流漏斗。此時(shí)段引航道內(nèi)水體向下游閘室方向流動(dòng),對(duì)鹽水上侵有一定的擬制作用。

圖5 充水時(shí)段閘首處流態(tài)及鹽度分布Fig.5 Flow pattern at lock head and salinity distribution in the period of water filling

圖6 閘門(mén)開(kāi)放時(shí)段閘首處流態(tài)及鹽度分布Fig.6 Flow pattern at lock head and salinity distribution in the period of gate opening

閘門(mén)開(kāi)啟時(shí)段的鹽度分布與流態(tài)見(jiàn)圖6。閘門(mén)開(kāi)放一次時(shí)間為31 min,閘門(mén)開(kāi)啟后閘室內(nèi)高鹽度水體與引航道內(nèi)低鹽度水體在密度差的驅(qū)動(dòng)下形成異重流,閘室內(nèi)高鹽度水體從底部潛入引航道,引航道內(nèi)低鹽度水體從上部進(jìn)入閘室,閘首附近呈現(xiàn)典型的鹽水楔,高鹽度水體與低鹽度水體間在縱剖面上形成回流。鹽分向引航道的輸運(yùn)依靠異重流在此時(shí)段內(nèi)完成。

鹽水自由上侵時(shí)段的鹽度分布與流態(tài)見(jiàn)圖7。此時(shí)段持續(xù)時(shí)間相對(duì)較長(zhǎng),閘首附近仍為異重流流態(tài),鹽分由高鹽度區(qū)向上游的低鹽度區(qū)輸運(yùn),閘首附近鹽度逐步降低,航道鹽度逐步升高。

圖7 鹽水自由上侵時(shí)段閘首處的流態(tài)及鹽度分布Fig.7 Flow pattern at lock head and salinity distribution in the period of free upward saltwater intrusion

5.2 斷面鹽度過(guò)程

在航道0～20 km范圍內(nèi)間隔2.5 km布置鹽度監(jiān)測(cè)斷面（斷面名稱為監(jiān)測(cè)斷面距第一閘首的距離),各斷面平均鹽度變化過(guò)程見(jiàn)圖8。

圖8 航道內(nèi)各斷面平均鹽度過(guò)程Fig.8 Variations of average salinity with time at different sections in the channel

從圖8可以看出:在船閘運(yùn)行的最初10 d內(nèi),閘首處（0 km斷面)鹽度迅速上升,10～40 d內(nèi)鹽度上升幅度逐漸變小,40 d后在13‰～16‰范圍內(nèi)變動(dòng)；2.5 km斷面鹽度在3～15 d內(nèi)鹽度迅速上升,15 d后增長(zhǎng)速度逐漸減慢,40 d后達(dá)到10‰左右；5 km斷面與7.5 km斷面在40 d鹽度分別達(dá)到7‰與3.5‰,隨后保持緩慢增長(zhǎng)；10 km及以上斷面鹽度在3～4 d后陸續(xù)開(kāi)始增長(zhǎng),隨后基本達(dá)到動(dòng)態(tài)平衡,鹽度在一個(gè)較小范圍內(nèi)波動(dòng),計(jì)算時(shí)段內(nèi)均未超過(guò)0.2‰。

5.3 縱向鹽度分布

圖9為40 d左右自由上侵時(shí)段某時(shí)刻航道的垂向鹽度分布圖。按照鹽度分布特征,可將上游航道分為異重流段、過(guò)渡段和擴(kuò)散段。

圖9 上游航道縱剖面鹽度分布Fig.9 Salinity distribution along longitudinal profile in upstream waterway

閘首至上游約9 km段為異重流段,從圖9（a)可以看出,鹽水楔由下游向上游延伸,鹽度等值線下起水面,上至航道底部,此段內(nèi)鹽分向上游的輸運(yùn)主要由異重流驅(qū)動(dòng),輸運(yùn)速度較快。

航道上游9 km附近為過(guò)渡段,從圖9（b)可以看出,鹽度的垂向分布由楔形分布過(guò)渡到成層分布,鹽分的輸運(yùn)由異重流驅(qū)動(dòng)轉(zhuǎn)向由擴(kuò)散驅(qū)動(dòng)為主,輸運(yùn)速度減緩。分析過(guò)渡段的鹽度等值線動(dòng)態(tài)過(guò)程可發(fā)現(xiàn),閘門(mén)開(kāi)放時(shí)段與鹽水自由上侵時(shí)段鹽水楔前鋒緩慢向上游移動(dòng),閘室充水時(shí)段鹽水楔前鋒短時(shí)向下游移動(dòng),鹽水楔前鋒基本穩(wěn)定在9 km附近。

上游約9 km至航道入湖口為擴(kuò)散段,從圖9（c)可以看出,鹽度等值線基本平行于水面,鹽度表現(xiàn)為成層分布規(guī)律,鹽分依靠濃度差由高鹽度水體向低鹽度水體擴(kuò)散,輸運(yùn)速度較慢。

圖10 鹽分入侵量過(guò)程線Fig.10 Variation of saline intrusion quantity with time

5.4 鹽分入侵量

以一個(gè)循環(huán)（4次上行和4次下行)為統(tǒng)計(jì)周期,得到圖10所示每個(gè)循環(huán)鹽分入侵量的變化過(guò)程。從圖10可以看出,船閘運(yùn)行初期閘門(mén)兩側(cè)水體密度差最大,入侵鹽量很大,每個(gè)循環(huán)可達(dá)15 000 t,隨著船閘運(yùn)行時(shí)間的增加,閘門(mén)兩側(cè)水體鹽度差減小,入侵鹽量亦逐漸減小,趨于恒定值,在船閘運(yùn)行45 d時(shí),每個(gè)循環(huán)入侵鹽量在1 800 t左右。

6 結(jié) 論

采用三維k-ε兩相混合流數(shù)值模型,模擬了海水經(jīng)船閘入侵淡水航道與尼加拉瓜湖的過(guò)程,主要研究結(jié)論如下:

（1)驗(yàn)證結(jié)果表明,數(shù)值模型入侵鹽量模擬精度較高,閘室內(nèi)鹽度的模擬值與實(shí)測(cè)值吻合較好；

（2)上游航道鹽度分布可分為異重流段、過(guò)渡段和擴(kuò)散段,各段鹽分輸運(yùn)速度依次減?。?/p>

（3)各斷面的鹽度呈周期性變化,在船閘運(yùn)行一段時(shí)間后逐步趨于動(dòng)態(tài)平衡；

（4)一個(gè)循環(huán)內(nèi)的鹽分入侵量在船閘運(yùn)行初期較大,隨著運(yùn)行時(shí)間的增加逐漸減小,并趨于恒定值。

［1］ TOM H G. Salt Water Intrusion Analysis Panama Canal Locks Future Situation: Post-Panamax Locks［R］. Delft, Netherlands: Delft University of Technology,2004.

［2］ ABRAHAM G, VAN DER BURGH P, DE VOS P. Pneumatic Barriers to Reduce Salt Intrusion Through Locks ［M］. The Hague, Netherlands: Government Publishing Office, 1973.

［3］ 周華興.海水以異重流方式入侵船閘時(shí)鹽量計(jì)算公式的初步探討［J］.海岸工程,1986,5（2): 24-31.

［4］ 竺存宏,李永振.船閘異重流交換系數(shù)的計(jì)算［J］.水運(yùn)工程,1986,（8): 33-37.

［5］ 周華興,孫玉萍.船舶過(guò)“海船閘”時(shí)鹽水入侵量的分析與計(jì)算［J］.海岸工程,1987,6（2): 48-56.

（編輯:劉運(yùn)飛)

本刊電子信箱:cjkb@ 163.net

Three-dimensional Numerical Simulation on Seawater Intrusion into Freshwater Canal Through Ship Lock

JIANG Zhi-bing,LU Hong, YANG Qing-yuan
（Hydraulics Department, Yangtze River Scientific Research Institute, Wuhan 430010, China)

Abstract:In order to study the law of salt transportation in the operation of ship lock and its influence on salinity in freshwater area, we set up a three-dimensional numerical model for k-ε two-phase mixed flow. In the model, control equation groups are dispersed by finite volume method, velocity-pressure coupling by SIMPLEC algorithm, time term by single-order implicit expression, and convection term by second-order upwind scheme, with computing area divided by hexahedral mesh. This model is verified by comparison with prototype experimental results of sea water intrusion into a ship lock, and simulated data are well in agreement with measured data of salinity in lock chamber. The simulation results about the process of seawater intrusion into freshwater zone show that 1) distribution of salinity in upstream waterway can be divided to 3 sections: density current section, transition section and diffusion section, whose salt transporting speed reduces in order；2) salinity periodically changes at each section and it gradually tends to be at a dynamic balance after ship lock operates for a period；3) as for a given cycle, saline intrusion quantity is big at early stage of ship lock operation, then, it gradually reduces and tends to be stable. The research results offer scientific reference for the design and operation of ship lock.

Key words:ship lock；seawater intrusion；numerical model；distribution of salinity；saline intrusion quantity

作者簡(jiǎn)介：姜治兵（1978-),男,湖北天門(mén)人,高級(jí)工程師,主要從事水工水力學(xué)與計(jì)算水力學(xué)方面的研究,（電話)027-82927197（電子信箱) zhibingj@163.com。

基金項(xiàng)目：國(guó)家自然科學(xué)基金項(xiàng)目（51209007)；國(guó)家科技支撐計(jì)劃課題資助項(xiàng)目（2012BAK10B04)

收稿日期：2015-10-12;修回日期：2015-11-03

doi:10.11988/ ckyyb.20150847

中圖分類號(hào)：TV135.4

文獻(xiàn)標(biāo)志碼：A

文章編號(hào)：1001-5485(2016)02-0052-05