葉明波，伍志元，楊明遠(yuǎn)，蔣昌波，3，陳 杰，3

北江是珠江水系第二大河流。為了適應(yīng)內(nèi)河水運(yùn)發(fā)展的需要，北江航道擬提升至III級航道標(biāo)準(zhǔn)，計(jì)劃實(shí)施北江（烏石至三水河口）217km河段III級航道整治工程。白石窯樞紐的一、二線船閘是在樞紐整體已建成條件下進(jìn)行的改擴(kuò)建工程，且均布置于同岸，存在著一、二線船閘水流條件和船舶進(jìn)出閘的相互干擾問題；白石窯壩下游的脫水段、水庫回水變動區(qū)河段和庫尾淺灘段需要同時進(jìn)行航道整治，對白石窯一、二線船閘的通航水位等亦將產(chǎn)生重要影響。船閘引航道口門區(qū)水流結(jié)構(gòu)是航電樞紐總平面布置方案論證及優(yōu)化的一個重要參考因素。按照內(nèi)河通航標(biāo)準(zhǔn)［1］要求，對Ⅰ～Ⅳ級船閘，平行于航線的縱向流速應(yīng)小于2.0m／s，垂直于航線的橫向流速應(yīng)小于0.3m／s，回流流速應(yīng)小于0.4m／s。

現(xiàn)階段絕大部分樞紐采用水工物理模型試驗(yàn)來研究船閘引航道口門區(qū)及連接段的通航水流條件，并以此為主要判別標(biāo)準(zhǔn)進(jìn)行樞紐的總平面布置論證。鄭寶友［2－3］對三峽船閘下游口門區(qū)和西江那吉航電樞紐上游口門區(qū)的水流條件的模型試驗(yàn)結(jié)果進(jìn)行了分析，闡述口門區(qū)的水流特性、回流機(jī)理以及影響回流強(qiáng)度的各種因素，并提出了改善口門區(qū)水流條件的有效措施；陳作強(qiáng)［4］通過物理模型、船舶模型和現(xiàn)場水流測驗(yàn)，研究了連接段布置形式對口門區(qū)水流條件的影響；李君濤［5］探討了導(dǎo)流墩削弱口門區(qū)內(nèi)斜流和回流、改善通航水流條件的機(jī)理。

由于物理模型受到模型場地、比尺效應(yīng)、研究周期及研究費(fèi)用等諸多限制，利用數(shù)學(xué)模型計(jì)算是解決這些問題的有效途徑。陳輝［6－7］等學(xué)者采用平面二維數(shù)值方法，對口門區(qū)水流條件進(jìn)行了分析，計(jì)算流速多為水深平均流速，但由于口門區(qū)水流運(yùn)動的復(fù)雜性，水流運(yùn)動三維特性明顯；且內(nèi)河通航標(biāo)準(zhǔn)［1］對于縱向流速、橫向流速及回流流速的要求均為直接影響船舶航行的表層流速，平面二維模型采用的水深平均流速不能準(zhǔn)確反映口門區(qū)不同水深位置的流速分布規(guī)律，故發(fā)展適用于復(fù)雜河段、復(fù)雜工況下的口門區(qū)三維水流數(shù)學(xué)模型，已成為工程應(yīng)用研究的一個重要趨勢。馮小香［8］基于平面曲線坐標(biāo)系、垂向σ坐標(biāo)系，建立了三維水流數(shù)學(xué)模型，研究了彎曲河段三維水流結(jié)構(gòu)特征和船閘引航道口門區(qū)的通航水流條件。

作者擬結(jié)合工程整治方案，應(yīng)用三維數(shù)學(xué)模型對白石窯樞紐下游引航道口門區(qū)水流條件進(jìn)行計(jì)算分析，開展白石窯樞紐一、二線船閘總體平面布置，對船閘通航水流條件的影響進(jìn)行研究。

1 三維數(shù)學(xué)模型計(jì)算方法

1.1 控制方程

在流體不可壓縮、Boussinesq假定下，采用笛卡爾坐標(biāo)系，給出基于RANS方程的口門區(qū)三維水動力控制方程組?？刂品匠探M由連續(xù)方程、動量方程和紊流方程組成。

連續(xù)性方程為：

X方向的水平動量方程為：

Y方向的水平動量方程為：

式中：t為時間；x，y，z均為笛卡爾坐標(biāo)系中的坐標(biāo)；d為靜止水深；h＝η＋d為總水深；η為水位；u，v，w 分別為流速在x，y，z方向上的分量；vt為垂向紊動流速。

