許虎，吳文勇,2*，王振華，喬長錄，王秋良

（1.石河子大學(xué)水利建筑工程學(xué)院，新疆石河子832000；2.中國水利水電科學(xué)研究院，北京100048；3.深圳市水務(wù)規(guī)劃設(shè)計院，廣東深圳518000）

0 引言

【研究意義】水資源優(yōu)化調(diào)配是農(nóng)業(yè)節(jié)水灌溉的基礎(chǔ)，做好灌區(qū)量水工作是保障灌區(qū)合理運營的關(guān)鍵。巴歇爾槽具有量測精度高、讀數(shù)方便、不易淤積等優(yōu)點，在我國灌區(qū)應(yīng)用廣泛[1]。水流由上游流經(jīng)進口收縮段時，水位升高、流速增大形成臨界流。此時上游水位與流量形成單一的函數(shù)關(guān)系，通過量測上游水位能夠推求渠道流量變化[2]。

隨著計算流體力學(xué)（CFD）的發(fā)展，數(shù)值模擬應(yīng)用于水利工程設(shè)施的相關(guān)研究[3]。通過建立數(shù)值模型，能夠?qū)⒘黧w運動軌跡、流速、壓強、溫度等物理量清晰地呈現(xiàn)出來。相較于傳統(tǒng)的模型試驗，數(shù)值模擬具有適用性強、試驗周期短、計算結(jié)果清晰等優(yōu)點[4]。

【研究進展】近年來，國內(nèi)外學(xué)者對各式量水裝置進行模型試驗與數(shù)值模擬分析，得出不同工況下槽內(nèi)水流水力特性和試驗結(jié)果。劉英等[5]、吉慶豐等[6]對不同體型圓柱形量水柱三維水流運動進行數(shù)值計算結(jié)果表明，水力參數(shù)的實測值與模擬值具有良好的一致性，在收縮比為0.63時，最大測流誤差為4.95%，三維水流流態(tài)與試驗結(jié)果一致。孫斌等[7]、宋金妍[8]、肖苡辀等[9]對不同收縮比的翼型量水槽進行流速場、水面線、水頭損失、佛汝德數(shù)水力計算，認為收縮比為0.50～0.60 比較合理，并對機翼形量水槽翼面進行優(yōu)化。張敏等[10]基于Flow-3D 軟件對不同喉道長度的短喉道量水槽進行數(shù)值模擬選型，得出臨界水深隨喉道長度的增加而增加，短喉道量水槽較長喉道和無喉道量水槽水頭損失小，且水流流態(tài)平穩(wěn)，形式最優(yōu)。侯瑩等[11]、冉聃頡等[12]、景志芳等[13]對無喉道量水槽過渡段形式水力特性分析，建立流量公式，得出直線形過渡段形式測流精度高但水頭損失大、圓形過渡段形式測流精度低但水頭損失小，湍動能耗散主要在變比和底部區(qū)域，在測流范圍內(nèi)水頭損失小于總水頭損失的10%，試驗值與模擬值相吻合。

數(shù)值模擬對量水堰的試驗研究也同樣適用，柳雙環(huán)等[14]基于FLUENT 軟件對U 形渠道三角形剖面堰數(shù)值模擬，得出不同底坡、渠道邊坡、堰高下的水頭損失，并對不同堰高擬合流量公式。葉青青[15]使用RNGk-ε湍流模型，VOF 方法對無喉堰在自由出流和淹沒出流情況下擬合流量公式，通過數(shù)值模擬得出自由液面、流速、壓力云圖分布情況?！厩腥朦c】在實際工程建設(shè)中，巴歇爾槽與渠道上游連接段并沒有明確的要求，目前常見的過渡段形式有圓弧段連接、斜面連接和無過渡連接。不同連接段形式對局部水頭損失、水面線波動、測流精度有很大影響，在實際的輸水過程中，連接段形式的選取決定了渠道輸水效率高低與水資源配置決策?！緮M解決的關(guān)鍵問題】本文采用數(shù)值模擬方法為巴歇爾槽連接段優(yōu)化選型，為灌區(qū)巴歇爾槽的流道設(shè)計提供理論依據(jù)。

