海洋工程數值風-浪-流水池模擬及分析*

劉紅兵 陳國明 呂 濤 黃 翱

(1. 中國石油大學(華東)海洋油氣裝備與安全技術研究中心 山東青島 266580； 2. 哈爾濱工程大學船舶工程學院 黑龍江哈爾濱 150001)

在復雜的海洋環境中，風、浪、流通常是共同存在且相互交替運動，從而產生類似于“Doppler”效應的風-浪-流耦合作用[1]，導致流場結構和波浪幅值、波高、周期等波浪要素發生改變，使得波浪發生折射乃至破碎[2]。因此，準確模擬海洋風-浪-流時空分布特征，探討風-浪-流耦合作用機理以及耦合作用下風-浪-流場分布特征，對于保證海洋工程結構物安全運行具有重要意義[3]。

目前，國內外關于海洋環境水動力性能的研究主要采用理論研究、實尺試驗、模型試驗和數值水池模擬等方法。基于無黏無旋假設的勢流理論可以快速計算波浪傳播及其對海洋工程結構物的作用，卻無法計算由于黏性效應引起的一些物理現象，如風-浪-流耦合遠程傳播和漩渦演變等，具有較大的局限性[4]。實尺試驗法受海洋環境影響因素較大，且試驗成本高并難以控制，在海洋工程應用中難以實現[5]。模型試驗法主要采用一定的縮尺比例轉化關系，通過縮尺模型水動力學性能反映整體結構水動力學性能，在海洋工程中應用較為廣泛，但存在試驗成本高、試驗周期長、試驗縮尺模型和原型之間準確轉換困難等缺陷，使其應用仍受到一定限制[6]。隨著近年計算機技術高速發展，數值波浪水池以其低成本、短周期和高精度等特點，在海洋工程波浪誘導運動和波浪載荷預測等方面得到了廣泛應用[7]。例如，Wang[8]和蘇曉杰 等[9]基于歐拉-拉格朗日算法追蹤流體瞬時自由液面的變化過程，建立了完全非線性的三維數值波浪水池。Alamian等[10]利用活塞運動造波法，在數值波浪水池中模擬了線性波和非線性波；Finnegan等[11]基于CFX軟件，建立了搖板式數值波浪水池；梁修峰 等[12]基于FLUENT軟件建立了一種面向海洋工程應用的二維數值波浪水池，并對非線性Jonswap波浪譜進行了模擬。但這些數值水池研究主要集中于單一波浪運動規律的研究，而對于風-浪-流三者之間相互耦合作用研究較少。本文基于FLUENT軟件，利用UDF二次開發，構建了一種適于三維海洋工程的數值風-浪-流水池，成功模擬了風-浪-流之間的相互耦合作用，并探討了風和海流對波浪的影響。

1 數值風-浪-流水池數學模型

1.1 控制方程

假定數值風-浪-流水池所涉及的流體(空氣和海水)均為不可壓縮黏性流體，則數值風-浪-流水池連續性方程和動量方程分別為

(1)

(2)

波浪傳播過程中，其氣-液交界面時刻發生改變，因而需要準確追蹤自由液面。本文采用VOF方法捕捉自由液面，即通過引入單元內各流體體積分數aq來捕捉氣液兩相自由液面。當aq=0時，表示該網格單元內沒有q相流體;當aq=1時，表示該網格單元內充滿q相流體;當0

(3)

式(3)中：aq為q相流體在網格內所占體積分數，當q=1、2時，aq分別表示空氣、水相體積分數。

1.2 造波方法

數值造波方法對于數值水池準確模擬波浪傳播特性具有重要作用，現階段主要的數值造波方法有仿物理造波法(推板造波法和搖板造波法)[13]、邊界造波法[14]和源項造波法[15]等,其中仿物理造波法原理簡單，造波功能完善，但對計算區域內網格質量要求較高，計算時間長；源項造波法產生的波浪效果好，波形穩定，但計算量大，且不易收斂；邊界造波法通過在波浪入射邊界給定波浪水質點速度和波面高度，從而生成數值波浪，具有較好的適用性并在海洋工程波浪領域應用較為廣泛。因此，采用邊界造波法將風-浪-流水池入口設置為風-浪-流相互作用的起始位置，各速度場方程可由Airy線性波浪和梯度風場理論描述，其中波面方程為

(4)

