李春良，朱 峰

(1.山東省交通規(guī)劃設計院有限公司，濟南 250031;2.中交第四航務工程勘察設計院有限公司，廣州 510230)

由于水路運輸?shù)谋憬菪院徒?jīng)濟性，內(nèi)河航運是內(nèi)陸地區(qū)大宗貨物運輸?shù)闹饕緩?。?nèi)河航道的水深條件與外海相比往往較淺[1]，而水流流速卻普遍較大，因此與遠洋船舶相比，內(nèi)河船舶更加注重船舶的可操縱性[2]。

圖1 典型大尺度床面水流Fig.1 Typical large scale bed surface flow

我國地勢東低西高，西部地區(qū)是我國主要水系的發(fā)源地，隨著水流從西部狹窄山谷到東部開闊平原的過渡變化，水流流速和水位變化也呈現(xiàn)出西側湍急雜亂，東側均勻平順的規(guī)律，因此中西部地區(qū)的惡劣水流條件給內(nèi)河航運帶來了較大挑戰(zhàn)[3]。

與水流條件變化相對應的，從東到西的過程中，水流對泥沙顆粒的分選效應逐漸顯著，中西部河流的河床組成顆粒直徑大于東部的河床組成顆粒。西部山區(qū)河流的河床主要由一定尺度的碎石和卵石組成，其中大型的塊石等效直徑可以達到數(shù)米，在枯水期常常露出水面[4]。

圖2 河工試驗人工加糙Fig.2 Artificial roughening in river experiment

眾所周知，當河床足夠光滑時，水流一般呈現(xiàn)自然平順狀態(tài)，河床糙率對水流狀態(tài)的影響完全可以采用明渠均勻流理論進行描述；而當河床的組成顆粒粒徑較大時，河床中的塊石會顯著影響水體的運動規(guī)律，從而導致表面水流狀態(tài)產(chǎn)生紊亂[5]，水流運動規(guī)律無法再用明渠均勻流理論進行描述。常見的大尺度粗糙河床表面形態(tài)如圖1所示。

粗糙床面對水流狀態(tài)的影響不僅體現(xiàn)在自然河流中。在物理模型試驗過程中，模型比尺選定為1:100時，為了準確模擬河床對水流的摩擦阻力，需要采用大顆粒對河床床面進行加糙[6]，導致加糙顆粒尺寸與試驗過程中的水深幾乎一致，如此加糙帶來的局部水流狀態(tài)不合理對實驗結果產(chǎn)生了不利影響[7]，加糙結果如圖2所示。

因此充分認識大顆粒粗糙床面對水流狀態(tài)的不利影響，是掌握河流航道演變規(guī)律所需要面對的關鍵性問題，本文基于模型試驗的方法研究了大顆粒粗糙床面對水流自由表面波動、流速分布及沿程阻力系數(shù)的影響，并基于計算流體力學分析和多元統(tǒng)計回歸相結合的方法，對試驗結果進行驗證和擴充，對山區(qū)河流的船舶航行以及航道整治工作具有顯著的參考價值[8]。

1 模型試驗

圖3 水流試驗現(xiàn)場照片F(xiàn)ig.3 Field photo of water flow test

模型試驗所采用的水槽寬度為0.5 m，總長度為28 m，高度達到0.5 m。通過自動化系統(tǒng)可以控制試驗水槽坡度在±7‰之間按照實際需求調(diào)整。水槽的最大循環(huán)造流量為120 L/s。在開始正式試驗前需要對實驗區(qū)域的水流狀態(tài)進行驗證，試驗布置方案為：從水槽中部向出流方向10 m范圍內(nèi)鋪設加糙塊體，試驗區(qū)段0～5 m范圍為水流過渡段，確保水流由進入試驗段之前的平順狀態(tài)完全發(fā)展為粗糙床面流態(tài)；試驗區(qū)段5～8 m范圍用于試驗數(shù)據(jù)采集；試驗區(qū)段8～10 m范圍用于消除出口段對粗糙床面流的影響。整個試驗過程中采用電磁流量計測定水流的流量，并采用超聲多普勒流速儀測量關心位置處的水流流速[9]。試驗現(xiàn)場片如圖3所示。

大直徑加糙顆粒選用邊長50 mm的混凝土立方體試件進行模擬，加糙塊體采用錯位間隔布置方案。加糙塊體沿著橫向的布置間距為50 mm，沿水槽縱向的布置間距根據(jù)不同試驗方案分別設置為25 mm、50 mm、100 mm、150 mm和200 mm。

由于水槽試驗沿著水槽縱向具有對稱性，且塊體布置具有一定周期性，因此水流流速測點僅布置在水槽的一側，根據(jù)塊體間距不同布置2～8組測速斷面，每個測速斷面包含5組垂向測線，每個測線在深度方向每隔1 cm布置一個測速點。典型的塊體布置和測速點布置如圖4所示。

