單樁波浪破碎荷載的數(shù)值模擬研究

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

（1.天津大學(xué)水利工程仿真與安全國家重點(diǎn)實(shí)驗(yàn)室，天津300072；2.交通運(yùn)輸部天津水運(yùn)工程科學(xué)研究所，天津300456）

單樁波浪破碎荷載的數(shù)值模擬研究

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

（1.天津大學(xué)水利工程仿真與安全國家重點(diǎn)實(shí)驗(yàn)室，天津300072；2.交通運(yùn)輸部天津水運(yùn)工程科學(xué)研究所，天津300456）

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

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

海洋石油開采平臺、海上風(fēng)力發(fā)電機(jī)、高樁碼頭等近海結(jié)構(gòu)物經(jīng)常采用樁式基礎(chǔ)。樁式基礎(chǔ)所受到的波浪荷載，特別是近岸波浪破碎荷載對樁式基礎(chǔ)的作用是工程設(shè)計(jì)所關(guān)心的重要問題。

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

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

1 波浪破碎數(shù)值模型的建立

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

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

造波源分布在一定寬度的區(qū)域內(nèi)，水槽兩端設(shè)置消波區(qū)域，如圖1所示。

1.1 控制方程

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

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

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

1.2 邊界條件和初始條件

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

計(jì)算初始條件設(shè)定全場速度和動壓為零。

1.3 計(jì)算時空步長

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

1.4 波浪破碎

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

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

2 波浪破碎荷載數(shù)值模擬的驗(yàn)證

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

2.1 數(shù)值模型驗(yàn)證水槽設(shè)置

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

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

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

表1 驗(yàn)證算例波高測點(diǎn)平面位置Tab.1 Horizontal coordinate of probes in verification simulationm

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

2.2 計(jì)算網(wǎng)格劃分

單樁底面中心坐標(biāo)為（19.1，0.24，1.0），直徑0.102m。網(wǎng)格均采用非結(jié)構(gòu)化的六面體單元。為準(zhǔn)確模擬單樁對波浪場的影響，在單樁的壁面附近加密網(wǎng)格；自由表面附近加密網(wǎng)格，以提高自由表面的模擬精度。水槽兩端的阻尼消波區(qū)內(nèi)，由計(jì)算域向水槽兩端網(wǎng)格稀疏漸變，以增加數(shù)值耗散，提高消波效果。

2.3 模擬結(jié)果分析

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

圖4 破波點(diǎn)附近波高測點(diǎn)（#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 單樁破碎波浪力包絡(luò)曲線Fig.7 Envelope of breaking wave forces acting on the pile

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

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

3 結(jié)論

本文基于開源程序OpenFOAM建立了三維波浪破碎數(shù)值模型，并將該模型應(yīng)用于單樁基礎(chǔ)破碎波浪荷載的研究，主要結(jié)論如下：

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

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

與Morison公式計(jì)算結(jié)果的比較表明，Morison公式不能直接應(yīng)用于淺水區(qū)單樁破碎波浪荷載的計(jì)算。

本文目前主要研究崩破波作用于單樁的破碎波浪荷載，模擬結(jié)果具有較好的精度。為了進(jìn)一步拓展該數(shù)值模型的應(yīng)用范圍，卷破波作用于單樁或群樁的破碎波浪荷載，將是進(jìn)一步研究的重點(diǎn)工作。

［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］馬汝建.淺海結(jié)構(gòu)物波浪力計(jì)算及實(shí)驗(yàn)研究［J］.石油大學(xué)學(xué)報,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開源程序的無反射數(shù)值波浪水槽［J］.中國港灣建設(shè),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］管寧.海上風(fēng)機(jī)基礎(chǔ)結(jié)構(gòu)波流荷載研究［D］.天津∶天津大學(xué),2014.

［18］查晶晶.基于OpenFOAM的數(shù)值造波與消波模型及其應(yīng)用［D］.上海∶上海交通大學(xué),2011.

青島迪拜環(huán)球碼頭工程初步設(shè)計(jì)通過部批

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

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

本刊從唐山港獲悉，唐山港口實(shí)業(yè)集團(tuán)有限公司委托環(huán)評機(jī)構(gòu)做的《京唐港區(qū)25萬t級航道工程環(huán)評》目前已正式向公眾公示。這意味著，京唐港區(qū)目前20萬t級航道將提升至25萬t級。25萬t級航道工程以現(xiàn)有第四港池20萬t級航道為基礎(chǔ)，通過浚深、拓寬擴(kuò)建成25萬t級航道。通航方式為單向通航，乘潮歷時為4 h，乘潮水位為1.10m；設(shè)計(jì)底標(biāo)高-22.2m；疏浚土方量1 687萬m3。工程建設(shè)期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

國家高技術(shù)研究發(fā)展863計(jì)劃基金資助項(xiàng)目（2012AA051702）；國家國際科技合作專項(xiàng)資助（2012DFA70490）

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

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