基于OpenFOAM的沉箱與梳式防波堤水動力特性比較*

萬 祥，李雪艷，謝 天，孫家文，解曉敏，程 志

(1.魯東大學(xué)，海岸研究所&山東省海上航天裝備技術(shù)創(chuàng)新中心，山東 煙臺 264025；2.國際海洋環(huán)境監(jiān)測中心，遼寧 大連 116085)

國家大力發(fā)展海洋強國戰(zhàn)略對海洋防護結(jié)構(gòu)的形式、功能等均提出了新的要求。因傳統(tǒng)重力式防波堤存在材料利用率低、水體交換能力差、使用功能單一、不利于海洋生態(tài)保護等缺點，功能多樣性的新型防波堤結(jié)構(gòu)成為研究熱點。

梳式防波堤是在直立沉箱防波堤的基礎(chǔ)上，將沉箱主體按適當比例取出并代之以沉箱翼板，在迎浪面上形成一種梳齒狀的防波堤形式[1]。方箱梳式防波堤具有沉箱主體和空腔處翼板兩個迎浪面，波浪作用在防波堤上時，會存在一個相位差，可有效減少最大波浪力[2]。方箱梳式防波堤已有較多的相關(guān)研究和工程使用[3]。Sun等[4]通過物理模型試驗，論證了方箱梳式防波堤能夠減小結(jié)構(gòu)所受波浪力的水力機制，給出了該結(jié)構(gòu)反射系數(shù)的經(jīng)驗公式。趙軍[5]通過物理模型試驗對方箱梳式透空式防波堤結(jié)構(gòu)進行研究，給出了不同波浪要素條件下防波堤的透射系數(shù)。房卓等[6-7]對波浪作用下的透空與非透空方箱梳式防波堤結(jié)構(gòu)進行了數(shù)值模型與物理試驗研究，探討波浪力特性和波浪反射系數(shù)。傅朝方等[8]在分析方箱梳式防波堤非線性耗能機理基礎(chǔ)上，采用阻抗分析方法給出堤前波浪反射系數(shù)的經(jīng)驗公式。朱大同[9]指出方箱梳式防波堤是一種兼有消波和聚能的海岸結(jié)構(gòu)，波浪進入空腔時局部水體受到的流動消減。王心玉等[10]基于線性勢流理論，建立了波浪與透空型方箱梳式防波堤相互作用的理論解，計算得到方箱梳式防波堤的反射和透射系數(shù)。目前方箱梳式防波堤已在大連、秦皇島等港區(qū)得到廣泛應(yīng)用。

沉箱式防波堤的消浪性能[11]、受力特性[12]及其周圍的流場特征[13]已有較多研究，而梳式防波堤的相關(guān)特性仍缺乏系統(tǒng)的分析，尤其是在相同條件下梳式透空堤與傳統(tǒng)沉箱防波堤水動力特性的差異方面。本文基于OpenFOAM構(gòu)建波浪與沉箱及梳式防波堤作用的數(shù)值模型，分析比較不同周期作用下兩者對波浪的反射系數(shù)、所受波浪壓強與結(jié)構(gòu)周圍流場分布特征的異同，以期為工程實際提供一定參考。

1 數(shù)值模型建立

1.1 控制方程

基于OpenFOAM建立波浪與防波堤作用數(shù)值模型所采用的流體控制方程為連續(xù)方程和雷諾平均納微斯托克斯方程?？刂品匠瘫磉_式如下：

(1)

Cκα-gXρ-Pρgh+S*

(2)

式中：U為速度；t為時間；ρ為流體密度；μ為動力黏度；ρτ為雷諾應(yīng)力項；τ為因脈動值引起的雷諾應(yīng)力張量；Cκα為表面張力項，C為表面張力系數(shù)，取0.07 kg/s2，κ為自由面的曲率，α為流體體積分數(shù)；g為重力加速度；X為位置矢量；S*為阻尼消波項；Pρgh為壓力項；P為總壓強，當結(jié)構(gòu)物壓力點位于自由表面時代表動水壓強，其表達式為：

