單樁波浪破碎荷載的數值模擬研究

王修亭,張慶河,陳同慶,管寧

（1.天津大學水利工程仿真與安全國家重點實驗室，天津300072；2.交通運輸部天津水運工程科學研究所，天津300456）

單樁波浪破碎荷載的數值模擬研究

王修亭1,張慶河1,陳同慶1,管寧2

（1.天津大學水利工程仿真與安全國家重點實驗室，天津300072；2.交通運輸部天津水運工程科學研究所，天津300456）

基于開源程序OpenFOAM建立了三維波浪破碎數值模型，模型控制方程采用三維不可壓縮粘性流體運動RANS方程，應用VOF方法捕捉自由表面，分別運用動量分布源方法、阻尼消波方法實現數值造波和消波，紊流模型選用RNGk?ε模型。利用所建立的三維波浪破碎數值模型模擬了破碎波對單樁結構物的作用力，并與物理模型實驗數據進行對比分析。結果證明該模型可以較好地模擬近岸單樁結構物所受破碎波浪荷載。

動量分布源造波；破碎波；波浪荷載；OpenFOAM；VOF

海洋石油開采平臺、海上風力發電機、高樁碼頭等近海結構物經常采用樁式基礎。樁式基礎所受到的波浪荷載，特別是近岸波浪破碎荷載對樁式基礎的作用是工程設計所關心的重要問題。

國內外早期對于破浪破碎荷載的研究主要采用物理模型試驗手段。Szwaragi和Nochino[1]、Apelt和Piore?wicz[2]、Swift[3]、Kyte和Torum[4]、馬汝建[5]等對于直立樁柱的波浪破碎荷載進行了一系列研究，基于實驗數據提出了破碎波浪荷載經驗計算公式。

近年來，隨著高性能計算機與并行計算技術發展迅速，數值模擬技術在波浪破碎研究中得到越來越廣泛的應用。Miyata[6]首次提出單相流模型研究破波區的流場特性，Lin和Liu[7]通過在單相流模型中引入非線性k?ε紊流模型，研究了崩破波的二維流場結構，Christensen[8]采用大渦模擬技術，應用MAC和VOF方法追蹤自由表面，揭示了崩破波、卷破波的一些重要的流場特性，Lubin等[9]提出兩相流模型研究波浪破碎問題，Hieu和Tanimoto[10]建立了基于VOF方法的兩相流模型，模擬了波浪與結構物的相互作用，以上數值模擬工作主要針對波浪破碎形態、破波帶內紊流特性和波浪破碎指標等方面進行了研究，而對結構物波浪破碎荷載進行數值模擬的研究還較為少見。Wu等[11]應用VOF方法模擬了破碎波對直立墻面的作用力；Zhang等[12]采用邊界元（BEM）方法研究了二維卷破波對直立固壁的作用力，采用破波指標來描述最大破碎波浪荷載；Mokrani等[13]應用VOF方法，通過求解Navier?Stokes方程模擬了卷破波對直立墻體的作用；Bredmose和Jacob?sen[14]基于開源程序OpenFOAM，采用波群聚焦技術生成極端大波，模擬了海上風機單樁基礎的破碎波浪荷載，但是模型中并沒有引入紊流模型，無法模擬波浪破碎過程中強烈的紊動效應，數值模型的模擬結果也沒有與物理實驗數據對比驗證，僅與Morison公式估算值進行了比較，而Morison公式在淺海波浪破碎條件下的適用性還有待驗證。本文將基于開源程序OpenFOAM建立三維波浪破碎數值模型，對直立樁柱所受的波浪破碎荷載進行數值模擬研究。

1 波浪破碎數值模型的建立

圖1 破碎波浪數值模型水槽示意圖Fig.1 Sketch of numerical wave tank for breaking wavemodel

本文基于OpenFOAM中的兩相流求解器in?terFoam，建立三維波浪破碎數值模型。通過求解RANS方程，模擬水體和空氣兩相不可壓縮流體運動，水體和空氣交界面通過VOF方法確定，紊流模型選取RNGk?ε模型。