式(4)中：η為波面高度，m；H為波高，m；k為在2π相位角內波的個數，rad/m；w為2π相位角內波的振動次數，rad/s；x為X方向位移,m。

風-浪-流水池入口速度包括自由液面以下水質點運動速度和自由液面以上空氣運動速度，其中入口水質點速度由波浪水質點速度和水流流動速度兩部分組成，可根據波流共同作用下流速場理論進行疊加[16]，即

ux=uwx+uc=

(5)

(6)

式(5)、(6)中：ux、uz分別為水質點X和Z方向速度，m/s；uwx為波浪水質點速度，m/s；z表示Z方向上的位移，m；uc為水流速度，m/s;T為波浪周期，s；d為水深，m。

自由液面之上的空氣流動速度可采用指數梯度風進行定義[9]，即

(7)

1.3 消波方法

采用數值風-浪-流水池進行模擬時，當波浪傳播至水池末端時，波浪會產生反射，對前進波浪產生干擾。因此，為準確獲得高質量、高精度的波浪，需對反射波進行消波處理。常用的消波方法主要有多孔介質消波法和動量阻尼消波法，其中多孔介質消波法通過在消波段設置一個多孔介質區域，實現消波目的，但關于多孔介質模型的選取及參數的設置較為復雜；動量阻尼消波法通過在消波段的流體動量方程中添加一個黏性阻尼衰減項，從而實現消波目的。消波段動量方程[10]為

(8)

(9)

式(8)、(9)中：θ為動量阻尼消波系數；x0和x1分別為消波段前端和尾端X方向坐標，m；θ1為動量阻尼消波經驗系數，取值10。

2 數值風-浪-流水池模擬及分析

2.1 網格模型及參數設置

采用GAMBIT建立三維數值風-浪-流水池幾何模型并劃分網格，其整體尺寸為：長50 m、寬20 m、高15 m、水深為10 m、右端消波段長度為8 m，如圖1所示。整個數值風-浪-流水池網格采用漸進式六面體網格，越靠近水池底端、右端和上端邊界處網格越稀疏，同時對自由液面附近兩倍波高區域進行加密處理，便于更好地捕捉波形變化。

圖1 數值風-浪-流水池幾何模型網格劃分

數值風-浪-流水池左邊界設置為速度入口，氣液兩相速度值和體積組分數可根據氣液兩相運動速度場和波面運動方程并采用UDF進行定義;右邊界設置為自由出流邊界，上邊界設置為對稱邊界，其余邊界為固壁邊界。造波過程中，采用三維非定常分離隱式求解器求解，其中湍流模型為RNGk-ε模型，壓力速度之間的耦合作用采用PISO算法，離散格式為中心差分格式，壓力方程采用加權體積力格式，動量方程、湍動能和湍動能耗散率方程均采用二階迎風格式。壓力參考值為一個標準大氣壓，重力加速度為-9.81 m/s2，時間步長為0.01 s。

2.2 模型驗證

為了驗證數值風-浪-流水池的準確性和適用性，采用線性Airy波浪對其進行驗證。設定模擬的線性Airy波浪參數為：波長L=8 m，波高H=0.40 m，波陡λ=0.05，周期T=2.263 6 s，波速c=3.534 2 m/s。圖2為該線性Airy波浪生成和傳播過程示意圖，上部藍色部分為空氣，下部紅色部分為水，中間綠色部分為波面，可以看出t=20 s時數值水池已形成穩定的Airy波面，同時水池右端采用動量阻尼消波效果較好，消波段波面基本保持不變。

圖2 Airy波浪生成和傳播過程

分別在數值風-浪-流水池中不同位置處設置浪高儀，實時監測不同位置處波面時程。x=2 m和x=18 m位置處波面時程如圖3所示,可以看出，兩個位置處波面時程均較為穩定，且與理論波面時程變化規律基本一致,其中x=2 m處數值波面時程與理論波面時程吻合度較高,而x=18 m處數值波面時程峰值略小于理論波面峰值,這主要是由于波浪傳播過程中由于黏滯耗能等影響而導致波浪最大幅值呈現一定程度的衰減。

圖3 波浪數值波形和理論波形比較

