呂超凡，趙西增,2，殷銘簡(jiǎn)

(1. 浙江大學(xué) 海洋學(xué)院，浙江 舟山 316021; 2. 浙江海洋大學(xué) 海洋工程裝備學(xué)院，浙江 舟山 316022)

重力式防波堤雖然能較好地起到擋浪作用，但它在一定程度上阻止港池內(nèi)外水體之間的自然循環(huán)和交換[1-2]，使得港內(nèi)水質(zhì)惡化，產(chǎn)生一系列的環(huán)境問題[3]。為在擋浪的同時(shí)，加強(qiáng)港池內(nèi)外水體交換，可在防波堤的水下部分開挖涵洞，形成涵洞式防波堤。利用波能主要集中在水面附近的特點(diǎn)，防波堤的上部結(jié)構(gòu)可以反射波浪，水下涵洞部分可進(jìn)行水體交換[4]。該措施既可滿足港池的泊穩(wěn)條件，又可提高港池的水質(zhì)；相較于其他常見的水體交換方式，具有結(jié)構(gòu)強(qiáng)度高、經(jīng)濟(jì)性好、施工容易等優(yōu)點(diǎn)[2]。

總體而言，針對(duì)涵洞式防波堤的消浪及透波特性研究己較為詳盡。但在水體交換更新方面，雖然已有基于淺水方程的大尺度數(shù)值模擬，但其針對(duì)的是潮流，波浪作用下的水體交換機(jī)理及其規(guī)律尚缺乏相關(guān)研究。為此，采用計(jì)算流體力學(xué)(CFD)數(shù)值模擬[11-13]結(jié)合物理模型試驗(yàn)[14-15]的方式精細(xì)化研究涵洞式直立堤在波浪作用下涵洞內(nèi)部水體的交換特性。數(shù)值模擬基于高精度的VPM (volume-average/point-value method)-THINC (tangent of hyperbola for interface capturing)/QQ (quadratic surface representation and Gaussian quadrature)模型[16-17]，該模型可以精確捕捉自由面，重現(xiàn)波浪與涵洞式直立堤的相互作用，為分析其運(yùn)動(dòng)機(jī)理提供依據(jù)[18]。此外，在其基礎(chǔ)上提出一種能夠標(biāo)記涵洞內(nèi)外水體變化的雙液相流體體積(VOF)方法[11]。該方法主要是把原有的計(jì)算水體分成多個(gè)子水體并分別標(biāo)記，而后計(jì)算其各自輸運(yùn)方程，進(jìn)而獲取不同水體的運(yùn)動(dòng)特征，以此來反映在水動(dòng)力的作用下，水體內(nèi)部的流動(dòng)規(guī)律。該方法可以同時(shí)獲取涵洞內(nèi)部的水體流動(dòng)以及整個(gè)水體的流動(dòng)，相對(duì)簡(jiǎn)單，計(jì)算成本低[19]。由此，文中利用上述方法定性以及定量描述涵洞內(nèi)外水體的交換特性，而后結(jié)合局部流場(chǎng)揭示水體交換機(jī)理，并分析不同涵洞參數(shù)以及波浪參數(shù)對(duì)水體交換的影響。

1 數(shù)值方法

1.1 基本方程

CFD數(shù)值模擬采用不可壓縮黏性流模型，其控制方程如下：

(1)

(2)

(3)

式中：U(x,z)為流體質(zhì)點(diǎn)速度，ρ為流體密度，p為相對(duì)動(dòng)壓力，m為動(dòng)力黏性系數(shù)，U(x,z)為加速度，U(x,z) 為坐標(biāo)點(diǎn)，F(xiàn)σ(x,z)為表面張力。Fσ表達(dá)式如下：

(4)

式中：σ為張力系數(shù)，κ為界面平均曲率，α為流體體積分?jǐn)?shù)。當(dāng)網(wǎng)格單元為空氣時(shí)，α=0；為水時(shí)，α=1；水汽摻混時(shí)，0<α<1。另外每個(gè)網(wǎng)格的流體特性定義如下：

λ=λ1α+λ2(1-α)

