蔣學(xué)煉，劉 暢，趙 悅，楊偉超*，柳淑學(xué)，朱福明

(1.天津城建大學(xué) 天津市軟土特性與工程環(huán)境重點實驗室，天津 300384；2.大連理工大學(xué) 海岸和近海工程國家重點實驗室，遼寧 大連 116024；3.天津市北洋水運水利勘察設(shè)計研究院有限公司，天津 300452)

作為兼具觀賞性和防護(hù)性以及水體交換功能的潛堤，近年來在海岸工程領(lǐng)域得到了廣泛的應(yīng)用，其水力性能的研究也成為熱點問題之一[1–2]。目前，隨著實驗技術(shù)和計算流體力學(xué)的發(fā)展，研究人員開始嘗試從流場結(jié)構(gòu)角度解析潛堤的水力性能。

Lin等[3]聯(lián)合運用激光誘導(dǎo)熒光（laser induced fluorescence，LIF）技術(shù)、粒子圖像測速（particle image velocimetry，PIV）技術(shù)、質(zhì)點追蹤（particle tracing，PT）技術(shù)考查了孤立波作用下矩形潛堤周圍的紊動射流現(xiàn)象。劉彥等[4]基于PIV圖像研究了不同流速下單體及雙體人工魚礁周圍的上升流和背渦流特征。Qie等[5]采用PIV技術(shù)在多種型式的潛堤周圍捕捉到交替瀉放的漩渦，并證實其是造成不同堤型反透射特性差異的重要原因。Zarruk等[6]運用染料示蹤法追蹤了孤立波作用下同軸雙圓柱潛體附近的渦流演變過程，認(rèn)為渦流的生成和泄放會造成動水壓力梯度的變化，從而影響建筑物所受到的波浪荷載。Athanassios等[7]基于2D–PIV量測數(shù)據(jù)分析了塊石護(hù)面陡坡的波生流場，探討了波浪爬坡過程中的流體分離和紊流場演變過程。Na等[8]聯(lián)合PIV技術(shù)和光纖反射儀（fiber optic reflectometry，F(xiàn)OR）研究了崩破波的破碎過程及能量轉(zhuǎn)換規(guī)律?？梢姡遮吘艿姆墙佑|式實驗技術(shù)已經(jīng)能夠在不干擾流場的前提下有效捕捉微觀流動結(jié)構(gòu)。

相比于模型實驗，計算流體力學(xué)（computational fluid dynamics，CFD）以其費用低、無比尺效應(yīng)、可知悉流動細(xì)節(jié)的優(yōu)點，成為波浪與建筑物相互作用研究領(lǐng)域的有力輔助工具，如相對干舷高度對矩形潛堤渦量場的影響[9]、斜坡堤后越浪水體的空間分布[10]、不同潛堤背面動水壓強的分布差異[11]、越浪水池式防護(hù)建筑物的越浪量[12]以及橢圓余弦波與弧形潛堤的相互作用[13]等。

上述研究多采用規(guī)則波、孤立波或橢圓余弦波，分別對應(yīng)深水和淺水工況，且未考慮拋石基床，而多數(shù)潛堤建設(shè)在過渡水深區(qū)域，由拋石基床支撐，其常規(guī)工況非線性的斯托克斯波形更適用。為此，本文執(zhí)行了斯托克斯2階波與拋石基床矩形潛堤相互作用的PIV試驗，獲取了迎浪區(qū)域和背浪區(qū)域的波面過程和速度場，展示了潛堤周圍不同相位時刻的渦流結(jié)構(gòu)。基于RANS–VOF格式建立了與試驗對應(yīng)的數(shù)值模型，探討了不同造波方法的適用性和海綿層的消能效果，進(jìn)一步分析了潛堤表面邊界層的流體分離現(xiàn)象。

1 試驗設(shè)計

