陳 剛，董增川，王海軍，顧世祥，楊紅宣

(1.河海大學(xué)水文水資源與水利工程科學(xué)國家重點(diǎn)實(shí)驗(yàn)室，江蘇 南京 210098；2.云南省水利水電勘測設(shè)計(jì)研究院，云南 昆明 650021；3.昆明理工大學(xué)電力工程學(xué)院，云南 昆明 650500)

河冰是陸地冰凍圈的重要組成要素[1-2]，盡管覆蓋范圍相對較小，但對人類社會的影響卻最為顯著[3]。河冰生消演變會改變河流中物質(zhì)、能量的輸移規(guī)律[4-5]，在特定的水文氣象、河勢和水力條件下甚至造成冰凌災(zāi)害，危及沿岸人民的生命財(cái)產(chǎn)安全[6]。中國陸地冰凍圈占陸地國土面積的78%[1]，冰情影響范圍廣，冰凌災(zāi)害頻發(fā)[7]。冰蓋是河冰的重要存在形式，結(jié)冰蓋輸水是中國長距離調(diào)水、灌溉渠系工程實(shí)現(xiàn)冬季安全輸水的常用運(yùn)行方式[8]。冰蓋形成附加的無滑移邊界對河渠水流紊動特性的影響是河渠冰凌災(zāi)害防治、輸水渠系冰期運(yùn)行、調(diào)水工程設(shè)計(jì)與運(yùn)行調(diào)度的重要理論基礎(chǔ)。因此，研究冰蓋下水流紊動特性對自然災(zāi)害防治和供水安全保障具有重要意義。

根據(jù)河冰生消演變對水流表面的影響，水流流態(tài)可分為明流、完全冰封和部分冰封[4，9]。部分冰封是指水流表面部分被冰覆蓋(稱為冰蓋區(qū))，部分為自由水面(明流區(qū))，是前2種流態(tài)相互轉(zhuǎn)換的過渡流態(tài)，岸冰是其主要存在形式。在水流紊動特性方面，相比明渠流已取得系統(tǒng)性研究成果，冰蓋下水流的研究相對較少，尤其是部分冰封。早期研究側(cè)重于揭示河冰生消的機(jī)理[10]，僅有少數(shù)研究者通過原型觀測分析冰蓋形成對流速分布和輸水能力的影響[11]。隨著先進(jìn)實(shí)驗(yàn)儀器的問世和室內(nèi)試驗(yàn)條件的改善，國內(nèi)外學(xué)者開始通過水槽試驗(yàn)研究冰蓋下水流的特性。Parthasarathy等[12]探究完全冰封水流縱向時(shí)均流速、雷諾應(yīng)力和紊動能的分布規(guī)律；陳建國等[13]測量不同冰蓋組合下的水流結(jié)構(gòu)，探討冰蓋流的輸水能力；拾兵等[14]研究不同流速、水深、冰蓋糙率條件下完全冰封水流的縱向紊動強(qiáng)度，發(fā)現(xiàn)冰蓋下水流縱向紊動強(qiáng)度比明渠流大；Peters等[15]測量岸冰下水流冰蓋區(qū)和明流區(qū)的縱向流速分布，發(fā)現(xiàn)岸冰形成使流動核心區(qū)向明流區(qū)收縮，導(dǎo)致明流區(qū)流速增加，冰蓋區(qū)流速減??；Nyantekyi-Kwakye等[16]采用PIV技術(shù)測量不同寬度岸冰下水流瞬時(shí)流場，用于檢驗(yàn)4種邊界切應(yīng)力計(jì)算方法的適用性；Wang等[17]測量岸冰下水流流速分布，用于校驗(yàn)冰蓋下水流橫向流速分布準(zhǔn)二維模型?？梢姡F(xiàn)有的研究多針對單一流態(tài)，對比研究3種流態(tài)紊動特性的成果還較少。