(5)

式中：λ可為網(wǎng)格內(nèi)流體的密度ρ或黏性系數(shù)μ。

其中，Navier-Stokes方程采用高階有限體積法VPM求解，自由面采用高精度的THINC/QQ算法[17]重構(gòu)。建立基于OpenFOAM底層函數(shù)庫(kù)的VPM-THINC/QQ模型，具體可參考文獻(xiàn)[18]。

為實(shí)現(xiàn)數(shù)值波浪水槽的造波與消波，采用速度邊界法造波[18]，以及 Jacobsen等[20]提出的松弛法來消除透射波和二次反射波。

1.2 雙液相算法

為直觀反映涵洞式直立堤的水體交換特性，在VPM-THINC/QQ模型的基礎(chǔ)上提出一種能夠模擬港池內(nèi)外水體交換的雙液相VOF方法[19]。該方法假設(shè)兩種性質(zhì)相同的液體分別為β1和β2，兩種液體均滿足VOF輸運(yùn)方程：

(6)

(7)

其中，

β1+β2=α

(8)

(9)

(10)

由此，港池外部以及涵洞內(nèi)部水體標(biāo)記為β1，港池內(nèi)部水體標(biāo)記為β2，開展基于雙液相VOF方法的涵洞式直立堤水體交換特性研究。

2 試驗(yàn)以及CFD數(shù)值模型布置

為定性以及定量分析涵洞式直立堤的水體交換特性，物理模型試驗(yàn)以及數(shù)值模擬布置如下。

2.1 試驗(yàn)設(shè)備及模型布置

物理模型試驗(yàn)在波浪水槽內(nèi)進(jìn)行，水槽尺寸為25 m×0.7 m×0.7 m(長(zhǎng)×寬×高)。水槽首端裝有主動(dòng)吸收式推板造波機(jī)，可產(chǎn)生單向二維正弦規(guī)則波，目標(biāo)波浪生成質(zhì)量不受造波時(shí)間影響，水槽尾端配有消波設(shè)施來消除波浪反射。模型試驗(yàn)布置見圖1，試驗(yàn)水槽中設(shè)置6根浪高儀，分別為EG1～EG6，所在位置依次為5 m、10 m、11 m、13 m、16 m、17 m，防波堤模型置于14 m處。模型采用亞克力板制作，防波堤長(zhǎng)度為B=0.5 m，涵洞高度及深度可調(diào)節(jié)，試驗(yàn)水深h=0.4 m。具體試驗(yàn)工況如表1所示，其中S為涵洞高度，dS為涵洞中軸線位置，Ti為試驗(yàn)周期，Hi為試驗(yàn)波高。

圖1 物理模型試驗(yàn)布置Fig. 1 Sketch of an experimental setup

表1 模型試驗(yàn)以及CFD模擬工況設(shè)置Tab. 1 Conditions of model experiments and CFD simulations

2.2 CFD數(shù)值模擬布置

CFD數(shù)值波浪水槽布置見圖2，水槽布置與模型試驗(yàn)布置基本一致。但數(shù)值水槽另外在14.25 m處設(shè)置流速測(cè)點(diǎn)FV1，用來監(jiān)測(cè)整個(gè)涵洞斷面的水平流速u(x,z)，進(jìn)而計(jì)算出整個(gè)斷面的振蕩流量Qn；同時(shí)，在0～4 m處設(shè)置波浪造波松弛區(qū)，用來吸收從結(jié)構(gòu)物傳來的反射波浪，18～26 m設(shè)置為消波區(qū)。數(shù)值模擬工況與模型試驗(yàn)一致，具體見表1。

圖2 數(shù)值波浪水槽示意Fig. 2 Schematic diagram of a numerical wave flume

2.3 模型驗(yàn)證