試驗在大連理工大學(xué)海岸和近海工程國家重點實驗室開展，波浪水槽長30.00 m，寬0.40 m，高0.65 m，一端裝有液壓推板式造波機，另一端設(shè)消能斜坡（圖1（a））。水槽底部和側(cè)面均為可透光的有機玻璃。模型堤體采用丙烯酸有機玻璃制成，放置于水槽中部拋石基床上，高度0.175 m，頂寬0.200 m，內(nèi)部鉛塊配重。在后續(xù)分析中，設(shè)定迎浪堤腳為坐標(biāo)軸原點(x= 0，z=0)，水平軸順波形傳播方向為正，垂直軸向上為正。

距離迎浪堤面和背浪堤面10 cm處分別放置2支浪高儀（WG–F、WG–R），用于波浪越堤過程中的波面變形分析，波面的采樣頻率為50 Hz。采用PIV技術(shù)捕捉潛堤迎浪面和背浪面的瞬時速度場，采樣區(qū)域設(shè)置在基肩上方見圖1（a）中的FOV–F和FOV–R（field of view，F(xiàn)OV），寬度為25.5 cm，高度為18.5 cm，對應(yīng)的分辨率為57 像素/cm。數(shù)字電荷耦合攝像機（CCD）安設(shè)于水槽側(cè)面采集PIV圖像，光源由安設(shè)于槽底下方儀器室的雙頭釔鋁柘榴石晶體諧波激光器提供，脈沖延遲400 ms，經(jīng)柱面鏡散射后透過玻璃槽底和基床狹縫形成雙倍頻片光點亮測速區(qū)域（見圖1（b））。激光器觸發(fā)和CCD相機采集由信號同步器控制同步。采用INSIGHT–3G軟件對連續(xù)采集的PIV圖像序列進(jìn)行多通道互相關(guān)分析，在x方向和z方向分別設(shè)置48個和35個速度矢量點，采樣精度為每厘米1.8個速度點，每幀圖像獲取1 680個網(wǎng)格點的瞬時速度矢量值。PIV采樣頻率為15 Hz，每次采樣的兩幀圖像之間的時間延遲為600 ms。

圖1 試驗布置Fig. 1 Experiment configuration

選擇一種代表性的潛堤工況，試驗水深h=0.333 m，波高H=0.068 m，周期T=0.996 s，波長L=1.4 m。根據(jù)le Méhauté[14]的波浪理論適用范圍判斷，試驗波形為斯托克斯2階波。參考流速取淹沒深度（d=0.083 m）處的斯托克斯2階波水質(zhì)點水平速度幅值，即：

基于可靠性考慮，試驗重復(fù)3次，自第6個穩(wěn)定波峰通過各采樣區(qū)域中心位置正上方時開始采樣，采集總時長10 s，包含10個完整波浪周期。每次試驗前在測速區(qū)域放置水下標(biāo)定尺，獲取相機成像像素和片光源位置平面尺寸之間的比例，用于PIV圖像分析中網(wǎng)格點瞬時速度值的參考基準(zhǔn)。示蹤粒子為平均粒徑10 μm的聚氯乙烯粉末，均勻摻混于水槽的模型段，以提高速度捕捉的靈敏度。

2 試驗結(jié)果

2.1 波面

圖2展示了有無潛堤情況下的波面過程。

圖2 波面的時間歷程Fig. 2 Temporal variation of surface elevations

由圖2（a）顯示空槽波面與入射波面吻合良好，空槽波高和波周期的均方根偏差均值分別為0.214 2和0.007 5，說明試驗波形穩(wěn)定。由圖2（b）可以看出，越堤前后的波面均發(fā)生了非線性變形。引入不對稱度（asymmetry）和偏度（skewness）衡量波面的變形程度[15]：