本文在國內(nèi)外已有研究成果的基礎(chǔ)上，通過水槽試驗(yàn)測量3種流態(tài)的恒定均勻流場，從縱向時(shí)均流速、雷諾應(yīng)力和紊動強(qiáng)度等方面分析冰蓋對水流紊動特性的影響。

1 試驗(yàn)裝置和研究方法

1.1 試驗(yàn)裝置

試驗(yàn)在昆明理工大學(xué)水流與水工結(jié)構(gòu)重點(diǎn)實(shí)驗(yàn)室的自循環(huán)玻璃水槽試驗(yàn)系統(tǒng)中進(jìn)行，該系統(tǒng)由地下水池、泵站、量水堰、穩(wěn)流池、水槽和回水池組成，見圖1。水槽采用有機(jī)玻璃制作，長15 m，寬0.80 m，深0.35 m，縱坡坡降固定為0.5‰。水流通過泵站從地下水池提升至進(jìn)水渠，經(jīng)量水堰進(jìn)入穩(wěn)流池，溢流后進(jìn)入玻璃水槽，從水槽尾門出流后經(jīng)回水池進(jìn)入地下水池，形成水流循環(huán)。為保證水槽的水流流態(tài)穩(wěn)定，水槽進(jìn)口端為深1.4 m、長4.5 m、寬6.5 m的穩(wěn)流池，通過調(diào)節(jié)水槽末端的尾門來控制水槽中的水位，使水槽保持恒定均勻流。量水堰為無側(cè)收縮矩形薄壁堰，用于測量流量，堰上水頭由精度為0.1 mm的水位測針測量得到。流速采用挪威諾泰克公司(Nortek AS)生產(chǎn)的聲學(xué)多普勒點(diǎn)式三維測速儀Vectrion測量得到。

圖1 試驗(yàn)裝置示意Fig.1 Sketch of experimental set-up

1.2 試驗(yàn)材料

試驗(yàn)探究冰蓋下水流的紊動特性，典型過流斷面見圖2(a)，其中，B為斷面寬度，b1為單側(cè)冰蓋的寬度，H為水深。以冰蓋邊緣為界劃分為冰蓋區(qū)和明流區(qū)，則斷面冰蓋覆蓋度定義為F=2b1/B，當(dāng)b1=0時(shí)，F(xiàn)=0；當(dāng)b1=B/2時(shí)，F(xiàn)=1，即明流、完全冰封可視為冰蓋覆蓋度分別為0、1的特例。試驗(yàn)以F為變量，共設(shè)計(jì)6種工況，F(xiàn)值分別為0、0.250、0.375、0.500、0.750、1.000，詳見表1。

在普通實(shí)驗(yàn)室中開展冰蓋下水流水力特性研究，通常采用模擬冰蓋替代真冰進(jìn)行水槽試驗(yàn)，利于形成穩(wěn)定流場。采用購置成本低、制作簡便的成型單塊泡沫板制作模擬冰蓋，使用火把T2型電熱絲泡沫切割機(jī)將泡沫板切割成長2.00 m，厚0.05 m，寬分別為0.10、0.15、0.20、0.30和0.40 m的條板。為確保水槽進(jìn)口水流流態(tài)平順穩(wěn)定，距進(jìn)口0.80 m段為槽內(nèi)穩(wěn)流段(圖1)，同時(shí)為避免尾門調(diào)節(jié)影響模擬冰蓋，距尾門0.20 m段不布置冰蓋，布置模擬冰蓋的渠段長14.0 m，每種工況各準(zhǔn)備14塊泡沫條板。

圖2 冰蓋下典型過流斷面和試驗(yàn)測量斷面布置示意Fig.2 Sketch of a typical under-ice section and measuring cross section

表1 試驗(yàn)工況

1.3 試驗(yàn)方案