水平力可以用梯度壓力關(guān)系表示為：

式中：A為水平紊動粘度。

u，v及w的表面和底層邊界條件為：

式中：（τsx，τsy）和 （τbx，τby）分別表示x，y方向上表面和底層的應(yīng)力分量。

本次計(jì)算紊動模型垂向采用k－ε模型，紊動粘度采用參數(shù)k和ε確定［9］：

水平方向采用Smagorinsky［10］提出的模型：

1.2 數(shù)值方法

采用分層網(wǎng)格處理，水平方向采用非結(jié)構(gòu)網(wǎng)格，垂直方向采用結(jié)構(gòu)網(wǎng)格。通過有限體積法對方程進(jìn)行空間離散，將計(jì)算域劃分成一系列互不交疊的控制單元，針對單個單元，運(yùn)用高斯定理可以將矢量方程改寫成積分形式。根據(jù)Roe格式［11］的近似黎曼求解方法，計(jì)算垂向界面的對流流量，采用線性梯度重建方法可以得到二階精度，其中平均梯度計(jì)算采用Jawahar和Kamath［12］的方法。為避免數(shù)值振蕩，采用二階TVD緩坡限制器，水平界面的對流流量使用一階迎風(fēng)格式求解。

1.3 干濕邊界處理

模型采用凍結(jié)法，設(shè)置3個特征水深，即干水深0.005m、淹沒水深0.01m和濕水深0.02m。當(dāng)某一單元的水深小于濕水深0.02m而大于干水深0.005m時，該單元只考慮質(zhì)量守恒，不考慮動量守恒；當(dāng)水深小于干水深0.005m時，單元被凍結(jié)，不參與計(jì)算；淹沒深度用來檢測單元是否已經(jīng)被淹沒。

2 模型的建立與驗(yàn)證

2.1 模型的建立

建立白石窯樞紐下游引航道口門區(qū)三維水流數(shù)學(xué)模型。該模型進(jìn)口位于白石窯樞紐，出口位于樞紐下游2.3km。計(jì)算模型平面方向采用三角形非結(jié)構(gòu)化網(wǎng)格系統(tǒng)，垂向采用σ坐標(biāo)，分為10層，最大網(wǎng)格長度約為50m。為適應(yīng)口門區(qū)、引航道及丁壩附近復(fù)雜水流情況，模型對樞紐、引航道和丁壩附近的網(wǎng)格進(jìn)行局部加密處理，最小網(wǎng)格長度約為3m，計(jì)算域內(nèi)共布置三角形網(wǎng)格節(jié)點(diǎn)16 725個，三角形單元32 427個。二線船閘修建前的河道計(jì)算范圍地形和數(shù)學(xué)模型計(jì)算網(wǎng)格如圖1所示。模型進(jìn)口為樞紐泄水閘和電站泄水口，計(jì)算過程中，模型可根據(jù)實(shí)際泄水情況調(diào)整入水位置，將入水口設(shè)置成泄水閘泄水、電站泄水或按流量比分別泄水，能較好模擬實(shí)際的水流過程。河道中丁壩及導(dǎo)流堤的水流、固體交界面上，采用固壁非穿越及無滑移邊界條件，邊界上的法向、切向流速為零。

圖1 數(shù)學(xué)模型計(jì)算地形及網(wǎng)格Fig.1 Topography and grid of numerical model

2.2 模型的驗(yàn)證

本次白石窯樞紐船閘口門區(qū)及連接段水流結(jié)構(gòu)特性數(shù)學(xué)模型計(jì)算采用的地形資料是由廣東省航道局提供的地形測圖，模型采用500，1 000和4 990m3／s三級流量的水面線和流速數(shù)據(jù)進(jìn)行驗(yàn)證。

2.2.1 水面線驗(yàn)證

2012年重慶西南水運(yùn)工程科學(xué)研究所對北江（烏石至三水河口）航道整治工程白石窯水利樞紐壩下河段進(jìn)行了水流泥沙物理模型試驗(yàn)研究。本次數(shù)學(xué)模型計(jì)算是根據(jù)其物理模型試驗(yàn)結(jié)果進(jìn)行糙率選取和模型驗(yàn)證的。根據(jù)物理模型試驗(yàn)在數(shù)學(xué)模型計(jì)算范圍內(nèi)的3個水尺點(diǎn)的水位觀測結(jié)果對數(shù)學(xué)模型進(jìn)行水面線驗(yàn)證，驗(yàn)證結(jié)果如圖2所示。驗(yàn)證結(jié)果表明：數(shù)學(xué)模型與物理模型的觀測水面線一致，該數(shù)學(xué)模型能夠達(dá)到本研究的計(jì)算要求，具有較高的精度。