式中： η (x,t)為x位置處t時刻的波面高程； η(x,t)為波面高程均值；H(·)表 示Hilbert變換；E(·)表示均值算子；σ為波面高程的標(biāo)準(zhǔn)差；不對稱度 A sym表征波浪形狀相對于垂直軸的不對稱程度，為負(fù)值時代表波面前傾（即波峰相對于垂直軸左偏斜），為正值時代表波面后傾（即波峰相對于垂直軸右偏斜）；偏度Skew表征波浪形狀相對于水平軸的不對稱程度，為正值時代表波峰上抬，波谷坦化，為負(fù)值時代表波峰下降，波谷尖化。

理論入射波面的不對稱度為 ? 3.29×10?5，偏度為0.23。與之相比，迎浪面的波形變化幅度較小，不對稱度和偏度分別為–0.21和0.04，由于潛堤的阻水效應(yīng)，其波面更為前傾。背浪面的波高減小幅度較大，說明部分波能在越堤過程中損耗，同時圖2中出現(xiàn)了次生峰谷，表明部分波能從低頻向高頻轉(zhuǎn)換，產(chǎn)生了高次約束諧波和自由諧波。背浪面波形的不對稱度和偏度分別為–0.39和0.99，相對于垂直軸和水平軸均表現(xiàn)出明顯的不對稱性。越堤過程中的波面非線性變化與潛堤附近的流場紊動有關(guān)，后續(xù)將從渦流結(jié)構(gòu)角度進(jìn)行分析。

2.2 速度場

基于PIV速度場分析波浪越堤過程中的渦流結(jié)構(gòu)。以波面通過各采樣區(qū)域正上方為準(zhǔn)，一個周期內(nèi)定義8個特征相位時刻（圖3），其中相位1對應(yīng)波峰時刻。

圖3 一個周期內(nèi)波面的相位時刻示意Fig. 3 Sketch of surface elevation phases within a wave period

圖4呈現(xiàn)了一個周期內(nèi)迎浪采樣區(qū)域（FOV–F）和背浪采樣區(qū)域（FOV–R）在不同相位時刻的速度場（矢量圖）和渦量場（云圖）。為了消除試驗過程中的波動，采用經(jīng)過相平均的速度場和渦量場，其計算式為：

由圖4可見：相位1時刻，波峰越堤，潛堤頂面射流通過背浪堤角發(fā)生流體分離，生成順時針負(fù)渦（圖4中A區(qū)），與正向水流從迎浪面攜帶來的逆時針正渦形成渦對（見相位8）。擴散過程中，順時針負(fù)渦向下游和底部對流，而逆時針正渦則向自由液面?zhèn)鬟f（相位2～相位3）。由于波浪的周期性運動，這一反向旋轉(zhuǎn)的渦對并不能充分發(fā)展，從相位4時刻開始，反向水流增強了逆時針正渦，順時針負(fù)渦開始拉伸分解，到相位7時刻基本耗散為小尺度渦旋。背浪區(qū)域的主渦對主要局限于距離堤背約1倍Keulegan–Carpenter數(shù)范圍（0.54h≈17.9 cm）內(nèi)。在整個周期過程中，背浪基肩上的渦量均很小，可視為紊流遮蔽區(qū)，工程中存在淤積的可能（圖4中B區(qū)）。

相位5時刻，波谷越堤，反向水流攜帶背浪側(cè)的逆時針正渦掃過迎浪堤角，疊加迎浪堤角處的剪切邊界層流體分離，生成了一個強逆時針正渦，并在迎浪堤面?zhèn)日T導(dǎo)出一個順時針負(fù)渦（圖4中C區(qū)）。這一渦對的運動方式不同：逆時針正渦基本維持在迎浪堤角附近；而順時針負(fù)渦則順著堤面向底床傳遞拉伸（相位6～相位7）；到相位8時刻，水流轉(zhuǎn)向，渦對隨正向水流向堤背傳遞。相對于背浪區(qū)域，迎浪堤角的主渦對運動范圍更窄，局限于距離堤面約0.5倍Keulegan–Carpenter數(shù)的范圍內(nèi)（ 0.5×0.54h≈9.0 cm）內(nèi)。與堤背不同的是，相位4時刻迎浪基肩上生成一個小型渦對（圖4中D區(qū)），隨后在上升水流和正向水流作用下向液面和下游對流擴散，逐漸分解破碎，影響范圍約兩個水質(zhì)點運動軌跡（2S=6.4 cm）。這一渦對的運動范圍較小，可能引起局部沖刷。