Pρgh=P-ρgX

(3)

1.2 邊界條件

1.2.1造波邊界

根據(jù)線性波浪理論，得到規(guī)則波波高與推板沖程之間的關(guān)系：

(4)

式中：Hs為波數(shù)傳遞函數(shù)；H為規(guī)則波波高；S為造波板沖程；k為波數(shù)；d為水深。假定造波板做小振幅簡諧運動，則求解的波面方程為：

(5)

式中：ω為目標波浪的圓頻率；Hsn為波數(shù)傳遞函數(shù)。

數(shù)值造波：基于OpenFOAM軟件采用主動吸收式造波技術(shù)，為了消除二次反射波影響，在推板原有位移的基礎(chǔ)上附加修正位移，該位移產(chǎn)生一個與反射波相位相反、幅值相等的波浪，抵消掉反射波，以達到吸收二次反射目的[14]。設(shè)推板原有位移為x0，修正位移為xa，主動吸收造波機推板的位移x和速度即u為：

(6)

(7)

式中：Sa為修正位移xa對應(yīng)的沖程；φ為相位，類比式(5)得到修正位移xa產(chǎn)生的波面ηa：

(8)

式中：φ為推板造波生產(chǎn)的波浪場速度勢。每個時刻的實測波面ηm由4部分組成：原有位移対應(yīng)的波面η、修正位移対應(yīng)的波面ηa、反射波対應(yīng)的波面ηr以及二次反射波対應(yīng)的波面ηrr。

ηm=η+ηa+ηr+Krrηr

(9)

(10)

將式(5)(6)(8)(10)代入式(9)并移項，得到：

(11)

1.2.2消波邊界

采用阻尼層消波法，即在動量方程中添加消波源項，使得消波區(qū)波浪的速度逐漸衰減，其S*表達式為：

S*=ρUχ

(12)

(13)

式中：χ為衰減函數(shù)；χ1為衰減系數(shù)；x為笛卡爾坐標系下的橫坐標；x0為阻尼消波區(qū)起點的橫坐標；Ls為阻尼消波區(qū)的長度，阻尼區(qū)長度宜為1～2倍波長，故本文衰減系數(shù)取5。

1.2.3自由表面條件

應(yīng)用流體體積函數(shù)法(VOF)追蹤流體自由表面，自由表面處以混合流體處理，可表示為：

(14)

式中：ρw和ρa分別為水和空氣的密度；μw和μa分別為水和空氣的動力黏度。

1.2.4固壁邊界條件

數(shù)值水槽的側(cè)壁和底部均采用不可滑移邊界條件，速度和壓力需滿足下列方程：

U|wall=0

(15)

(16)

2 數(shù)值模型驗證

依據(jù)模型尺寸，建立長30 m、寬0.4 m、高0.8 m的數(shù)值水槽，空水槽數(shù)值模型見圖1。利用此模型進行水深h=0.45 m、波高H=0.1 m、周期T=1.4 s的水槽數(shù)值驗證，合理布置水槽波浪運動曲線的浪高儀，即1#～4#。水槽方向自左向右，左邊為造波區(qū)，取x=21 m；右邊為消波區(qū)，取9 m。

圖1 空水槽數(shù)值模型

依據(jù)OpenFOAM數(shù)值模擬的波面歷時曲線與理論解[15]數(shù)據(jù)結(jié)果進行對比，見圖2。在消波段前，1#、2#兩條浪高儀數(shù)值模擬的波面歷時曲線波高幾乎與理論波一致；在消波段，3#波面歷時曲線有明顯減小、4#波面歷時曲線趨于水平線，說明數(shù)值水槽的消浪區(qū)消波效果顯著，即本文所建立數(shù)值水槽模型可用于有結(jié)構(gòu)物的水槽數(shù)值計算。