為驗(yàn)證數(shù)值模擬的正確性，需將數(shù)值模擬結(jié)果與試驗(yàn)結(jié)果進(jìn)行對(duì)比。圖3表示在涵洞深度dS=0.177 m，涵洞高度S=0.096 m 的條件下，CFD結(jié)果與試驗(yàn)結(jié)果的波面曲線對(duì)比，其中圖3(a)表示波高Hi=0.08 m,波周期Ti=1.57 s的波面對(duì)比，圖3(b)表示Hi=0.02 m,Ti=1.98 s的波面對(duì)比。從圖3中可以看出，圖3(a)的波浪透射試驗(yàn)結(jié)果以及CFD模擬結(jié)果幾乎一致，圖3(b)的透射結(jié)果略微偏大主要是由于試驗(yàn)消波區(qū)略短，存在部分反射波與透射波的疊加，故隨著時(shí)間的增長(zhǎng)，透射值逐漸變大；反射結(jié)果兩者略有差異，但整體差異不大。上述結(jié)果表明雙液相VPM-THINC/QQ模型能較好重現(xiàn)實(shí)際波浪與涵洞式防波堤的相互作用，具有一定可靠性。

圖3 試驗(yàn)及CFD波面對(duì)比(dS=0.177 m,S=0.096 m)Fig. 3 Comparisons of experimental results and CFD wave surface (dS=0.177 m,S=0.096 m)

為分析涵洞內(nèi)部水體輸移規(guī)律，試驗(yàn)中，涵洞內(nèi)部水體用染色劑進(jìn)行染色。該部分染色水體在一段時(shí)間內(nèi)，能定性的反映出波浪作用下的水體運(yùn)動(dòng)規(guī)律。圖4表示在波高Hi=0.06 m，周期Ti=2.4 s，涵洞深度dS=0.177 m，涵洞高度S=0.046 m條件下試驗(yàn)與數(shù)值模擬的涵洞內(nèi)部水體輸移軌跡對(duì)比。圖4(a)、4(b)分別表示初始時(shí)刻以及波浪運(yùn)動(dòng)初期數(shù)值模擬條件下的雙液相水體分布情況；圖4(c)、4(d)分別表示初始時(shí)刻以及波浪運(yùn)動(dòng)初期試驗(yàn)條件下的染色水體分布情況。從圖4中可以看出，當(dāng)波浪傳到結(jié)構(gòu)物時(shí)，涵洞內(nèi)部水體在波浪的作用下會(huì)向防波堤內(nèi)側(cè)輸移，并在涵洞口形成近乎對(duì)稱的兩個(gè)渦旋，逐漸把涵洞內(nèi)部水體帶出。對(duì)比試驗(yàn)與數(shù)值結(jié)果發(fā)現(xiàn)，兩者現(xiàn)象幾乎一致，渦旋軌跡及大小也近乎相同，表明雙液相模型能夠較好地模擬出波浪作用下的水體輸移規(guī)律。

圖4 雙液相驗(yàn)證(Hi=0.06 m,Ti=2.4 s ,dS=0.177 m,S=0.046 m)Fig. 4 Verifications of two-liquid phase (Hi=0.06 m,Ti=2.4 s ,dS=0.177 m,S=0.046 m)

模擬與試驗(yàn)結(jié)果對(duì)比表明，文中所提出的雙液相VPM-THINC/QQ模型能夠真實(shí)反映波浪作用下涵洞式直立堤的水體交換特性。

3 結(jié)果與討論

由于已經(jīng)驗(yàn)證文中數(shù)值模擬能夠較好重現(xiàn)波浪與涵洞式直立堤的相互作用，故下文結(jié)果均為數(shù)值所得。為探究涵洞式直立堤的水體交換特性，定量分析涵洞內(nèi)部的流體運(yùn)動(dòng)規(guī)律，定義涵洞內(nèi)14.25 m處的水平振蕩流量Qn如下：

(11)

3.1 涵洞內(nèi)水體交換的機(jī)理分析

為分析涵洞式直立堤的水體交換機(jī)理，通過分析涵洞附近的雙液相分布、渦量分布以及流線分布規(guī)律對(duì)其展開研究。其中雙液相為了揭示涵洞式直立堤的水體交換作用，渦量分布用來說明波浪與防波堤作用機(jī)理，流線分布則是渦量的一個(gè)補(bǔ)充，用來說明涵洞內(nèi)部水體的振蕩規(guī)律。圖5、6表示波浪與結(jié)構(gòu)物作用穩(wěn)定后，不同周期的4個(gè)典型相位下(時(shí)間間隔為Ti/4)雙液相分布、渦量分布以及流線分布。