圖4 不同相位時刻的相平均速度場和渦量場Fig. 4 Phase-averaged velocity and vorticity fields at different phases

3 數(shù)值模型

為了理解渦流的生成機理及其相互作用，采用數(shù)值模擬進(jìn)行更深入的研究。

3.1 RANS–VOF架構(gòu)

數(shù)值波浪水槽采用RANS–VOF架構(gòu)，水流運動控制方程為描述不可壓縮黏性流動的雷諾時均納維爾–斯托克斯方程（Reynolds–averaged-Navier–Stokes equations，RANS）：

拋石基床內(nèi)的流體運動采用空間平均的納維爾–斯托克斯方程控制[17]：

波浪運動過程中的自由液面采用流體體積分?jǐn)?shù)法（volume of fluid，VOF）追蹤[18]：

式中，F(xiàn)為各計算單元內(nèi)的流體體積分?jǐn)?shù)。

上述RANS–VOF架構(gòu)采用有限差分離散、兩步映射法求解，詳細(xì)過程參考文獻(xiàn)[19–20]等。這一數(shù)值模式已應(yīng)用于眾多波浪與建筑物相互作用問題中，如護(hù)面塊體消浪能力[21]、弧形堤波浪荷載[22]、直立堤前沖刷[23]等。

3.2 數(shù)值造波

常用的數(shù)值造波方法包括邊界造波、動板造波和源函數(shù)造波3種。邊界造波法[24]根據(jù)目標(biāo)波形指定入流邊界的速度和波面時間歷程，通過動量傳遞實現(xiàn)數(shù)值造波。動板造波法[25]根據(jù)目標(biāo)波浪要素確定動板的速度時程，通過造波板的往復(fù)強迫運動引起水體波動。上述2種方法適用于弱反射或反射波能夠事前估計的情況，當(dāng)RANS水槽內(nèi)有建筑物且二次反射較強時，即使引入吸收式造波技術(shù)，長時間模擬后也可能累積較大的數(shù)值誤差。源函數(shù)造波法[26]是在水流連續(xù)性方程或動量方程中增加質(zhì)量源項或動量源項，通過源域往復(fù)釋放吸收質(zhì)量或動量形成水體波動，由于入射波和反射波分別從兩端開放邊界流出，可有效消除二次反射影響。

本文選用質(zhì)量源函數(shù)造波法。源域內(nèi)的連續(xù)性方程（6）改寫為：

式中，S(x,z,t)為 質(zhì)量源函數(shù)，?代表源域。

當(dāng)源域?qū)挾认鄬τ谀繕?biāo)波長很小時，質(zhì)量源函數(shù)S(x,z,t) 與源域正上方的目標(biāo)波形 η(t)之間的關(guān)系為：

式中，系數(shù)2表示質(zhì)量源同時產(chǎn)生傳播方向相反的兩列波，C為目標(biāo)波速， η (t)為目標(biāo)波面方程。

由斯托克斯2階波的波面方程可推導(dǎo)出質(zhì)量源函數(shù)的表達(dá)式：

丁陽等[27]建議將質(zhì)量源函數(shù)乘以一個遞增的包絡(luò)系數(shù)，使前3個周期的波面逐漸增大，以避免初始計算階段源域上方產(chǎn)生過大波動影響波面的穩(wěn)定性，包絡(luò)系數(shù)R的表達(dá)式為：

對于非線性波，在開放邊界處僅采用輻射邊界難以完全透射波能。為了減小邊界反射，基于人工衰減消波的思路，本文在開放邊界前部設(shè)置海綿層消浪。在海綿層內(nèi)，水平向動量方程（7）修正為[28]：