1 材料與方法

1.1 物理結(jié)構(gòu)模型

采用標準設(shè)計尺寸（表1）和SolidWorks 建模軟件建立喉道寬度為0.25m 巴歇爾槽三維物理模型。巴歇爾槽由上游收縮段、喉道段、下游擴散段組成，結(jié)構(gòu)尺寸如圖1、圖2 所示。Z1 斷面為巴歇爾槽俯視圖（圖1），巴歇爾槽上下游0.25m 處X1、X2 斷面作為水力特性計算斷面（圖2）。置于寬度為1.10 m 的矩形渠道內(nèi)，在上游收縮段和喉道段設(shè)置水位觀測點，記錄水位信息。

圖1 巴歇爾槽俯視Fig.1 Parshall flume vertical view

圖2 巴歇爾槽左視Fig.2 Parshall flume left view

上游收縮段與渠道采用4 種形式連接（圖3），分別為方案1（內(nèi)接圓弧過渡段）、方案2（外接圓弧過渡段）、方案3（直面過渡段）、方案4（無過渡段）。方案1、方案2 和方案3 進口連接段L1 長度為0.4m（3 種方案巴歇爾槽長度尺寸相同，故不在表1 中詳細描述，統(tǒng)一由方案1 表示。），方案4 無連接段L1長度為0。

圖3 巴歇爾槽上游連接段形式Fig.3 Parshall flumeupstream inlet connection section

表1 標準巴歇爾槽結(jié)構(gòu)尺寸Fig.1 StandardParshall flumestructural dimensions

1.2 數(shù)學(xué)計算模型

明渠流動涉及氣-液兩相流動，VOF（流體體積）模型在處理多種互不交融的交界面時有很好的應(yīng)用。在VOF 模型中，不同流體組分共用一套動量方程，對每一個網(wǎng)格流場計算時都記錄下各組分所占有的體積率，跟蹤自由流體表面[16]。在處理自由液面、分層流動、水壩決堤時有比較精確的效果。對于水氣二相流場，假設(shè)水和空氣具有相同的速度，在每個網(wǎng)格單元中，水和空氣的體積分數(shù)之和為1，即：

式中：aw、aa分別表示計算域中每一個控制單元內(nèi)水和空氣的體積分數(shù)（aw=1 時，表示控制體積內(nèi)全部充滿水；當aw=0 時，表示控制體積內(nèi)全部充滿空氣）。

水氣交界面使用連續(xù)性方程進行追蹤，即：

式中：ui為速度分量；xi為坐標分量。

流體運動遵循三大守恒定律，即：質(zhì)量守恒定律、能量守恒定律、動量守恒定律。明渠內(nèi)水流為不可壓縮流且黏性恒定，不涉及能量交換。對質(zhì)量守恒定律和動量守恒定律的連續(xù)性方程和動量方程進行描述[17]。

連續(xù)性方程：

式中：ρ為密度；t為時間；u為速度矢量。

動量方程：

式中：ui、uj分別為x方向流速分量（i、j=1、2、3）；p為靜壓力；t為時間；μ為黏滯系數(shù)；ρ為流體密度；g 為重力加速度。

明渠水流為充分發(fā)展的湍流，需要滿足湍流輸運方程。RNGk-ε湍流模型考慮了平均流動中的旋轉(zhuǎn)流動，能夠很好地處理明渠水流流線彎曲較大的情況。

湍動能k表示為：

耗散率ε表示為：

式中：k為湍動能（m2/s2）；ε為湍動能耗散率kg/（m2·s2）；μ為動力黏滯系數(shù)N/（s·m2）；μeff為有效動力黏滯系數(shù)N/（s·m2），（μ與μt之和，其中αk=αs=1.39），Gk為平均流速梯度引起的湍動能k的產(chǎn)生項C2ε=1.68；i、j=1、2、3）。

1.3 網(wǎng)格劃分與無關(guān)性分析

