波浪在潛堤和潛堤地形上的傳播變形

孫亮， 那鑫宇， 陳柯宇， 李芬， 陳明勝

(1.武漢理工大學 交通學院，湖北 武漢 430063； 2.大連理工大學 海岸和近海工程國家重點實驗室，遼寧 大連 116024； 3.河海大學 港口海岸與近海工程學院，江蘇 南京 210098)

波浪是設計近岸結構物需要考慮的主要環境因素，波浪從深水傳播到近岸區域會由于水下結構物或者地形的影響而發生形態的改變。潛堤是一種建造在海岸線外、完全淹沒在水面下的護岸建筑物[1]，其主要功能是阻擋和減弱外海波浪的作用，為港口、防波堤、閘門等提供掩護作用[2]。波浪在潛堤上傳播時，不同頻率的波浪會發生相互作用，產生高階諧波，也可能發生破碎等現象。分析波浪在潛堤和潛堤地形上的傳播變形的方法主要有理論分析、物理模型試驗和數值模擬3種研究方法。理論分析方法對試驗研究及數值模擬具有指導意義，但在解決實際工程問題時具有較大的局限性。物理模型試驗可通過試驗波浪水槽進行準確的現場觀測，得到較為可靠的實驗數據，但成本偏高，流程復雜。隨著數值分析相關理論的發展和計算機硬件性能的提高，數值模擬方法成本低，重復性強，相比物理模型試驗可以節約大量人力、資源和時間。

分析近岸波浪傳播變形的數學模型主要包括：Boussinesq方程、緩坡方程、N-S方程等。Beji等[3]基于一維的Boussinesq方程模擬了非破碎波在潛提上的傳播并與試驗結果相比較。Losada等[4]使用緩坡方程模擬了可滲透梯形潛堤上的波浪演變與附近流場情況。滕愛國等[5]以完整N-S方程為基礎，改進湍流模型，較為精確地模擬出橢圓余弦波在潛堤上的傳播過程。周援衡等[6]同樣地基于N-S方程，采用VOF法追蹤自由液面和湍流模型，對橢圓余弦波在滲透潛堤上的傳播變形以及潛堤前后的波形變化進行數值模擬。劉誠等[7]建立了比前人更為精確的PLIC-VOF數值方法并使用有限差分法求解N-S方程，研究了波浪越過斜坡潛堤和垂直潛堤時的波浪場。張婷婷[8]利用Boussinesq方程模擬了規則波和自由長波在潛堤地形上的傳播變形。吳迪[9]采用Boussinesq方程，模擬了波浪在非滲透潛堤和滲透潛堤上的傳播并進行對比，顯示出滲透作用對波浪的巨大消波作用。馬小舟等[10]利用完全非線性Boussinesq方程數值模型研究了不規則波浪在緩坡潛堤上傳播時非線性參數的變化。廖斌等[11]利用五階WENO有限差分格式以及三階TVD Runge-Kutta格式離散求解RANS方程分析了潛堤上波浪的傳播問題。

本文基于開源計算流體力學工具OpenFOAM建立了二維數值波浪水槽[12]，模擬規則波和橢圓余弦波在潛堤及潛堤地形上的傳播變形并與已發表的結果對比驗證，為進一步分析波浪對海岸結構物的作用提供了堅實基礎。

1 數值波浪水槽的建立

OpenFOAM利用有限體積法求解N-S方程分析可壓縮以及不可壓縮流體問題。對于不可壓縮流體，N-S方程由連續性方程和動量守恒方程組成[13]：

·U=0

(1)

?ρU/?t+·(ρUU)-·(μU)-ρg=

(2)

式中：ρ為流體密度；p為流體壓力；μ為粘滯系數；U為速度矢量；fσ為表面張力。

在建立數值波浪水槽的左端邊界，根據線性波浪理論計算波面、水質點速度和壓力分布從而產生入射波浪。在左端入射邊界和右端邊界都設置了相應的消波海綿層[14]。

2 規則波在潛堤上的傳播

2.1 物理模型試驗及數值波浪水槽的布置

為了驗證數值模型的精度，本文首先使用Gobbi等[15]研究波浪在淹沒潛堤上的傳播問題時所采用的代爾伏特物理模型試驗數據，入射波浪的參數如表1所示，潛堤的形狀和本文對比中所選取的浪高儀位置如圖1所示。

表1 代爾伏特物理模型試驗的入射波浪參數

圖1 代爾伏特物理模型試驗中潛堤與浪高儀的位置(單位：m)

本文建立了如圖2所示的二維數值波浪水槽，長30 m，高0.5 m，水深設為0.4 m。數值波浪水槽的左端為造波邊界，波浪由左向右傳播，右端為出流邊界。在造波邊界后設置1.0 m長的海綿層，避免由地形產生的反射波到達左端造波邊界。在右端出流邊界前設置6.0 m長的海綿層，避免越過潛堤的波浪到達右端邊界在數值水槽內形成反射。對于如圖2所示的計算區域，沿水槽寬度方向只劃分了1個網格。在水槽不同垂直截面處劃分了相同的單元數，這樣在潛堤斜坡上和坡頂網格會自然加密，有利于精確模擬潛提附近區域復雜的物理過程。初始網格的控制尺寸取為0.01 m，每個垂直截面上有50個單元，水平方向有3 000個單元，網格1在計算區域的總單元數為15萬。網格2在每個垂直截面上有100個單元，計算區域單元總數30萬。網格3在每個垂直截面上有200個單元，計算區域單元總數為60萬。