式中， ?Cdu為人工阻尼項，其中，Cd為阻尼系數(shù)，順海綿層由0線性增長到1.0。

3.3 模型驗證

3.3.1 波面

圖5比較了試驗波面和3種造波方法得到的數(shù)值波面。分別采用不對稱度、偏度、相關(guān)系數(shù)、歐式距離衡量波形的相似度和波面過程的累積偏差，如表1所示。其中，歐式距離公式見式（18），用于衡量試驗波面與數(shù)值波面的累計偏差。綜合比較可知，源項造波的數(shù)值波面與試驗波面的相似度最高，累積偏差最小。

表1 數(shù)值波面與試驗波面的波形參數(shù)Tab. 1 Parameters of numerical and experimental wave profiles

圖5 數(shù)值波面與試驗波面比較Fig. 5 Comparison of numerical and experimental surface elevations

式中，η1i為 試驗波面，η2i為數(shù)值波面。

3.3.2 速度場

圖6比較了迎浪面2個特征時刻瞬時速度場的PIV實測數(shù)據(jù)和數(shù)值模擬結(jié)果。從圖6整體來看，數(shù)值結(jié)果與PIV實測結(jié)果吻合良好。由于分辨率的限制和數(shù)據(jù)平滑處理，少量區(qū)域的PIV測速結(jié)果不如數(shù)值結(jié)果平順，表明數(shù)值模型可以更精細(xì)化地模擬建筑物附近的流場。

圖6 瞬時速度場PIV實測和數(shù)值模擬結(jié)果比較Fig. 6 Comparison of instantaneous velocity fields results captured by PIV and numerical simulation

圖7比較了迎浪堤腳和背浪堤腳2個特征時刻的水平速度幅值沿水深的分布。由圖7可知，數(shù)值結(jié)果與實測數(shù)據(jù)的分布趨勢和量值均吻合較好。此外，可以看出，由于渦旋引起的動量摻混，流速的分布并不符合理論的拋物形，而是趨向于均勻分布。

圖7 水平速度幅值的垂向分布Fig. 7 Vertical distribution of horizontal velocity amplitudes

上述驗證表明，本文構(gòu)建的RANS-VOF數(shù)值波浪水槽可以較準(zhǔn)確地模擬非線性波浪與矩形潛堤的相互作用。

3.4 流體分離分析

基于數(shù)值速度場分析波浪越堤過程中邊壁附近的流體分離現(xiàn)象（圖8）。以波面通過潛堤質(zhì)心正上方為準(zhǔn)，參照圖3定義特征相位時刻。 圖8中，渦量色條下為自由液面的渦量等值線，由于波浪越堤過程中的淺化效應(yīng)，發(fā)生了部分破碎，故波面始終存在小尺度渦旋。其下域中的符號A、B、C、D分別指代背浪堤角和壁面、后趾、迎浪堤角和壁面、前趾4個結(jié)構(gòu)區(qū)域。

圖8 近壁區(qū)渦量場的演變Fig. 8 Evolution of vorticity field in near-wall regions

由圖8可知：相位6～相位8時刻，一個新的波到達(dá)堤頂，C區(qū)的剪切邊界層內(nèi)開始生成負(fù)渦量，并向上部發(fā)展。相位1～相位2時刻，波峰通過潛堤上方，C區(qū)剪切邊界層內(nèi)的順時針負(fù)渦從結(jié)構(gòu)表面泄放，向自由液面和下游拉伸、對流，迎浪堤角邊界層內(nèi)開始誘導(dǎo)出反向旋轉(zhuǎn)的正渦量，同時，堤頂負(fù)渦邊界層逐漸發(fā)展并向背浪堤角傳遞。相位3時刻，水流開始反向，結(jié)構(gòu)表面正渦剪切層逐漸發(fā)展。到相位5時刻已占滿堤頂表面，切斷了順時針負(fù)渦的渦量供應(yīng)，并在迎浪堤角從結(jié)構(gòu)表面泄放出逆時針正渦，向底床和上游拉伸、對流，但這一過程并沒有充分發(fā)展，隨著正向水流的到來，逆時針正渦轉(zhuǎn)而向自由液面運動，并與之相互作用。同時，迎浪堤角泄放出的順時針負(fù)渦進(jìn)一步向自由液面?zhèn)鬟f，與之發(fā)生強烈的相互作用而耗散；背浪面A區(qū)的順時針負(fù)渦也從結(jié)構(gòu)表面泄放，向底床拉伸彎曲，隨后在反向水流的作用下與潛堤表面的正渦邊界層相互作用而耗散。在一個周期循環(huán)中，前趾D區(qū)的渦旋變化過程與C區(qū)類似，存在一個小型環(huán)流系統(tǒng)，后趾B區(qū)的渦旋強度則始終較弱。