圖5 一個(gè)周期內(nèi)涵洞內(nèi)部水體交換特性(Ti=2.40 s,Hi=0.08 m,dS=0.177 m,S=0.096 m)Fig. 5 Water exchange characteristics inside the culvert in a period (Ti=2.40 s,Hi=0.08 m,dS=0.177 m,S=0.096 m)

圖5的試驗(yàn)條件為波周期Ti=2.40 s，波高Hi=0.08 m，涵洞深度dS=0.177 m，涵洞高度S=0.096 m。由圖5可看出，自入射波從波谷到達(dá)波峰的半個(gè)周期內(nèi)(圖5(a))，港池內(nèi)側(cè)出現(xiàn)了未染色水體而港池外側(cè)出現(xiàn)部分染色水體，表明涵洞內(nèi)部以及港池內(nèi)側(cè)(透射波方向)水體逐漸被帶入港池外部(來波方向)，期間水體大量摻混，涵洞內(nèi)外水體發(fā)生交換。對(duì)比渦旋軌跡(圖5(b))能發(fā)現(xiàn)，當(dāng)?shù)糖安ɡ颂幱诓ü葧r(shí)，涵洞入口(來波方向)附近會(huì)形成上下兩個(gè)渦，涵洞出口(透射波方向)水面以及水底附近也會(huì)形成渦旋。此后，在波浪從波谷到波峰的運(yùn)動(dòng)過程中，涵洞入口附近的兩個(gè)渦逐步被帶到堤前水面附近，并逐漸消散，同時(shí)帶走了大量涵洞內(nèi)部的水體；而防波堤外側(cè)水面及水底附近的渦逐漸運(yùn)動(dòng)到涵洞出口處，并進(jìn)入涵洞。觀察流線圖能發(fā)現(xiàn)，涵洞內(nèi)部流線較為平順，僅在波浪轉(zhuǎn)向(波峰、波谷)時(shí)內(nèi)部流線發(fā)生改變，表明在波浪與其作用的過程中，整個(gè)水體交換在半個(gè)周期內(nèi)(波峰到波谷，或者相反)是穩(wěn)定而平順的，在整個(gè)周期內(nèi)則是呈現(xiàn)出振蕩的特性。當(dāng)入射波從波峰到達(dá)波谷的半個(gè)周期內(nèi)，涵洞內(nèi)部以及港池外側(cè)水體逐漸被帶入港池內(nèi)部，對(duì)比渦旋軌跡能發(fā)現(xiàn)，此時(shí)渦旋運(yùn)動(dòng)規(guī)律與波谷到波峰時(shí)的規(guī)律相反。堤前水面以及水底形成的渦逐漸運(yùn)動(dòng)到涵洞內(nèi)部，而涵洞出口附近的渦則逐漸被帶入到水面以及水底附近。上述結(jié)果表明，涵洞內(nèi)部水體在波浪的驅(qū)動(dòng)下，往復(fù)環(huán)流，形成一股振蕩的水流，進(jìn)而完成水體交換。在此過程中，雖然單次波浪周期所完成的水體交換量不算太大，如圖5的雙液相可以發(fā)現(xiàn)港池外側(cè)逐漸出現(xiàn)港池內(nèi)測(cè)的水體，但在長(zhǎng)時(shí)間的波浪作用下，水體交換范圍將逐漸擴(kuò)大，進(jìn)而可以較好地實(shí)現(xiàn)港池內(nèi)外水體交換的目的。