2.2.2 流速驗(yàn)證

圖2 水面線驗(yàn)證結(jié)果Fig.2 Verification results of water surface profile

流速分布驗(yàn)證同樣利用2012年重慶西南水運(yùn)工程科學(xué)研究所對白石窯水利樞紐壩下河段進(jìn)行的水流泥沙物理模型試驗(yàn)結(jié)果。試驗(yàn)進(jìn)行了天然狀態(tài)下 （Q＝98，500，1 000，2 000，3 000，4 990和6 740m3／s共7級流量）的流速、水位和流態(tài)觀測，根據(jù)流量大小、電站和閘門工作與否選取Q＝500，1 000和4 990m3／s三種典型工況進(jìn)行流速分布驗(yàn)證，結(jié)果如圖3所示。從圖3中可以看出，除少數(shù)測點(diǎn)略有差異外，各斷面的流速分布規(guī)律趨勢與物理模型試驗(yàn)一致，表明所建立的三維水流數(shù)學(xué)模型中各參數(shù)取值合理，較真實(shí)地反映了研究河段的水流特性，可用于三維水流分析計(jì)算。

圖3 斷面流速驗(yàn)證結(jié)果Fig.3 Verification results of velocity

3 計(jì)算結(jié)果與分析

樞紐二線船閘修建后，對船閘下游引航道口門區(qū)通航水流條件數(shù)學(xué)模型進(jìn)行研究。分別對98，500，1 000，2 000，3 000，4 990和6 740m3／s等7種流量下的9種工況進(jìn)行了計(jì)算分析。以船閘下引航道平面布置采用重慶西南水運(yùn)工程科學(xué)研究所提出的設(shè)計(jì)方案二、流量為1 000m3／s為例，對船閘口門區(qū)水流條件進(jìn)行了討論，在下游引航道口門區(qū)附近布置6個斷面共36個表面流速測量點(diǎn)，分別測量表面縱向、橫向、回流流速以及計(jì)算區(qū)域三維流速分布和口門區(qū)流場分布。

3.1 表面流速分布規(guī)律

當(dāng)流量為1 000m3／s時，泄水閘全部關(guān)閉，1?！?＃機(jī)組各引用200m3／s。由36個表面流速測量點(diǎn)計(jì)算得出口門區(qū)表層縱向流速、橫向流速及回流流速，引航道口門區(qū)最大縱向流速為0.935m／s，小 于 2.0m／s；最 大 橫 向 流 速 為0.173m／s，小 于 0.3m／s；最 大 回 流 流 速 為0.095m／s，小于0.4m／s。因此，下游引航道口門區(qū)能滿足通航水流條件要求。

各級流量下，口門區(qū)表面最大縱向流速、最大橫向流速及最大回流流速的計(jì)算結(jié)果見表1。從表1中可以看出，在兩年一遇（Q＝4 990m3／s）和五年一遇 （Q＝6 740m3／s）洪水下泄時，引航道口門區(qū)出現(xiàn)橫向流速稍大于0.3m／s的情況，應(yīng)注意行船安全，必要情況下考慮限制通航或者禁止通航；其余各流量條件下，引航道口門區(qū)縱向流速均小于2.0m／s、橫向流速均小于0.3m／s、回流流速均小于0.4m／s，引航道口門區(qū)滿足通航水流條件要求。

3.2 流場特征

通過計(jì)算，得到了引航道口門區(qū)流場分布。以流量1 000m3／s時為例，給出的流速矢量分布如圖4所示。從圖4中可以看出，水體表層流態(tài)與底層一致，表層流速大于底層流速。下泄水流在坑口咀附近受礁石洲影響分成左、右兩汊，左汊為主流，流向與河流左岸成45°左右夾角下行，至下游引航道出口處，主流逐漸向右擴(kuò)散，在口門區(qū)及連接段形成與航道中心線存在一定夾角的斜流，關(guān)閉，電站泄水區(qū)流速較大，泄水閘泄水區(qū)流速較?。欢问饔脤?dǎo)致河道變窄，導(dǎo)致丁壩壩頭附近流速較大；引航道口門區(qū)流速較小，表明導(dǎo)航墻作用明顯，能較好地滿足船舶通航的需求。與平面二維模型采用的水深平均流速相比［13］，三維模型可以給出水體各層的流速分布，能更加準(zhǔn)確地反映引航道口門區(qū)流速對船舶航行的影響。從圖5中還可以看出，表面流速大于中層流速和底部流速，即相對二維模型而言，三維模型的表面流速較水深平均流速危險(xiǎn)，故船閘口門區(qū)水流條件與通航條件分析采用三維模型更為安全。斜流效應(yīng)明顯，斜流夾角為10°～30°，但橫向流速較小，不影響口門區(qū)通航水流條件；同時在口門區(qū)右側(cè)形成一個逆時針方向的回流區(qū)，回流強(qiáng)度較弱，引航道出口下150m斷面下游河道水流全斷面順流，引航道內(nèi)僅出口處為弱回流。