圖8表明，隨著波浪越堤過程中的水流周期性變化，迎浪面C區(qū)和背浪面A區(qū)的剪切邊界層內(nèi)交替生成負(fù)渦量和正渦量，其生成機理可用圖9解釋。

圖9 剪切邊界層內(nèi)反向渦量的生成機理（相位8）Fig. 9 Formation mechanism of reverse vorticity in shear boundary layer (phase 8)

以相位8時刻迎浪堤面生成負(fù)渦量為例，堤角處邊界層流體分離泄放出的逆時針正渦隨著正向水流的到來向結(jié)構(gòu)物附近傳遞，使得迎浪面近壁區(qū)流體減速，動水壓力增大，形成自下而上的逆壓梯度和抽吸效應(yīng)，將水體加速推向堤角上部，在高低壓之間的剪切邊界層內(nèi)誘導(dǎo)出負(fù)渦量。

4 結(jié) 論

采用PIV技術(shù)和數(shù)值模擬研究了非線性波與矩形潛堤相互作用過程中的渦流結(jié)構(gòu)和流體分離，得到如下結(jié)論：

1）越堤過程中的波面非線性變化與潛堤附近的流場紊動有關(guān)。迎浪區(qū)域的主渦對運動范圍較窄，波能損耗小，波形變化幅度較小，背浪區(qū)域的主渦對擴散范圍較大，波高減小幅度較大，部分波能在越堤過程中損耗，部分波能從低頻向高頻轉(zhuǎn)換，產(chǎn)生了高次約束諧波和自由諧波，波形相對于垂直軸和水平軸均表現(xiàn)出明顯的不對稱性。

2）在波浪作用下，反向旋轉(zhuǎn)的渦對在迎浪面和背浪面周期性地生成、泄放和傳遞，但均未能充分發(fā)展，由隨后到來的變向水流攜帶，與自由液面或結(jié)構(gòu)表面相互作用而耗散。數(shù)值模擬結(jié)果表明，渦旋的渦量供應(yīng)主要來自結(jié)構(gòu)表面的剪切邊界層，而這些反旋渦量是由分離渦運動引起的逆壓梯度誘導(dǎo)生成的。

3）目前的工程設(shè)計主要基于無黏勢流理論，未考慮渦旋引起的復(fù)雜流動效應(yīng)，但研究結(jié)果表明渦旋的生成、泄放、拉伸、對流、耗散會顯著改變局部流場，進(jìn)而影響結(jié)構(gòu)物的局部沖淤和受力。RANS–VOF數(shù)值波浪水槽為潛堤周圍渦流效應(yīng)的精細(xì)化研究提供了一種有效工具，可用于波浪與潛堤相互作用的參數(shù)分析和工程設(shè)計方案優(yōu)化。[致謝] 衷心感謝赫瑞瓦特大學(xué)（Heriot–Watt University）的Zou Qingping教授為本文工作提供的技術(shù)指導(dǎo)，衷心感謝四川大學(xué)的林鵬智教授為本文的數(shù)值模擬工作提供的RANS–VOF源代碼。