圖2 波面歷時曲線

3 計算參數(shù)

3.1 結(jié)構(gòu)模型

本文設(shè)計的水槽中結(jié)構(gòu)物是梳式與沉箱兩種防波堤，沉箱與梳式結(jié)構(gòu)物長為0.4 m、高為0.6 m、寬為0.4 m，但梳式防波堤從中心對稱處摳除0.35 m×0.60 m×0.10 m及0.05 m×0.15 m×0.10 m長方體。采用結(jié)構(gòu)化網(wǎng)格對計算區(qū)域進行剖分，在自由表面處和防波堤結(jié)構(gòu)附近進行網(wǎng)格加密，以獲取更精確的數(shù)值結(jié)果。水平方向為x軸，沿波浪傳播為正；垂直為z軸，向上為正；經(jīng)過加密數(shù)值水槽2次，梳式與沉箱兩種結(jié)構(gòu)局部網(wǎng)格剖分三視圖見圖3。

圖3 沉箱與梳式結(jié)構(gòu)物局部網(wǎng)格剖分三視圖(單位：m)

數(shù)值水槽中沉箱防波堤的壓力測點主要分布于迎浪方向上的正面中心對稱處，從上到下依次布置①～⑥號測點；梳式防波堤壓力測點分別沿迎浪向正面背板與梳齒的中心對稱處，從上到下依次布置①～⑥號與⑦～號測點。在沉箱防波堤迎浪面與梳式防波堤0.05 m厚的背板及梳式防波堤梳齒上沿水深方向從箱體高度為Hw=0.48 m處向下每隔0.06 m布置1個壓力測點，依次布置6個。結(jié)構(gòu)及壓力測點布置見圖4。

圖4 防波堤結(jié)構(gòu)及壓力測點布置(單位：m)

3.2 數(shù)值水槽

基于波浪與防波堤作用的水槽模型，在水深h=0.45 m，波高為H=0.1 m，周期T分別為1.2、1.4、1.6、1.8、2.0 s的規(guī)則波下進行沉箱與梳式防波堤的水動力特性研究對比。以沉箱防波堤為例的數(shù)值水槽模型見圖5。數(shù)值模型按照Goda兩點法[16]布置浪高儀位置，其他波浪要素不變的情況下，選取周期T=1.2 s，兩根浪高儀固定位置取0.6 m，為了滿足數(shù)值模擬要求，將沉箱防波堤放置于距離數(shù)值水槽造波邊界15 m處，進行沉箱防波堤數(shù)值水槽模擬。

圖5 沉箱防波堤數(shù)值模型(單位：m)

4 數(shù)值結(jié)果分析與討論

4.1 波面高程

水深h=0.45 m，周期T=1.2、1.6 s，波高H=0.1 m情況下，沉箱與梳式防波堤迎浪向1#與2#測點處波面歷時曲線見圖6。結(jié)果表明，兩種防波堤不同測點處的波面高程差別較為顯著，在迎浪側(cè)1#測點處，梳式防波堤波面高程低于沉箱防波堤；在2#測點處，梳式防波堤波面高程高于沉箱防波堤。

圖6 沉箱與梳式防波堤迎浪面的波面歷時曲線比較

4.2 消浪特性

在不同周期下(T=1.2、1.4、1.6、1.8、2.0 s)，波高H=0.1 m、波長為L、堤寬為W，波浪反射系數(shù)隨相對堤寬W/L(0.10、0.12、0.14、0.16、0.20)的變化見圖7。根據(jù)沉箱與梳式防波堤迎浪向1#與2#浪高儀監(jiān)測的波高歷時曲線，由Goda兩點法分離出反射波高Hr與入射波高Hi，兩者的比值即為所計算的反射系數(shù)Kr。

圖7 H=0.1 m反射系數(shù)隨相對堤寬變化