試驗(yàn)時(shí)貼近水槽邊壁布置切割好的泡沫條板，并利用水槽橫梁上的旋轉(zhuǎn)螺絲桿進(jìn)行固定，使其貼近水槽邊壁，不能隨流波動，但能上下浮動。同時(shí)，由于真冰和泡沫板密度差異較大，通過旋轉(zhuǎn)螺絲桿施加垂向壓力，使泡沫板浸沒在水中的深度與同體積真冰的浸沒深度相等。

為確保試驗(yàn)在準(zhǔn)恒定均勻流條件下進(jìn)行，需要調(diào)整尾門使水槽內(nèi)水面近似平行于床面，為此布置5個水深測量斷面，分別距水槽進(jìn)口5.0、6.5、7.5、8.5和10.5 m，見圖2(b)。各斷面都安裝有水位測針，其零刻度線與水槽底板齊平，精度為±0.1 mm。在冰蓋寬度和流量固定時(shí)，通過調(diào)整尾門高度，使各斷面的水深相等。為確保水流充分發(fā)展和避免尾門的影響，流速測量斷面為水槽中段的C-C斷面，結(jié)合斷面的對稱性，僅測量半寬的斷面。自水槽右邊壁起，每間隔5 cm布置1條測線，共布置8條測線，依次記為Y05、Y10、Y15、Y20、Y25、Y30、Y35和Y40，每條垂線上均布置10個測點(diǎn)，見圖3。對于部分冰封工況，均布置1條測線與冰蓋邊緣垂線重合。采用標(biāo)準(zhǔn)Vectrion測量各測點(diǎn)的瞬時(shí)流速，采樣頻率設(shè)為25 Hz，采樣時(shí)間不少于60 s，即每個測點(diǎn)至少采樣1 500次。根據(jù)輸出數(shù)據(jù)的振幅數(shù)、信噪比和相關(guān)數(shù)，對采樣的瞬時(shí)流速數(shù)據(jù)點(diǎn)進(jìn)行質(zhì)量控制[18]，符合要求的用于計(jì)算時(shí)均流速與脈動流速，進(jìn)而計(jì)算得到雷諾應(yīng)力和紊動強(qiáng)度，用于揭示冰蓋對水流紊動特性的影響。

圖3 C-C測量斷面上測點(diǎn)布置示意Fig.3 Sketch of measuring points on the C-C measuring section

2 結(jié)果與分析

2.1 縱向時(shí)均流速分布

時(shí)均流速計(jì)算公式如下：

(1)

式(1)表明時(shí)均流速受采樣時(shí)間的影響。以B20工況Y20測線上靠近冰蓋的測點(diǎn)為例，其縱向時(shí)均流速與采樣時(shí)間的關(guān)系見圖4，其中u′、v′、w′表示x、y、z方向的脈動流速，為瞬時(shí)流速與時(shí)均流速之差。

圖4 采樣時(shí)間對時(shí)均流速的影響及其脈動流速分布Fig.4 Influence of sampling time on time average velocity and its fluctuating velocity distribution

由圖4可見，當(dāng)采樣時(shí)間較短時(shí)，縱向時(shí)均流速波動較大；當(dāng)采樣時(shí)間為30 s時(shí)，縱向時(shí)均流速趨于穩(wěn)定值，因此，試驗(yàn)確定采樣時(shí)間為60 s是合理的。同時(shí)，橫向、垂向時(shí)均流速均不等于0，這表明水流除縱向的主流動外，存在疊加于主流之上的二次流；橫向時(shí)均流速小于0，表明該點(diǎn)的流量總體上是從明流區(qū)流向冰蓋區(qū)；垂向時(shí)均流速大于0，表明該點(diǎn)流量總體上是流向無滑移邊界的。在脈動流速分布上，3個方向的瞬時(shí)流速在時(shí)均流速附近波動，脈動流速均滿足正態(tài)分布。