4-a阻流塊體擺放密度e=0.114-b阻流塊體擺放密度e=0.42圖4 加糙塊體布置方案(單位:mm)Fig.4Layoutplanofrougheningblock

5-a 低速水流狀態(tài) 5-b 高速水流狀態(tài)圖5 不同流速條件下水流狀態(tài)試驗結果Fig.5 Test results of flow state under different velocity conditions

分別針對不同的塊體擺放密度、不同的底坡坡度和不同水深條件的組合工況開展系列試驗，并基于試驗結果采集得到了自由表面波動、流速分布和沿程水頭損失等試驗結果。圖5展示了試驗過程中的典型自由表面波動結果。

從高、低流速條件下的典型試驗結果可以看到：在水流低速運動狀態(tài)下，試驗水面相對平穩(wěn)，水面線沿程均勻分布；在水流高速運動狀態(tài)下，加糙塊體對水流的阻滯作用逐漸加強，使得塊體前側水位發(fā)生壅高；加糙塊后方由于水流有效通過面積增加，因此后側局部水位降低，并由此引起液面發(fā)生一定程度的波動，波動在水槽內(nèi)傳播過程中發(fā)生不斷的反射疊加，最終使得水流液面出現(xiàn)不規(guī)律形態(tài)。

2 數(shù)值計算

2.1 計算模型

圖6 CFD數(shù)值計算模型Fig.6 CFD numerical model

盡管物理模型試驗具有天然的真實性，但是試驗結果的重復性往往難以保障[10]，且模型試驗的成本通常較高，試驗周期往往也較長，因此數(shù)值仿真技術一般用于對模型試驗的驗證和補充[11]。本文采用CFD軟件建立與本試驗相一致的計算模型，對試驗中所涉及的不同加糙條件下的水流運動特性進行模擬分析。

數(shù)值計算模型與試驗模型的尺寸相一致，為確保水流狀態(tài)得到充分發(fā)展，數(shù)值計算中設置了與試驗類似的水流過渡段。計算模型中設置模擬計算的總時長為100 s。為了確保數(shù)值計算的精度，計算單元尺度設置為5 mm，并啟用k-e紊流模型計算水流運動中的粘性效應。數(shù)值計算模型如圖6所示。

2.2 基本方程

質(zhì)量守恒方程和動量守恒方程是計算流體力學的基本控制方程分，根據(jù)計算單元中的水體質(zhì)量的保守特性，質(zhì)量守恒方程可以表達為如下形式

(1)

式中：VF是單元中水的體積分數(shù)；RSOR是質(zhì)量源項；Ax、Ay、Az是單元體中水的流動面積；u、v、w是單元中水在三個方向的流動速度。

根據(jù)水體所受外力沖量和自身動量改變之間的關系，水體的動量守恒方程可以表達為如下形式

(2)

式中：Gx、Gy、Gz是作用在流體上三個方向的體力；fx、fy、fz是作用在三個方向上的面力；bx、by、bz是水流在多孔介質(zhì)中三個運動方向上的水頭損失。

水流的紊動是一種涉及多尺度、多學科的復雜運動狀態(tài)，目前采用質(zhì)量守恒和動量守恒方程直接計算水流的紊動特性面臨極大的計算量和極大的收斂困難。因此工程上往往用紊動動能來表達水流的紊動狀態(tài)

(3)

紊動動能在水流中的產(chǎn)生、輸運和擴散可以基于如下方程進行描述

(4)

式中：PT是由于速度梯度所產(chǎn)生的紊流動能；GT是由于浮力所產(chǎn)生的紊流動能；εT是紊流動能消散項。

2.3 計算結果

通過數(shù)值計算可以準確獲取水流的運動情況。以塊體加糙密度為0.21工況為例，圖7展示了不同流速下的水流自由液面狀態(tài)。從圖中可以看到：低速水流作用下，自由表面形態(tài)呈現(xiàn)出與塊體擺放形態(tài)相適應的橫向條帶狀變形，液位在塊體的前后兩側出現(xiàn)規(guī)律性的高低變化；隨著來流量的不斷增加，水體質(zhì)點的慣性效應逐漸加強，原自由液面的橫向條帶狀分布受到破壞，逐漸過渡為具有短峰特性的自由液面波動。

7-a e=0.21,H=0.11m,Q=2.5L/s7-b e=0.21,H=0.11m,Q=22.5L/s圖7 不同流速下水流自由表面結果Fig.7Freesurfaceresultsofflowatdifferentvelocities

8-a t=3.4,e=0.21,H=0.107m,Q=2.5L/s8-b t=11.0,e=0.21,H=0.107m,Q=2.5L/s