利用ANSYS18.0 軟件進行網(wǎng)格劃分，結(jié)構(gòu)性網(wǎng)格具有網(wǎng)格節(jié)點分布規(guī)則、占用計算資源較小等優(yōu)點，本次計算模型采用六面體結(jié)構(gòu)網(wǎng)格劃分。為了提高計算精度，更好地模擬近壁面水體流動特性，對壁面網(wǎng)格自下而上進行整體加密。為了確定網(wǎng)格單元個數(shù)及網(wǎng)格整體質(zhì)量，將全局網(wǎng)格最大尺寸分別設(shè)定為0.05、0.04、0.03、0.02、0.01 劃分網(wǎng)格。為了探究不同網(wǎng)格單元大小對于模擬結(jié)果的影響，對5 種尺寸的網(wǎng)格文件進行數(shù)值模擬計算。設(shè)置進口流量0.1m3/s，流速0.232m/s,渠道進口水深0.392m，求解方法與邊界條件設(shè)置按照1.4 部分內(nèi)容。分析計算結(jié)果：當網(wǎng)格最大單元尺寸小于0.03 以后，上游水深不再隨著網(wǎng)格加密發(fā)生變化（表2）。因此選用最大網(wǎng)格尺寸大小為0.03 對模型網(wǎng)格劃分（圖4），網(wǎng)格整體質(zhì)量在0.7 以上。4 種方案下網(wǎng)格單元個數(shù)分別為303452、303918、303726、301794。

表2 網(wǎng)格無關(guān)性分析Fig.2 Analysis of meshing independence

圖4 六面體網(wǎng)格劃分Fig.4 Hexahedral meshing

1.4 數(shù)值求解方法與邊界條件

1.4.1 求解方法

數(shù)值模型基于三維穩(wěn)態(tài)壓力法求解器，VOF 方法，對Y 方向設(shè)置重力-9.81m/s2。基本控制方程為N-S 方程，采用RNGk-ε湍流模型計算。采用一階迎風(fēng)格式離散，PISO 算法進行流場計算，對全局區(qū)域進行初始化，設(shè)置殘差收斂值為0.0001，迭代步數(shù)為1000 步。

1.4.2 邊界條件

在自由出流情況下分別設(shè)置來流流量：0.01、0.05、0.10、0.15、0.20m3/s。渠道進口分別設(shè)置氣液兩相，空氣相采用壓強入口，默認為標準大氣壓強、水流相采用速度入口，分別設(shè)置進口流速和水力直徑、出口設(shè)置為自由出流、壁面設(shè)置為標準無滑移壁面。

2 結(jié)果與分析

2.1 模型合理性驗證

為了驗證數(shù)值模型和邊界條件能夠滿足實際工況，對方案1 模型進行數(shù)值模擬并與文獻試驗結(jié)果進行對比分析（文獻試驗?zāi)Ｐ蜑楹淼缹挾葹?.25m 的標準巴歇爾槽，與方案1 模型結(jié)構(gòu)尺寸一致），比較二者誤差[18]。上、下游水深與流量通過對計算結(jié)果進行數(shù)據(jù)處理得到。如表3 所示，數(shù)值模擬結(jié)果與模型試驗結(jié)果最大誤差小于5.22%，表明此次模型選取較合理，能夠滿足數(shù)值模擬計算。

表3 數(shù)值模擬結(jié)果與試驗結(jié)果誤差分析Fig.3 Error analysis of simulation results and test results

2.2 水力特性分析

2.2.1 水面線與流線

運用Tecplot 軟件對模擬結(jié)果進行后處理分析，取Z1 斷面作為研究斷面，圖5 為過流流量為0.1m3/s時，4 種方案的巴歇爾槽水面線模擬情況。水流由渠道進口流經(jīng)至上游收縮段，由于流道束窄形成雍水水位升高，此時水流流速增大形成臨界流。水流由喉道段流經(jīng)至下游擴散段，流道擴寬形成跌水水位降低，此時水流流速減小流向下游。4 種方案整體而言，喉道段前后水位落差較大，流線曲率變化在此處達到極值。水流方向發(fā)生較大的改變，分子間相互作用加劇，動能轉(zhuǎn)化為內(nèi)能，水頭損失較大。方案1、方案2、方案3 由于過渡段的存在水流流線較為平緩；方案4由于沒有過渡段連接，水流流線在渠道與巴歇爾槽交界處有較大的起伏。