常用于描述冰蓋下縱向時(shí)均流速的分布律有雙對數(shù)律和雙冪律。雙對數(shù)律存在最大流速處不連續(xù)、當(dāng)量粗糙高度不易獲取、紊流核心區(qū)垂向流速分布不符合雙對數(shù)律[19]等問題，而雙冪律能描述冰蓋下水流垂向流速光滑連續(xù)的特性，與實(shí)測資料也較吻合[12，20]，故采用雙冪律描述連續(xù)冰蓋下水流縱向時(shí)均流速的垂向分布：

(2)

式中：k0為常數(shù)；mb、mi分別為與床面、冰蓋粗糙度有關(guān)的冪指數(shù)。特別地，當(dāng)mi→∞時(shí)，k0→umax，umax為最大流速。式(2)簡化為常用于描述明流縱向時(shí)均流速垂向分布的單冪律。采用平均百分比誤差和擬合度指標(biāo)評價(jià)雙冪律的模擬性能，其中平均百分比誤差(EMAP)為相對誤差絕對值的算術(shù)平均值，計(jì)算公式為

(3)

式中：i表示數(shù)據(jù)點(diǎn)數(shù)；n為數(shù)據(jù)量；Y為模型預(yù)報(bào)值；X為實(shí)測值。EMAP最小值為0，其值越小說明模型模擬值越接近參考數(shù)據(jù)。擬合度指標(biāo)是常用于判定非線性回歸方程擬合度的統(tǒng)計(jì)參數(shù)，有別于線性回歸的確定系數(shù)，故記為IFO，其計(jì)算公式為

(4)

式中：IFO最大值為1，其值越接近1，說明模型預(yù)報(bào)值與參考值的吻合程度越好。對于明流O00工況，各垂線的EMAP(u)均小于5.0%，IFO(u)均大于0.90，這表明明流工況縱向流速的垂向分布符合指數(shù)律，見圖5(a)。該工況的寬深比為4.17，受二次流影響最大流速點(diǎn)偏離自由水面[21]，特別是靠近邊壁的Y05測線。

圖5 O00工況和C40工況縱向時(shí)均流速垂向分布Fig.5 Vertical profiles of streamwise velocity for Run O00 and Run C40

完全冰封C40工況的實(shí)測流速與式(2)計(jì)算值進(jìn)行對比見圖5(b)。當(dāng)mb=9.69、mi=7.99時(shí)，各條垂線的EMAP(u)均小于5.0%，IFO(u)均大于0.90，這表明雙冪律可用于描述連續(xù)冰蓋下水流縱向時(shí)均流速的垂向分布。受冰蓋附加無滑移邊界的影響，各垂線最大流速點(diǎn)位置不同于明流工況，位于垂線中部區(qū)域。理論上，最大流速點(diǎn)的速度梯度為0，故將雙冪律對z求導(dǎo)，令導(dǎo)數(shù)為0得到最大流速點(diǎn)的相對水深為

(5)

部分冰蓋形成導(dǎo)致冰蓋區(qū)上邊界為無滑移邊界，明流區(qū)上邊界為自由水面，水流結(jié)構(gòu)較明流工況和完全冰封工況都復(fù)雜[15]。目前，針對部分冰蓋下水流流速垂向分布的理論分析較少。圖6給出B10、B15、B20和B30工況各垂線上的縱向時(shí)均流速的垂向分布。由圖6可見，部分冰封工況縱向時(shí)均流速垂向分布兼具明流和完全冰封工況的特征：① 在冰蓋區(qū)，所有垂線均類似于完全冰封工況，采用C40工況mi和mb的值，由式(2)得到冰蓋區(qū)縱向時(shí)均流速的計(jì)算值，并同實(shí)測值對比，結(jié)果表明除B20工況Y15垂線外，各垂線的EMAP(u)均小于5.0%，IFO(u)均大于0.90，表明部分冰蓋下水流冰蓋區(qū)縱向流速分布可采用雙冪律描述；② 明流區(qū)部分垂線類似于明流工況，采用明流O00工況冪指數(shù)mb的值，由式(2)得到明流區(qū)部分垂線縱向流速垂向分布，各垂線EMAP(u)均小于5.0%，IFO(u)均大于0.90，表明部分冰蓋下水流明流區(qū)靠近紊流核心區(qū)的流速垂向分布符合單冪律；③ 在兩者間存在過渡區(qū)，縱向時(shí)均流速垂向分布既不同于明流工況，也不同于完全冰封工況。根據(jù)實(shí)測值按雙冪律進(jìn)行擬合，得到的mb均小于C40工況的值，但相反地mi均大于C40工況的值。