表1 各流量下，口門區(qū)表面最大縱向、橫向及回流流速Table 1 Surface maximum longitudinal，lateral and backward flow velocity in each discharge

圖4 口門區(qū)流場分布Fig.4 Velocity field on entrance area

3.3 工程前、后水流條件比較

流量為1 000m3／s情況下二線船閘修建前、后樞紐下游的三維流速分布云圖如圖5所示，分別給出了上、中及下3層的流速分布規(guī)律。從圖5中可以看出，由于該工況下采用電站發(fā)電、泄水閘

圖5 樞紐下游三維流速分布云圖Fig.5 Cloud picture of 3Dvelocity field

3.4 渦量分布特征

渦量是流體的基本物理量，用來描寫水體的旋渦運(yùn)動，水體中渦量決定了水體輸運(yùn)特性和能量耗散的程度。口門區(qū)表層與底層的渦量分布如圖6所示。從水平方向可以看出，渦量峰值均出現(xiàn)在導(dǎo)航墻墻頭附近，峰值出現(xiàn)在口門區(qū)外側(cè)，對航道影響較小；從垂直方向可以看出，表層水體的渦量峰值比底層水體渦量峰值大，且峰值分布區(qū)域更廣。總體來講，口門區(qū)渦量較小，表明水體旋渦運(yùn)動微弱，水流穩(wěn)定，利于通航。

圖6 口門區(qū)渦量等值線Fig.6 Contour of vorticity magnitude

4 結(jié)論

利用有限體積法，對分層網(wǎng)格進(jìn)行空間離散，建立了口門區(qū)三維水流數(shù)學(xué)模型。對白石窯樞紐下引航道口門區(qū)的流速分布情況、流場特征及渦量特征進(jìn)行了數(shù)值分析，討論了通航水流條件，得到結(jié)論：

1）本模型水平采用三角形無結(jié)構(gòu)網(wǎng)格，垂向采用σ坐標(biāo)分層，較好地?cái)M合了工程區(qū)復(fù)雜的地形特征。通過驗(yàn)證，模型能用于樞紐引航道口門區(qū)水流條件的計(jì)算分析。計(jì)算結(jié)果較二維模型更為準(zhǔn)確，可以更科學(xué)合理地判斷引航道口門區(qū)通航水流條件。

2）下泄主流在口門區(qū)及連接段形成與航道中心線存在斜流效應(yīng)，在兩年一遇（Q＝4 990m3／s）和五年一遇（Q＝6 740m3／s）洪水下泄時，口門區(qū)局部出現(xiàn)橫向流速稍大于0.3m／s的情況，應(yīng)注意行船安全，必要情況下考慮限制通航或者禁止通航；其余各流量條件下，橫向流速均小于0.3m／s；引航道口門區(qū)存在水體回流，但回流強(qiáng)度較小，回流流速均小于0.4m／s，引航道口門區(qū)滿足通航水流條件要求。

3）口門區(qū)水體表層流態(tài)與底層流態(tài)一致，表層流速大于底層流速；渦量峰值出現(xiàn)在導(dǎo)航墻墻頭附近、口門區(qū)外側(cè)，表層水體渦量峰值比底層水體渦量峰值大，且峰值分布區(qū)域更廣，口門區(qū)內(nèi)渦量較小，表明水體旋渦運(yùn)動微弱，利于船舶下行。

4）很多情況下，引航道中水流具有非恒定性，雖然本研究中采用的方程為非恒定流方程，但現(xiàn)階段計(jì)算的是恒定流量級，因而對水流條件分析具有一定的局限性?？紤]水流非恒定特性研究，引航道口門區(qū)通航的水流條件是今后研究的方向。

（