2.2.2 局部水頭損失

依據(jù)不可壓縮流體恒定總流的能量方程，對巴歇爾槽上下游水頭損失計算。由于沿程渠道較短，沿程水頭損失相對于整體水頭損失比值過小，不予考慮。取上下游斷面X1、X2 斷面作為計算斷面，其能量守恒方程為[19]：

式中：z1、z2分別為上下游液位高度（m）；p1、p2分別為上下游靜水壓強（kN/m2）；γ為流體重度（N/m3）；α1、α2分別為動能修正系數(shù)，取值為1，v1、v2分別為斷面平均流速（m/s）；g 為重力加速度，取值為9.8 m2/s，hw1-2為上下游水頭損失（m）。

圖5 4 種方案水面線Fig.5 Four schemes for surface lines

由表4 所示，4 種方案下的巴歇爾槽局部水頭損失與流量之間關(guān)系，局部水頭損失隨著來流流量的增大而增加。4 種方案平均水頭損失：方案4＞方案2＞方案3＞方案1。方案1、方案2 采用弧面過渡，水流流經(jīng)圓弧壁面形成2 個分力。一部分分力垂直于壁面，產(chǎn)生損耗；另一部分分力平行于壁面，匯入主流。方案4 無過渡段直接采用邊墻連接，水流流動方向垂直于邊墻發(fā)生撞擊，動能轉(zhuǎn)化為內(nèi)能，水力損失最大。

表4 不同流量下4 種方案水頭損失Fig.4 Head loss of four schemes under different flow rates

2.2.3 測流精度分析

依據(jù)《JJG（水利）004—2015 明渠堰槽流量計計量檢定規(guī)程》流量計算公式，對4 種方案巴歇爾槽水流流量進行計算。自由出流情況下，流量與上游水頭呈現(xiàn)單一的因變關(guān)系，流量計算式為：

式中：Q為流量（m3/s）；C為流量系數(shù)；h為上游實測水頭；n為指數(shù)。

由表5 所示，4 種方案巴歇爾槽測流精度相對誤差分析，測流精度相對誤差隨來流流量的增大而減小。

當來流流量較小時，上游液位較低，液位高度讀數(shù)誤差相較于液位高度占比較大，因而產(chǎn)生較大的流量計算誤差；當來流流量增大時，上游液位較高，液位高度讀數(shù)誤差相較于液位高度占比較小，流量計算誤差變小。方案1、方案2、方案3 由于過渡段的存在，能夠發(fā)揮引導(dǎo)來流水流的功能，使上游液位波動范圍較小，讀數(shù)精度較方案4 更加準確。通過對4 種方案巴歇爾槽測流精度平均誤差計算得出結(jié)論：方案4＞方案2＞方案1＞方案3，方案3 連接形式的巴歇爾槽有較高的測流精度。

表5 不同流量下4 種方案測流精度相對誤差Fig.5 The relative errors of measuring accuracy of four schems under different flow raters

2.3 流場分析

2.3.1 速度

取Z1 斷面作為研究斷面，圖6 為過流流量為0.1m3/s 時，4 種方案巴歇爾槽縱剖面速度云圖。水流流速沿渠道進口至渠道出口方向依次增大，流速梯度變化明顯。水流流速大小與流道寬度成反比，流速變化引起水流流態(tài)由緩流向急流過渡，在喉道附近形成臨界流，達到測流目的[20]。在水流與空氣交界處，水面波動引起空氣局部范圍內(nèi)的擾動，產(chǎn)生較小的流速變化。

圖6 4 種方案速度Fig.6 Four schemes forvelocity

2.3.2 壓強