圖6 部分冰封工況縱向時(shí)均流速垂向分布Fig.6 Vertical profiles of streamwise velocity in open channels with partial ice cover

2.2 雷諾應(yīng)力垂向分布

聯(lián)立連續(xù)性方程和動量方程，得到恒定均勻流的控制方程為[23]

(6)

(7)

(8)

(9)

(10)

(11)

式中：u*為斷面摩阻流速，這表明紊流核心區(qū)雷諾應(yīng)力在水深方向上呈線性分布。圖7(a)對比了明流O00工況Y40垂線上雷諾應(yīng)力實(shí)測值和計(jì)算值，結(jié)果表明，明流工況紊流核心區(qū)的雷諾應(yīng)力在水深方向上呈線性分布，雷諾應(yīng)力的最大值出現(xiàn)在床面附近。但是，對應(yīng)于最大流速點(diǎn)因二次流影響偏離自由水面，雷諾應(yīng)力零值點(diǎn)也偏離自由水面。

(12)

采用斷面摩阻流速對切應(yīng)力進(jìn)行量綱一化。采用試驗(yàn)測量的斷面尺寸和比降，通過下列公式計(jì)算得到u*：

(13)

(14)

式中：hb為床面層的水深；c1、c2為回歸系數(shù)。顯然，切應(yīng)力零值點(diǎn)相對水深不同于最大流速點(diǎn)，兩者相對位置與床面、冰蓋的相對粗糙程度有關(guān)。以C40工況為例，由式(14)計(jì)算的切應(yīng)力零值點(diǎn)相對水深為0.399，位于最大流速點(diǎn)和相對光滑的固壁之間。

圖7 O00工況和C40工況雷諾應(yīng)力垂向分布Fig.7 Vertical profiles of shear stress for Run O00 and Run C40

對于部分冰封工況，由于冰蓋區(qū)、明流區(qū)的垂向流速分布差異較大，存在橫向動量交換，由式(7)推求雷諾應(yīng)力時(shí)不能忽略二次流的影響。目前尚未建立雷諾應(yīng)力的數(shù)學(xué)模型，僅能從測量數(shù)據(jù)探究其垂向分布規(guī)律。以部分冰封B20、B30工況為例，冰蓋邊緣垂線上雷諾應(yīng)力分布見圖8。受橫向動量交換影響，雷諾應(yīng)力在水深方向上散亂分布，不滿足線性分布的規(guī)律。

圖8 部分冰封工況雷諾應(yīng)力垂向分布Fig.8 Vertical profiles of shear stress in open channels with partial ice cover

2.3 紊動強(qiáng)度垂向分布

紊動強(qiáng)度是流場內(nèi)任一空間點(diǎn)脈動流速的二階中心矩，表示瞬時(shí)流速在其均值附近的分散程度，反映水體流速脈動的強(qiáng)弱程度。對于明流工況，Nezu等[25]采用指數(shù)律描述其縱向紊動強(qiáng)度垂向分布：

(15)