造波源分布在一定寬度的區域內，水槽兩端設置消波區域，如圖1所示。

1.1 控制方程

三維波浪破碎數值模型通過在RANS方程的動量方程中分別添加動量分布源項、阻尼項，實現數值造波和消波[15]。模型控制方程采用有限體積法離散，變量均定義在單元中心。模型控制方程如下

式中：γ為流體體積函數，是一個定義在每個單元中心的標量，取值范圍在0～1之間：等于0表示該單元完全被空氣充滿，等于1則表示單元完全被水占據，介于0～1表示該單元內同時存在空氣和水。ρS為造波源項；ρrU為消波項，U、Ur分別為速度矢量和相對速度矢量，r為阻尼消波系數；μeff為考慮了分子動力粘性和紊動作用的有效黏滯系數；p_rgh為修正壓力，p_rgh=p-ρg·X，X為位置矢量。坐標系統為笛卡爾坐標系，波浪沿x軸方向傳播，垂向為y軸，水槽寬度方向為z軸。S和r表達式[15-16]如下

式中：β為源寬系數；ω為目標波圓頻率；D為源強幅值；k=(kx,ky,kz)為波數矢量；θ為波浪入射角；xs、xe分別是阻尼消波區起點和終點的x坐標；h為水槽靜水深；αc、M是待定系數，推薦值分別為20和2[16]。

1.2 邊界條件和初始條件

水槽底部設為滑移邊界；頂部位于空氣中，采用壓力出口邊界條件，設定為標準大氣壓；水槽入口、出口為零梯度條件，即速度、壓力和體積函數等變量在邊界處梯度為零；三維水槽的兩個側邊界為對稱邊界。

計算初始條件設定全場速度和動壓為零。

1.3 計算時空步長

空間與時間步長對波浪模擬結果的精度與計算速度有較大影響[17-18]，在本文的計算中，時間步長與空間步長采用管寧[17]的推薦值：時間步長Δt=T/100；縱向空間步長Δx=L/100；垂向波面附近一個波高范圍內取10個網格，即Δy=H/10，向兩側逐漸稀疏；橫向網格中間密兩邊疏。

1.4 波浪破碎

波浪破碎屬于強非線性的波浪運動，水體紊動強烈，自由表面變化劇烈。本文建立的三維波浪破碎數值模型中，通過引入RNGk?ε紊流模型使三維不可壓縮粘性流體控制方程封閉，模擬水體的紊動；同時應用VOF方法精確捕捉自由表面。為進一步提高模擬精度，劃分網格時在自由表面附近、單樁壁面附近加密網格。RNGk?ε紊流模型的控制方程如下

式中：k、ε分別為紊動動能和紊動動能耗散率；μ、μturb分別為水體動力黏滯系數和紊動黏滯系數；τij、Sij分別為應力張力和變形張量；δij為克羅內克爾符號；η為應變參數，模型中其他經驗參數取值如下：Cμ=0.084 5，αk=αε=1.39，C1ε=1.42，C2ε=1.68，η0=4.38，β=0.012。

2 波浪破碎荷載數值模擬的驗證

本文通過Apelt[2]的物理模型實驗數據對波浪破碎荷載的數值模擬結果進行驗證。物理模型實驗水槽長15m，寬3m，高0.5m，水深0.35m，水槽末端底部坡度為1∶15，實驗樁徑為0.102m，實驗深水波波高為0.066m，周期0.94 s，波長1.291m，波浪破碎類型為崩破波。

2.1 數值模型驗證水槽設置

波浪破碎數值模型的水槽設置如圖2所示。動量分布源和阻尼消波區占用一定范圍的計算域，為保證數值模型有效計算域與物理模型實驗吻合，在數值模型水槽入口、出口分別增加動量分布源和阻尼消波區的長度，數值水槽總長24m。沿水槽中軸線布置12個波高測點，用于測定入射波高和確定破波點位置，其平面坐標見表1。