結(jié)果表明：在其他波浪要素相同的情況下，本次數(shù)值試驗中，沉箱防波堤對波浪的反射系數(shù)大于梳式防波堤。隨著相對堤寬的增加，沉箱與梳式防波堤的反射系數(shù)隨之增大。沉箱防波堤的反射系數(shù)最大可達0.864，梳式防波堤的反射系數(shù)最大可達0.683。沉箱防波堤對波浪的反射較大，梳式防波堤因中間部位存在空腔，減弱了波浪的反射。

4.3 波浪壓強

在不同周期下(T=1.2、1.4、1.6、1.8、2.0 s)，波高H=0.1 m時沉箱與梳式防波堤迎浪面不同測點處的壓強包絡(luò)圖見圖8，P代表各工況下兩種防波堤上的總壓強，無量綱化處理為P/(ρgh)?？梢钥闯觯?)靜水面以上，梳式防波堤背板①測點處的正壓均大于其梳齒及沉箱迎浪面對應(yīng)位置處的正壓。靜水面附近，梳式防波堤背板②號測點處的正壓較梳式防波堤背板③～⑥、梳齒⑧～及沉箱防波堤迎浪面②～⑥測點處的正壓都大，最大值達0.931 kPa；其負壓與梳式防波堤梳齒⑧測點及沉箱防波堤迎浪面②測點處的負壓相等，且都小于其它測點處的負壓，最小值為0.304 kPa。2)靜水面以下，梳式防波堤梳齒⑧～及沉箱防波堤迎浪面③～⑥測點處的正負壓強均隨相對堤寬的增大而增大，且均呈近似對稱分布；梳式防波堤背板③～⑥測點處的正負壓強值，均隨水深的增大而減小。

圖8 各測點的壓強隨相對堤寬變化的包絡(luò)圖

4.4 流場

水深h=0.45 m、周期T=1.2 s、相對波高H=0.1 m時，梳式防波堤背板中心位置處的二維剖面流場與沉箱防波堤對應(yīng)位置處的流場比較見圖9。t=1/4T時，梳式防波堤空腔內(nèi)水體運動速度較大，且伴有漩渦產(chǎn)生；沉箱防波堤相應(yīng)位置處的水體運動速度較??；t=1/2T時，梳式防波堤空腔內(nèi)水體向開口外側(cè)持續(xù)流動，產(chǎn)生更加顯著的漩渦；沉箱防波堤相應(yīng)位置處的水體運動方向與波浪傳播方向相反；t=3/4T時，梳式防波堤空腔內(nèi)水體處在波谷位置，其水體運動方向同波浪傳播方向相反，開口處水體碰撞劇烈，有明顯的紊動現(xiàn)象；沉箱防波堤相應(yīng)位置處的水體出現(xiàn)回流現(xiàn)象；t=T時，梳式防波堤空腔內(nèi)水體運動產(chǎn)生的漩渦明顯減弱；沉箱防波堤相應(yīng)位置處的水體運動方向與波浪傳播方向相向。

圖9 梳式與沉箱防波堤一個周期內(nèi)流場變化的對比

5 結(jié)論

1)波浪要素相同的情況下，隨著相對堤寬的增加，沉箱與梳式防波堤的反射系數(shù)隨之增大，且沉箱防波堤對波浪的反射較大，梳式防波堤因中間部位存在空腔減弱了波浪的反射。

2)沉箱與梳式防波堤不同測點處所受的波浪正、負壓強值大小差異顯著。靜水面附近，梳式防波堤②測點處的正壓強較其它測點處壓強均大；靜水面以下，梳式防波堤梳齒及沉箱防波堤迎浪面測點處的正、負壓強均隨相對堤寬的增大而增大，且均呈近似對稱分布。

3)在一個波浪周期內(nèi)，梳式防波堤對其周圍流場影響顯著。梳式防波堤背板下方的開口處，形成較大漩渦，水體紊流運動劇烈，加快傳播迎浪面部分水體的能量。