圖6試驗(yàn)條件與圖5類似，只是波周期變?yōu)?.18 s，可以看出，其分布規(guī)律與長(zhǎng)周期大體相似，但仍有一定差異，主要體現(xiàn)在水體交換強(qiáng)度弱于長(zhǎng)周期波浪，同時(shí)，涵洞附近的渦脫軌跡也有一些不同。主要由于，在短周期波浪作用下，水體垂向能量部分不均勻，水面附近能量遠(yuǎn)遠(yuǎn)大于水底，故涵洞入口處僅在涵洞上部形成渦，而后逐漸被帶入到水體表面，如此來回振蕩；涵洞出口處與長(zhǎng)波類似，形成上下兩個(gè)渦，但由于透過波能較小，其運(yùn)動(dòng)軌跡遠(yuǎn)遠(yuǎn)小于長(zhǎng)波，僅在涵洞出口附近來回振蕩。波浪作用下涵洞內(nèi)部的流線分布規(guī)律幾乎一致。

圖6 一個(gè)周期內(nèi)涵洞內(nèi)部水體交換特性(Ti=1.18 s,Hi=0.08 m,dS=0.177 m,S=0.096 m)Fig. 6 Water exchange characteristics inside the culvert in a period (Ti=1.18 s,Hi=0.08 m,dS=0.177 m,S=0.096 m)

3.2 涵洞參數(shù)對(duì)水體交換的影響

圖7表示當(dāng)波高Hi=0.02 m，涵洞深度dS=0.177 m時(shí)不同涵洞高度S下的涵洞內(nèi)部振蕩流Qn隨時(shí)間t的變化曲線，圖7(a)、7(b)分別表示波周期Ti=1.00 s、Ti=2.40 s時(shí)的Qn變化曲線。從圖7中可以看出，Qn隨t正弦變化；隨著S的增大，Qn在短周期下隨著高度的增長(zhǎng)近乎同比例增大，在長(zhǎng)周期下，增大倍數(shù)略小于高度增長(zhǎng)率。

圖7 涵洞高度對(duì)涵洞內(nèi)部振蕩流Qn的影響(Hi=0.02 m,dS=0.177 m)Fig. 7 The effects of the culvert heights on the oscillatory flow Qn inside the culvert (Hi=0.02 m,dS=0.177 m)

圖8為當(dāng)波高Hi=0.02 m，涵洞高度S=0.096 m時(shí)不同涵洞深度dS下的涵洞內(nèi)部振蕩流Qn隨時(shí)間t的變化曲線，圖8(a)、8(b)分別表示波周期Ti=1.00 s、Ti=2.40 s時(shí)的Qn變化曲線。從圖8中可以看出，隨著dS的增大，短周期下的Qn不斷變?。婚L(zhǎng)周期下的Qn也在變小，但變化幅度很小。主要由于短周期波浪能量集中于水面附近，越靠近水面，涵洞所能透過的波能越多，進(jìn)而Qn變大；而長(zhǎng)周期下的波浪水質(zhì)點(diǎn)速度整個(gè)垂向分布較為均勻，改變涵洞位置對(duì)Qn的影響不大。

圖8 涵洞深度對(duì)涵洞內(nèi)部振蕩流Qn的影響(Hi=0.02 m,S=0.096 m)Fig. 8 The effects of the culvert depths on the oscillatory flow Qn inside the culvert (Hi=0.02 m,S=0.096 m)

上述結(jié)果表明，增大涵洞高度能夠明顯地提高涵洞內(nèi)部的流量，促進(jìn)水體交換；減小涵洞深度對(duì)短周期波浪能夠一定程度上提高涵洞內(nèi)部的流量，但對(duì)長(zhǎng)周期波浪作用很小。

3.3 波浪參數(shù)對(duì)水體交換的影響

由前文可知，在小波高情況下，涵洞內(nèi)部振蕩流隨時(shí)間呈現(xiàn)出正弦變化，可以將其表示為Qmaxcos(σt)，其中Qmax為振蕩流振幅，σ為波浪圓頻率，t為時(shí)間。故在一個(gè)周期內(nèi)涵洞內(nèi)部水體的總輸移量(周期積分)等于0。所以應(yīng)該考慮的是半個(gè)周期的平均凈流出，而半個(gè)周期的平均凈流出流量Qave與振蕩流振幅存在線性關(guān)系：

(12)