算例模擬時間為100 s。圖3為斜坡下部平臺上3個測點（#1、#2、#3）的波面歷時曲線，測點#1、#2、#3處的入射波高均可達到0.066m。圖4為破波點附近測點（#7、#8、#9）的波面歷時曲線，測點# 7～#8之間波高迅速減小，表明破波點的x坐標在19～19.2m。圖5為88.065 s時刻斜坡上的波面空間分布曲線，波浪向19m處傳播過程中非線性逐漸增強，波峰變陡，波谷變坦，隨后發生破碎，由波面空間分布曲線推測破波點的x坐標在19～19.5m范圍內。綜合分析，破波點的x坐標位于19～19.2m，計算取值為19.1m。

圖2 破碎波數值水槽示意圖(單位∶m)Fig.2 Sketch of numerical tank for breaking wave(unit∶m)

表1 驗證算例波高測點平面位置Tab.1 Horizontal coordinate of probes in verification simulationm

圖3 斜坡下部平臺波高測點（#1、#2、#3）波面歷時曲線Fig.3 Surface elevation for probes(#1、#2、#3)at the platform downside the slope

2.2 計算網格劃分

單樁底面中心坐標為（19.1，0.24，1.0），直徑0.102m。網格均采用非結構化的六面體單元。為準確模擬單樁對波浪場的影響，在單樁的壁面附近加密網格；自由表面附近加密網格，以提高自由表面的模擬精度。水槽兩端的阻尼消波區內，由計算域向水槽兩端網格稀疏漸變，以增加數值耗散，提高消波效果。

2.3 模擬結果分析

圖6為三維波浪破碎數值模型模擬出的單樁周圍波浪破碎過程。由圖6可以看出，波浪向單樁傳播過程中非線性逐漸增強，波峰逐漸變陡，樁前發生破碎增水，樁后破碎波面擾動劇烈且波高迅速減小；圖6依次為波浪破碎過程中4個時刻的自由水面，波浪首先在波峰頂部出現浪花并向波浪前沿面擴散聚集，波浪在樁前發生破碎，水體撞擊樁前產生增水現象，破碎波面紊亂、變化劇烈并伴有氣液摻混現象，波浪破碎后繼續向樁后傳播，樁前逐漸恢復波浪破碎前的狀態。從圖6可以看出，模擬出的破碎波符合崩破波的特點，證明本模型可以準確模擬波浪破碎前后的非線性變化過程和自由表面的復雜變化。

圖4 破波點附近波高測點（#7、#8、#9）波面歷時曲線Fig.4 Surface elevation for probes(#7、#8、#9)around the breaking point on the slope

圖5 88.065 s波面空間分布曲線Fig.5 Spatial distribution of surface elevation at 88.065 s

圖6 單樁周圍波浪破碎過程Fig.6 Snapshots of wave break around the pile

圖7 單樁破碎波浪力包絡曲線Fig.7 Envelope of breaking wave forces acting on the pile

圖7為波浪破碎數值模型模擬出的單樁破碎波浪力x方向歷時曲線及其包絡線。受波浪破碎水流的影響，實際上樁柱周圍的波浪不是完全恒定的，因此破碎波波浪力也是變動的。x方向破碎波浪力極大值范圍為4.83～7.85 N，均值為6.38 N，最大破碎波浪力發生在57.589 s；Apelt在實驗中連續測量了100個波的最大波浪力，測量值范圍為4.14～7.99 N，均值為5.64 N。可以認為，本文建立的波浪破碎數值模型能較為準確地計算單樁基礎所受的破碎波浪荷載。比較《海港水文規范》（JTS 145-2-2013）單樁破碎波作用力計算公式，同等波浪條件下計算得最大破碎波浪力為8.96 N，可以認為規范結果略大于實驗和數值模擬結果，作為外包絡值能夠保證荷載計算的設計安全性。