8-c t=40,e=0.21,H=0.107m,Q=22.5L/s圖8 不同流速下的平面流速分布Fig.8Planevelocitydistributionunderdifferentvelocity

圖9 模型試驗與數(shù)值計算結果對比Fig.9 Comparison of model test and numerical calculation results

水流的縱向流速分布是研究塊體阻流效應的重要參考依據(jù)。因此文章重點研究了不同塊體擺放密度下水流的縱向流速分布。圖8展示了在塊體擺放密度為0.21條件下的縱向流速平面分布。從結果等值線云圖中可以看到：在入口處來流量較小條件下，縱向流速分布與塊體的橫向布置具有較大的相關性，高流速區(qū)出現(xiàn)在塊體頂部；隨著來流量的逐漸加大，水質(zhì)點的慣性力與粘性力相比大幅度增加，使得縱向流速的分布特性逐漸轉(zhuǎn)變?yōu)橐钥v向為主的分布特性。在水槽邊緣部位，由于水槽橫向?qū)挾扔邢逕o法形成完整的橫向周期性，從而導致在水槽邊緣區(qū)域的塊體阻流效應較弱，高速水流在慣性作用驅(qū)使下，主要從阻力較低的區(qū)域通過，因此縱向流速在模型邊界處有所增加。

通過基于本文所建立數(shù)值模型計算得到的不同組合工況下的沿程水頭損失結果與物理模型試驗結果對比如圖9所示。對比試驗結果可以發(fā)現(xiàn)：數(shù)值模型的計算結果與物理模型試驗結果之間的共線性達到0.94，說明數(shù)值計算和模型試驗結果是相匹配的；除此以外，計算結果和試驗結果之間的相關度高達0.99，由此可以判斷：數(shù)值計算的結果是合理的，可以充分準確表達水流與塊體相互作用的實際情況。

3 回歸建模分析

文章列舉了不同塊體擺放密度、不同底坡、不同水深組合條件下20種工況的試驗結果如表1所示，為了使得研究成果可以為實際工程提供參考依據(jù)，文章基于試驗結果，開展經(jīng)驗公式修正分析。

針對天然河床在水流運動過程中的水頭損失問題，Bathurth等人通過對40多個河流區(qū)段進行觀察統(tǒng)計，提出了縱向沿程水頭損失的經(jīng)驗公式如下

(5)

基于上述經(jīng)驗公式，對本文試驗所得到的數(shù)據(jù)進行驗證后可以得到如圖10所示結果。從結果對比中可以看到：盡管試驗結果和經(jīng)驗公式結果之間的相關度很高，但是兩者間的共線性較差(僅達到0.1)。該結果是由于經(jīng)驗公式僅僅考慮了粗糙顆粒的尺度和水深這兩個因素，尚無法全面客觀的描述水流的運動規(guī)律。

因此，文章進一步結合試驗結果和上述經(jīng)驗公式的基本形式，增加考慮了河床底坡和塊體密度兩個影響因素，對經(jīng)驗公式進行了修正

(6)

式中：A=λB，B=a1i+a2，λ=a3exp(a4·e)。

針對上述修正經(jīng)驗公式形式，本文基于非線性多元回歸方法，確定了經(jīng)驗公式中的四個系數(shù)

采用經(jīng)過修正的經(jīng)驗公式再次對試驗和計算結果進行驗證，結果如圖11所示。增加考慮了水槽縱向坡度和底坡加糙密度這兩個參數(shù)后，經(jīng)驗公式與試驗結果符合度有所改善。

圖10 Bathurth公式結果對比(修正前)Fig.10ComparisonoftheresultsofBathurthformula(BeforeModification)圖11 Bathurth公式結果對比(修正后)Fig.11Comparisonoftheresultsofbathurthformula(AfterModification)

4 結論

文章基于物理模型試驗、計算流體力學模擬和數(shù)據(jù)建模分析的方法，研究了人工加糙河床上的水流流動特征，研究所得到的主要結論如下：

(1)通過一系列工況組合下的水槽試驗，可以準確記錄并得到水流流動過程中，自由表面變化、縱向流速分布和沿程水頭損失等參數(shù)；(2)對文章所涉及試驗的流體力學建模與計算可以發(fā)現(xiàn)：通過設置可靠的計算網(wǎng)格尺寸、施加準確的進出口和壁面條件以及選用合理的紊流模型后，通過計算可以獲得與試驗結果相一致的水流自由表面形態(tài)、水流流速分布以及水頭損失系數(shù)等；(3)文章通過增加考慮水槽坡度和塊體擺放間距兩個影響因素對傳統(tǒng)經(jīng)驗公式開展修正，并采用非線性多元回歸法對修正公式進行參數(shù)估計，使得修正后的公式結果更符合試驗結果。