圖9為在波高Hi=0.02 m，涵洞高度dS=0.177 m，涵洞深度S=0.096 m條件下，不同波周期Ti下涵洞內(nèi)部振蕩流頻譜。

圖9 不同波周期下涵洞內(nèi)部振蕩流頻譜(Hi=0.02 m,dS=0.177 m,S=0.096 m)Fig. 9 Spectrum of the oscillatory flow in different wave periods (Hi=0.02 m,dS=0.177 m,S=0.096 m)

可看出，隨著Ti的增大，振蕩流基頻振幅AQ不斷增大；除了基頻振幅外，涵洞內(nèi)部存在部分倍頻項(xiàng)，但倍頻振幅較小，與基頻振幅相比幾乎可忽略。同時(shí)，發(fā)現(xiàn)波周期Ti=2.40 s (kh=0.55,k為波數(shù),h為水深) 所對(duì)應(yīng)的基頻振幅約為Ti=0.80 s (kh=2.55)的7倍。上述結(jié)果表明，波周期對(duì)水體交換有較大影響，長(zhǎng)周期波浪的振蕩流量遠(yuǎn)遠(yuǎn)大于短周期波浪。

圖10 波高對(duì)涵洞內(nèi)部振蕩流的影響(Ti=0.80 s,dS=0.177 m,S=0.096 m)Fig. 10 The effects of wave heights on the oscillatory flow inside the culvert (Ti=0.80 s,dS=0.177 m,S=0.096 m)

圖11為50 s后不同波高下的雙液相水體分布范圍?？梢婋S著波浪波高的增大，水體交換范圍明顯擴(kuò)大；當(dāng)作用時(shí)間足夠時(shí)，水體交換能夠逐漸擴(kuò)散到整個(gè)港區(qū)，進(jìn)而提高水體質(zhì)量。上述結(jié)果表明，涵洞式直立堤在波浪作用下能夠有效地實(shí)現(xiàn)港池內(nèi)外的水體交換。

圖11 50 s后不同波高的水體交換對(duì)比(dS=0.177 m,S=0.096 m)Fig. 11 Comparisons of the water exchange at different wave heights after 50 s (dS=0.177 m,S=0.096 m)

4 結(jié) 語

在VPM-THINC/QQ模型的基礎(chǔ)上，提出一種能夠標(biāo)記涵洞內(nèi)外水體變化的歐拉兩相VOF方法，并結(jié)合物理模型試驗(yàn)，開展了波浪作用下涵洞式直立堤水體交換特性的研究，得到如下結(jié)論：

1) 基于VPM-THINC/QQ算法的雙液相VOF方法能夠較好反映真實(shí)流體的內(nèi)部流動(dòng)，具有較高的精度。

2) 涵洞內(nèi)部水體交換機(jī)理如下：當(dāng)防波堤堤前波浪從波谷向波峰運(yùn)動(dòng)時(shí)，涵洞入口處和涵洞出口的水面以及水底附近形成渦旋，帶動(dòng)涵洞內(nèi)部水體向堤前運(yùn)動(dòng)；當(dāng)?shù)糖安ɡ藦牟ǚ逑虿ü冗\(yùn)動(dòng)時(shí)，涵洞入口的水面附近和涵洞出口處形成渦旋，帶動(dòng)涵洞內(nèi)部水體向堤后運(yùn)動(dòng)。故在波浪作用下，涵洞內(nèi)部會(huì)形成往復(fù)的振蕩流，且該振蕩流主導(dǎo)著港池內(nèi)外的水體交換。

3) 涵洞高度對(duì)涵洞內(nèi)部振蕩流有較大影響；涵洞深度對(duì)短周期波浪有較大影響，而對(duì)長(zhǎng)周期波浪影響不大；隨著波浪非線性的增強(qiáng)或波周期的增大，涵洞內(nèi)外水體交換作用加強(qiáng)。同時(shí)，隨著非線性的增強(qiáng)，涵洞內(nèi)部存在一個(gè)類似于斯托克斯漂移的凈輸移。在波浪長(zhǎng)時(shí)間的作用下，涵洞式直立堤能夠有效地促進(jìn)港池內(nèi)外水體交換。