樁前波高為0.084m、水深為0.113m，波長為0.904m，樁徑波長比0.11，滿足Morison公式適用條件（樁徑波長比小于等于0.2），按照Morison公式計算的最大波浪力為5.50 N，小于數值模型和物理實驗給出的破碎波浪力極大值的均值，因此，Morison公式不能準確估算淺水區單樁破碎波浪荷載，需要針對波浪破碎特性做相應的修正。

3 結論

本文基于開源程序OpenFOAM建立了三維波浪破碎數值模型，并將該模型應用于單樁基礎破碎波浪荷載的研究，主要結論如下：

該數值模型通過引入RNGk?ε紊流模型來封閉求解三維不可壓縮粘性流體控制方程，模擬水體的運動；采用VOF方法精確捕捉自由表面。數值模擬結果證明，該模型可以準確模擬波浪破碎前后的非線性變化過程以及自由表面的復雜變化。

與物理模型實驗數據比較，該數值模型計算出的最大破碎波浪力為7.85 N，僅比物理實驗值小1.7%。與《海港水文規范》（JTS 145-2-2013）公式計算結果相比較，該數值模型給出的最大破碎波浪力比公式計算值小12.3%，規范公式作為外包絡值能夠保證荷載計算的設計安全性。以上分析證明，本文建立的三維波浪破碎數值模型可以應用于近岸破碎波浪荷載的模擬計算。

與Morison公式計算結果的比較表明，Morison公式不能直接應用于淺水區單樁破碎波浪荷載的計算。

本文目前主要研究崩破波作用于單樁的破碎波浪荷載，模擬結果具有較好的精度。為了進一步拓展該數值模型的應用范圍，卷破波作用于單樁或群樁的破碎波浪荷載，將是進一步研究的重點工作。

［1］Sawaragi T,Nochinom.Impact forces of nearly breaking waves on a vertical circular cylinder［J］.Coastal Engineering in Japan, 1984,27∶249-263.

［2］Apelt C J,Piorewicz J.Laboratory studies of breaking wave forces acting on vertical cylinders in shallow water［J］.Coastal Engi?neering,1987,11（3）∶263-282.

［3］Swift R H.Prediction of breaking wave forces on vertical cylinders［J］.Coastal Engineering,1989,13（2）∶97-116.

［4］Kyte A,Torum A.Wave forces on vertical cylinders upon shoals［J］.Coastal Engineering,1996,27（3）∶263-286.

［5］馬汝建.淺海結構物波浪力計算及實驗研究［J］.石油大學學報,2000,24（2）∶79-81.mA R J.Experimental study and calculation of wave force acting on the offshore structures in shallow sea water［J］.Journal of the University of Petroleum,2000,24（2）∶79-81.

［6］Miyata H.Finite?difference simulation of breaking waves［J］.Journal of Computational Physics,1986,65（1）∶179-214.

［7］Lin P Z,Liu P L F.A numerical study of breaking waves in the surf zone［J］.Journal of Fluidmechanics,1998,359∶239-264.

［8］Christensen E D.Large eddy simulation of spilling and plunging breakers［J］.Coastal Engineering,2006,53（5）∶463-485.

［9］Lubin P,Vincent S,Abadie S,et al.Three?dimensional large eddy simulation of air entrainment under plunging breaking waves［J］.Coastal Engineering,2006,53（8）∶631-655.

［10］Hieu P D,Tanimoto K.Verification of a VOF?based two?phase flowmodel for wave breaking and wave?structure interactions［J］. Ocean Engineering,2006,33（11-12）∶1 565-1 588.

［11］Wu N T,Oumeraci H,Partenscky H W.Numericalmodeling of breaking wave impacts on a wall［C］//American Society of Civil Engineers.Proceedings of the 24th International Conference on Coastal Engineering.Kobe∶American Society of Civil Engineers, 1994,1∶1 672-1 686.

［12］Zhang S,Yue D K P,Tanizawa K.Simulation of plunging wave impact on a vertical wall［J］.Journal of Fluidmechanics,1996, 327∶221-254.