進一步統計該數值模擬波浪的波長、波幅和周期等特征，并與理論值進行比較分析(表1)，可以看出：數值波長平均值為8.301 m，與理論波長相對誤差(取絕對值，余同)為3.76%；x=2 m處，模擬波幅為0.196 m，與理論波幅相對誤差為2.00%；隨著x取值逐漸遠離造波邊界，模擬波幅與理論波幅相對誤差絕對值逐漸增大，這主要是由于波浪傳播過程中會存在逐步的衰減作用，但波幅相對誤差均小于10%，且不同位置處波浪周期與理論周期相對誤差均小于1%，充分驗證了該數值風-浪-流水池模擬波浪的可靠性。

表1 波長、波幅和周期模擬值與理論值相對誤差統計

2.3 浪-流耦合作用模擬及分析

根據式(5)～(10)，分別設定浪-流耦合作用下造波和消波方法，其中速度入口邊界通過UDF進行定義(包括流速場和風場)，消波段通過增加UDF添加動量阻尼源項進行消波，然后利用數值風-浪-流水池模擬浪-流耦合作用。為了便于比較，設定浪-流耦合作用幾何模型、邊界條件、物理屬性以及波浪參數等均與模型驗證時保持一致；海流速度共設置4組，分別為uc=0 、-0.1 、 0.1 、0.3 m/s。

圖4為4種不同海流流速下x=2 m處波面時程曲線，可以看出，不同流速下波面時程曲線形狀基本保持一致，但流速越大波面時程形成穩定波形所需時間越短，即uc=0.3 m/s時波形最先達到穩定，uc=-0.1 m/s時波浪最后達到穩定狀態。這主要是由于波流同向時，海流對波浪傳播具有推動作用，導致波浪傳播速度加快；而波流逆向時，海流對波浪傳播具有阻礙作用。

圖4 浪-流耦合作用下x=2 m處波面時程對比

分別提取4種不同流速下的波長和波幅，并同波流耦合作用下波長和波幅理論計算公式[16]進行比較分析，結果見表2。由表2知，4種流速下數值波長和波幅與理論值吻合較好，其中數值波長和理論波長最大相對誤差為3.76%，數值波幅和理論波幅最大相對誤差為5.35%。對比4種流速下波浪峰值大小可知，uc=-0.1 m/s時波浪峰值最大，uc=0.3 m/s時波浪峰值最小，這主要是由于波流同向時，波浪傳播速度增大，波形被拉長，導致波浪峰值減小；而波流逆向時，波浪傳播速度減小，波形被壓縮，導致波浪峰值增大所致。

表2 浪-流耦合作用下波長和波幅模擬值與理論值

2.4 風-浪-流耦合作用模擬及分析

實際海洋環境中通常風、浪、流三者是共同存在的，進一步基于該數值風-浪-流水池模擬風-浪-流耦合作用。根據式(7)在造波邊界自由液面上方空氣區域施加梯度風風場，同時海流速度uc取值0.1 m/s，波浪及其他物理參數與上述完全相同。圖5為考慮波浪、浪-流耦合以及風-浪-流耦合等3種情況下，x=2 m位置處波面時程對比，可以看出，風-浪-流耦合作用下波面達到穩定時間和浪-流耦合作用下達到穩定時間差不多，但比單純波浪達到穩定時間快；風-浪-流耦合作用下波面峰值附近呈現一定的傾角，這主要是由于海風使得波面峰值附近波浪運動速度加快，導致波峰發生傾斜，加速波浪發生碎浪現象。進一步統計風-浪-流耦合作用下波浪波長和波幅，其平均值分別為8.836 m和0.185 m，相比于浪-流耦合作用下的波長和波幅分別增大2.45%和減小3.14%。

圖5 風-浪-流耦合作用下x=2 m處波面時程對比

3 結論

1) 基于FLUENT軟件，通過UDF二次開發，構建了適于海洋工程的數值風-浪-流水池，實現了數值風-浪-流水池造波和消波功能，成功模擬了浪-流以及風-浪-流相互耦合作用，模擬結果與理論結果吻合度較好，驗證了該數值風-浪-流水池具有較好的精度。

2) 浪-流耦合作用數值模擬結果表明，海流對于波浪運動參數具有較大影響。波流同向時，海流導致波浪傳播速度和波長增大，幅值減小；波流逆向時，波浪傳播速度和波長均減小，波幅增大，且海流速度越大，波浪參數增大或減小的幅度越明顯。

3) 風-浪-流耦合作用數值模擬結果表明，風會導致波浪傳播過程中波長增大，波幅減小，同時會導致波峰呈現一定傾角，加速波浪發生破碎。