式中：urms為縱向紊動強(qiáng)度；Du和ξu為經(jīng)驗(yàn)系數(shù)。根據(jù)明流O00工況實(shí)測的瞬時(shí)流速，計(jì)算得到各垂線的縱向紊動強(qiáng)度，見圖9(a)。從紊動強(qiáng)度垂向分布來看，在床面附近存在最大值，隨水深增加逐漸減小。當(dāng)相對水深z/H>0.5后，紊動強(qiáng)度總體上趨于穩(wěn)定，在自由水面附近達(dá)到最小。由圖9(a)可見，明流O00工況，各垂線紊動強(qiáng)度垂向分布符合指數(shù)律，EMAP(urms)均值為9.43%，IFO(urms)均值為0.895 6。從橫向分布來看，總體上越靠近邊壁，紊動強(qiáng)度越大，越靠近紊流核心區(qū)，紊動強(qiáng)度越小。結(jié)合縱向流速的垂向分布，紊動強(qiáng)度在最大流速附近達(dá)到最小值，說明最大流速處流速脈動較小，相鄰流層間的動量交換最弱。

圖9 O00工況和C40工況紊動強(qiáng)度垂向分布Fig.9 Vertical profiles of turbulence intensity for Run O00 and Run C40

對于完全冰封工況，按照雙層假定將過流斷面在垂向上劃分為床面層和冰蓋層，在各流層中紊動強(qiáng)度垂向分布均符合指數(shù)律[26]，則

(16)

式中：Du，b、Du，i、ξu，b和ξu，i為經(jīng)驗(yàn)系數(shù)。根據(jù)C40工況實(shí)測的瞬時(shí)流速，計(jì)算得到各垂線的紊動強(qiáng)度，見圖9(b)。在垂向上，紊動強(qiáng)度在床面和冰蓋附近存在最大值，隨著到無滑移邊界距離的增加而逐漸減小，在最大流速附近達(dá)到最小值。經(jīng)回歸分析，C40工況各垂線上床面層和冰蓋層的縱向紊動強(qiáng)度均滿足指數(shù)分布的規(guī)律，EMAP(urms)=6.94%，IFO(urms)=0.920 3。

對于部分冰封工況，目前尚未提出描述其紊動強(qiáng)度的分布律，僅能通過分析試驗(yàn)測量數(shù)據(jù)探究其垂向分布規(guī)律。圖10給出了B10、B15、B20和B30工況各測線上紊動強(qiáng)度垂向分布，總體上無滑移邊界(床面、冰蓋)附近的紊動強(qiáng)度相對較大，且邊界相對粗糙的冰蓋區(qū)紊動強(qiáng)度大于邊界相對光滑的床面區(qū)。部分冰封4種工況冰蓋區(qū)urms/u*均值分別為1.68、1.87、1.90和1.79，明流區(qū)urms/u*均值分別為1.41、1.66、1.78和1.62，冰蓋區(qū)紊動強(qiáng)度大于明流區(qū)紊動強(qiáng)度。當(dāng)F≤0.5時(shí)，冰蓋區(qū)、明流區(qū)urms/u*均隨F的增加而增大；當(dāng)F>0.5時(shí)，相比B20工況，B30工況冰蓋區(qū)、明流區(qū)urms/u*均減小，紊動強(qiáng)度減弱。

3 結(jié) 論

在普通實(shí)驗(yàn)室循環(huán)水槽系統(tǒng)中開展模型試驗(yàn)，采用三維聲學(xué)多普勒測速儀測量不同冰蓋覆蓋度下恒定均勻流場，從縱向流速、雷諾應(yīng)力和紊動強(qiáng)度等方面分析冰蓋下水流紊動特性，主要得到以下結(jié)論：

(2) 明流工況雷諾應(yīng)力在水深方向上符合線性分布，最大值出現(xiàn)在床面附近；完全冰封工況雷諾應(yīng)力在水深方向也符合線性分布，切應(yīng)力零值點(diǎn)位于最大流速點(diǎn)和相對光滑邊界之間；部分冰封工況受橫向動量交換的影響，雷諾應(yīng)力在水深方向上不符合線性分布。

(3) 明流工況紊動強(qiáng)度垂向分布符合指數(shù)律，完全冰封工況基于雙層假定劃分的冰蓋層和床面層中紊動強(qiáng)度可分別采用指數(shù)律描述，部分冰封工況紊動強(qiáng)度受冰蓋影響在橫向上存在明顯的差異。