2.2 波面驗證與精度分析

圖3給出了10個斷面處(如圖3所示)的波面隨時間變化的計算和試驗測量值的比較。圖3(a)和(b)中x=2.0 m和x=4.0 m處的結果用于研究入射波浪遇到潛堤之前的變化情況，由圖3(a)和(b)中的對比可以看出基于3種網格的數值模擬結果和測量值吻合較好，數值波浪水槽可以較好地重復試驗中的入射波。

圖2 數值波浪水槽的設置(單位：m)

圖3 代爾伏特試驗中波面測量值與數值結果

圖3(c)～(g)中斷面x=10.5，12.5，13.5，14.5，15.7 m位于潛堤斜坡上和坡頂。在x=10.5 m(位于潛堤前坡)，x=12.5和13.5 m(位于潛堤頂)處基于3種網格的數值結果幾乎沒有區別并且與測量值吻合得較好，表明目前的3種單元劃分都滿足了精度要求。在x=14.5，15.7 m(位于潛堤后坡)處基于3種網格的數值結果出現了差別，網格2和3的結果更加接近測量值。正如Gobbi等[15]所分析所得結論：波浪經過潛堤后高階諧波發生分離使得波浪場變得更加復雜。所以更加精細的網格才能滿足計算精度的要求。圖3(h)～(j)中斷面x=17.3，19.0，21.0 m完全位于潛堤后，基于網格3的數值結果與測量值吻合得更好。

3 橢圓余弦波在潛堤地形上的傳播

3.1 數值波浪水槽的布置

張婷婷[8]進行了規則波在潛提地形上的試驗研究,潛堤地形由橢圓型曲線和指數型曲線組成，其數學表達式為：

z+h=

(3)

式中：h是常水深h=0.5 m。hs是潛堤最高點的高度hs=0.4 m，a1=4.0 m，a2=2.5 m，θ1=arctan(a1z1/hs)，z1=0.125，θ2=arctan(a2z2/hs)，z2=0.2。如圖4所示曲線地形分為4個部分，其中x0=20.2 m，為潛堤最高點對應的橫坐標，此處z+h=hs=0.4 m。x1、x2為潛堤最大坡度處對應的橫坐標，x1=17.077 m，x2=22.152 m。張婷婷[8]將基于Boussinesq方程的數值結果與模型試驗數據比較，驗證了數值結果的可靠性。

張婷婷[8]進一步利用Boussinesq方程模擬了橢圓余弦波在潛堤地形上的傳播變形，分析了如圖4所示x=10.0，16.0，18.0，19.2，20.2，22.5，24.8 m處的波面變化。本研究中所采用的橢圓余弦波的波浪參數如表2所示。

圖4 大連理工大學試驗水槽中的潛堤地形

表2 橢圓余弦波的波浪參數

本文中建立了與圖2類似的總長130 m、高0.6 m的二維數值波浪水槽，左端吸波海綿層長40 m，右端吸波海綿層長45 m。二維數值波浪水槽中潛堤地形頂部的水平坐標由圖4中的20.2 m變為65.2 m，數值浪高儀的位置變為：x=55.0，61.0，63.0，64.2，65.2，67.5，69.8 m。計算區域中每個垂直截面上有300個單元，水平方向有1 300個單元，共計39萬個單元。

3.2 波面驗證與精度分析

圖5給出了數值波浪水槽中7個位置處基于Boussinesq方程和N-S方程的波面計算結果。

圖5 橢圓余弦波在潛堤地形的傳播過程

圖5(a)中x=55.0 m(對應圖4中x=10.0 m)位于潛堤地形前，用于研究入射波浪遇到潛堤之前的變化情況；圖5(b)中斷面x=61.0 m(對應圖4中x=16.0 m)，圖5(c)中x=63.0 m(對應圖4中x=18.0 m)，圖5(d)中x=64.2 m(對應圖4中x=19.2 m)位于潛堤前坡處，由圖5(a)～(d)中的對比可以看出：基于Boussinesq的數值結果和基于N-S方程的數值結果吻合較好。基于Boussinesq的數值結果和基于N-S方程的數值結果吻合較好。圖5(e)中x=65.2 m(對應圖4中x=20.2 m)位于潛堤坡頂，圖5(f)中x=67.5 m(對應圖4中x=22.5 m)和圖5(g)中x=69.8 m(對應圖4中x=24.8 m)位于潛堤后，OpenFOAM數值模擬結果和Boussinesq方程數值模擬結果出現微小差別。

值得注意的是在圖5(g)中斷面x=69.8 m(對應圖4中x=24.8 m)處，在基于Boussinesq方程的數值模擬結果中發現了毛刺波[8]。張婷婷[8]認為：自由長波本身含有高階諧波，波浪越過潛堤后發生反淺水效應，原來被約束的諧波發生分離成為自由波。然而基于OpenFOAM的數值模擬結果中并沒有體現這一現象，需要進一步分析。

4 結論

1)本文建立了基于OpenFOAM的二維數值波浪水槽，分別模擬了規則波在潛堤上的傳播、橢圓余弦波在潛堤地形上的傳播，所得的波面結果與試驗數據、基于Boussinesq方程的數值結果都吻合較為理想。

2)本文的研究驗證了基于OpenFOAM的數值波浪水槽在預測波浪與近岸結構物相互作用時的可靠性，為近岸結構物的相關研究及設計提供可靠依據。