）：

［1］ 中華人民共和國交通部.GBJ50139－2004，內(nèi)河通航標(biāo)準(zhǔn)［S］.北京：人民交通出版社，2004.（Ministry of Communications of the People’s Republic of China.GBJ50139－2004，Navigation standard of inland waterway［S］.Beijing：China Communications Press，2004.（in Chinese））

［2］ 鄭寶友.三峽船閘下游口門區(qū)流態(tài)及其改善措施的研究［J］.水道港口，1990（3）：21－28.（ZHENG Baoyou.Study on the flow regime at the downstrem entrance area of the Three Gorge Lock and its improving measures［J］.Journal of Waterway and Harbor，1990（3）：21－28.（in Chinese））

［3］ 鄭寶友，陳波.船閘上游口門區(qū)連接段為彎道的通航水流條件［J］.水道港口，2005（2）：99－102，108.（ZHENG Bao－you，CHEN Bo.Test on navigation flow condition of curving connect section of upstream lock entrance［J］.Journal of Waterway and Harbor，2005（2）：99－102，108.（in Chinese））

［4］ 陳作強(qiáng).通航建筑物口門區(qū)及連接段通航水流條件研究［D］.成都：四川大學(xué)，2006.（CHEN Zuo－qiang.Research on flow conditions for navigation at the entrance area and the transitional reach of the navigation structures［D］.Chengdu：Sichuan University，2006.（in Chinese））

［5］ 李君濤，張公略，馮小香.導(dǎo)流墩改善口門區(qū)水流條件機(jī)理研究［J］.中國港灣建設(shè)，2011（2）：1－3.（LI Juntao，ZHANG Gong－lve，F(xiàn)ENG Xiao－xiang.Improvement mechanism of flow conditions at the entrance of approach channel by diversion pier［J］.China Harbour Engineering，2011（2）：1－3.（in Chinese））

［6］ 陳輝，劉志雄，江耀祖.引航道通航水流條件數(shù)值模擬［J］.水利水運(yùn)工程學(xué)報(bào)，2012（4）：13－18.（CHEN Hui，LIU Zhi－xiong，JIANG Yao－zu.Numerical simulation and analysis of navigation conditions in a shiplock approach channel［J］.Hydro－Science and Engineering，2012（4）：13－18.（in Chinese））

［7］ 黃明海，黃國兵.三峽工程引航道通航水流數(shù)值模擬應(yīng)用研究［J］.湖北水力發(fā)電，2007（3）：36－39.（HUANG Ming－h(huán)ai，HUANG Guo－bing.Numerical simulation on navigable flow in approach channel of Three Gorges Project［J］.Hubei Water Power，2007（3）：36－39.（in Chinese））

［8］ 馮小香，李丹勛，張明.樞紐船閘引航道口門區(qū)三維水流數(shù)值模擬應(yīng)用研究［J］.水運(yùn)工程，2012（1）：122－126.（FENG Xiao－xiang，LI Dan－xun，ZHANG Ming.Application research on 3Dflow numerical model at the entrance area of ship lock approach channel［J］.Port ＆ Waterway Engineering，2012（1）：122－126.（in Chinese））

［9］ Li Y S，Zhan J M.An efficient three－dimensional semi－implicit finite element scheme for simulation of free surface flows［J］.International Journal for Numerical Methods in Fluids，1993，16：187－198.

［10］ Smagorinsky J S.General circulation experiments with the primitive equations［J］.Monthly Weather Review，1963，91：99－164.

［11］ Brufau P，Garcia P，Vazquez－cendon M E.Zeromass error using unsteady wetting－drying conditions in shallow flow sover dry irregular topography［J］.International Journal for Numerical Methods in Fluids，2004，45（10）：1047－1082.

［12］ Jawahar P，Kamath H.A high－resolution procedure for euler and navier－stokes computations on unstructured grids［J］.Journal of Computational Physics，2000，164（1）：165－203.

［13］ 楊明遠(yuǎn)，楊武，郁達(dá)，等.北江千噸級航道整治研究：Ⅱ二維數(shù)學(xué)模型與整治河段通航水流條件分析［J］.交通科學(xué)與工程，2013，29（4）：66－72.（YANG Ming－yuan，YANG Wu，YU Da，et al.Study on the Beijiang River 1 000twaterway regulation：Ⅱtwodimensional numerical model development and study on the navigation flow conditions in the regulation river［J］.Journal of Transportation Science and Engineering，2013，29（4）：66－72.（in Chinese））