取Z1 斷面作為研究斷面，圖7 為過流流量0.1m3/s 時，4 種方案巴歇爾槽縱剖面壓強云圖。靜水壓強隨著水深的增加而增大，壓強梯度變化明顯。渠道上游雍水，上游連接段處靜水壓強最大。在巴歇爾槽的安裝施工過程中，應(yīng)當在上游連接段處襯砌加固，防止因靜水壓強過大導(dǎo)致巴歇爾槽物理結(jié)構(gòu)變形而產(chǎn)生測流誤差。

圖7 4 種方案壓強Fig.7 Four schemes forpressure

2.3.3 湍動能

取Z1 斷面作為研究斷面，圖8 為過流流量0.1m3/s 時，4 種方案巴歇爾槽縱剖面湍動能云圖。湍動能數(shù)值較大處皆位于氣相所分布的區(qū)域，因為氣體相較于液相有更大的流動性與擴散性，分子間的能量交換與傳遞更強烈。湍動能數(shù)值最大處位于水流與空氣交界處，由于水體運動帶動水氣交界處空氣做相對運動，氣體分子間做剪切運動產(chǎn)生摩擦力，消耗內(nèi)能。水相區(qū)域的湍動能分布沿著水流方向有增大的趨勢，流道的突擴和突縮改變水流的流動方向使水體內(nèi)產(chǎn)生旋渦，液體分子之間發(fā)生大量無規(guī)律的碰撞與摩擦，消耗水流能量。湍動能數(shù)值較大處一般位于分子運動活躍的區(qū)域。

圖8 4 種方案湍動能Fig.8 Four schemes for turbulent kinetic energy

3 討論

巴歇爾槽普遍應(yīng)用于我國大中型灌區(qū)，上游連接段作為流道的一部分，在施工過程中考慮其最優(yōu)形式尤為重要。巴歇爾槽進口連接段形式一方面要考慮測流精度、水頭損失的因素，另一方面也要根據(jù)渠道斷面形狀、施工條件來決定。內(nèi)接圓弧過渡段形式水頭損失小，直面過渡段形式測流精度高[11]，無過渡段形式水頭損失較大、測流精度低。

本文研究發(fā)現(xiàn)巴歇爾槽在來流流量小的情況下，測流精度較低，為了保證較高的測流精度應(yīng)在來流流量較大的情況下進行流量量測工作。局部水頭損失隨著來流流量的增加而變大，流速云圖與壓強云圖分布符合水力學(xué)定律，水流在喉道段由于流道縮窄形成臨界流此時流速最大，靜水壓強與淹沒水深呈線性關(guān)系，隨著水深的增加而變大。湍動能大小衡量分子間能量、動量的傳遞關(guān)系，氣體相較于液體有更強的流動性與擴散性。氣體分子的活躍度與無規(guī)律運動更加強烈，因此耗散性更強。

從工程建設(shè)方面考慮：由于流道上游雍水的緣故，上游連接段處靜水壓強最大，在安裝過程中對此處進行加固，防止流道變形沉降。數(shù)值模擬具有成本低、計算結(jié)果準確、后處理可視化等優(yōu)點，能夠應(yīng)用于水利工程的建設(shè)。文中僅對自由出流情況下的巴歇爾槽進行水力特性數(shù)值模擬分析，后期將對淹沒出流情況下的巴歇爾槽進行水力特性數(shù)值計算，并比較二者之間的差異。

4 結(jié)論

1）采用進口連接段過渡的巴歇爾槽相比無連接段過渡的巴歇爾槽水流流線更加平緩、起伏較小。

2）內(nèi)接圓弧過渡段形式巴歇爾槽局部水頭損失最小、無連接段過渡的巴歇爾槽局部水頭損失最大。

3）為提高測流精度，巴歇爾槽上游連接段宜采用直面過渡形式。

4）巴歇爾槽縱剖面速度云圖與壓強云圖變化梯度明顯，流速最大處位于喉道段、靜水壓強最大處位于上游雍水段，湍動能數(shù)值較大處位于流體運動活躍的區(qū)域。