［13］Mokrani C,Abadie S,Grilli S,et al.Numerical simulation of the impact of a plunging breaker on a vertical structure and subse?quent overtopping event using a Navier?Stokes VOFmodel［C］//International Society of Offshore and Polar Engineers.Proceed?ings of the 20th International Conference on Offshore and Polar Engineering.Beijing∶International Society of Offshore and Polar Engineers,2010∶729-736.

［14］Bredmose H,Jacobsen N G.Breaking wave impacts on offshore wind turbine foundation∶Focused wave group and CFD［C］// American Society ofmechanical Engineers.Proceedings of the 29th International Conference on Ocean,Offshore and Arctic En?gineering.Shanghai∶American Society ofmechanical Engineers,2010,3∶397-404.

［15］Choi J,Yoon S B.Numerical simulations usingmomentum source wave?maker applied to RANS equationmodel［J］.Coastal Engi?neering,2009,56（10）∶1 043-1 060.

［16］張博杰，張慶河.基于OpenFOAM開源程序的無反射數值波浪水槽［J］.中國港灣建設,2012（5）∶1-3. ZHANG B J，ZHANG Q H.2D numerical wave flume without reflection based on OpenFOAM［J］.China Harbour Engineering, 2012（5）∶1-3.

［17］管寧.海上風機基礎結構波流荷載研究［D］.天津∶天津大學,2014.

［18］查晶晶.基于OpenFOAM的數值造波與消波模型及其應用［D］.上海∶上海交通大學,2011.

青島迪拜環球碼頭工程初步設計通過部批

本刊從交通運輸部獲悉，2015年4月，交通運輸部批準了青島港前灣港區迪拜環球碼頭工程初步設計。本工程主要建設2個10萬t級和2個3萬t級集裝箱泊位（水工結構均按靠泊10萬t級集裝箱船舶設計）及配套設施，碼頭泊位長度1 320m，設計年通過能力220萬TEU，工程總概算約41.54億元。（殷缶，梅深）

京唐港主航道將升至25萬t級

本刊從唐山港獲悉，唐山港口實業集團有限公司委托環評機構做的《京唐港區25萬t級航道工程環評》目前已正式向公眾公示。這意味著，京唐港區目前20萬t級航道將提升至25萬t級。25萬t級航道工程以現有第四港池20萬t級航道為基礎，通過浚深、拓寬擴建成25萬t級航道。通航方式為單向通航，乘潮歷時為4 h，乘潮水位為1.10m；設計底標高-22.2m；疏浚土方量1 687萬m3。工程建設期12個月。（殷缶，梅深）

Numerical simulation of breaking wave force on a pile

WANG Xiu?ting1,ZHANG Qing?he1,CHEN Tong?qing1,GUAN Ning2
（1.State Key Laboratory of Hydraulic Engineering Simulation and Safety,Tianjin University,Tianjin 300072, China；2.Tianjin Research Institute for Water Transport Engineering,Tianjin 300456,China）

A 3?D numericalmodel for breaking wave was developed based on the open source code Open?FOAM.The governing equations of themodel are the RANS equations of the three?dimensional incompressible vis?cous fluid.The VOFmethod was used to capture the interface between air and fluid.The wave?maker of distributedmomentum source and the damping wave absorption were applied in themomentum equations for generating and ab?sorbing waves respectively.The RNG k?εturbulencemodel was adopted.The presentedmodel was used to simu?late the breaking wave force acting on the pile on a beach.The numerical results were compared with the previous experimental data.It is shown that themodel is of sufficient accuracy to simulate the breaking wave force onmono?pile in shallow water.

wave?maker of distributedmomentum source;breaking wave;wave force;OpenFOAM;VOF

TV 142；O 242.1

A

1005-8443（2015）03-0191-06

2013-09-22；

2014-11-17

國家高技術研究發展863計劃基金資助項目（2012AA051702）；國家國際科技合作專項資助（2012DFA70490）

王修亭（1989-），男，安徽省宣城市人，碩士研究生，主要從事波浪數值模擬研究。

Biography：WANG Xiu?ting（1989-），male，master student.