波浪對(duì)淹沒(méi)水平圓柱作用的非線性數(shù)值模擬＊

嚴(yán) 斌，林紅星，寧德志＊

（1.大連理工大學(xué) 海岸和近海國(guó)家重點(diǎn)實(shí)驗(yàn)室，遼寧 大連116024；2.中國(guó)科學(xué)院 可再生能源與天然氣水合物重點(diǎn)實(shí)驗(yàn)室，廣東 廣州510640）

波浪與淹沒(méi)潛體相互作用，在潛體上由于水深變淺，波浪非線性增強(qiáng)，會(huì)導(dǎo)致大部分鎖定波轉(zhuǎn)換為自由波的形態(tài)存在，波浪場(chǎng)的變化會(huì)對(duì)潛堤本身及潛堤后工程結(jié)構(gòu)有很大影響。而水平圓柱是眾多淹沒(méi)潛體中應(yīng)用較多的一種，譬如海底石油、天然氣輸送管道和圓柱式防波堤等，開(kāi)展波浪與淹沒(méi)水平圓柱式結(jié)構(gòu)相互作用的研究，對(duì)準(zhǔn)確了解結(jié)構(gòu)后的波浪場(chǎng)及對(duì)結(jié)構(gòu)的作用具有重要意義。

國(guó)內(nèi)外許多學(xué)者已經(jīng)對(duì)波浪與淹沒(méi)潛體作用展開(kāi)了廣泛的研究，如 Lee［1］，Vada［2］，F(xiàn)riis、Grue和Palm［3］應(yīng)用二階勢(shì)流理論計(jì)算二階自由諧波幅值；Cointe［4］采用完全非線性數(shù)值方法模擬了較大幅值波浪對(duì)淹沒(méi)水平圓柱作用下的二階自由諧波；Grue等［5－6］對(duì)深水波浪通過(guò)淹沒(méi)潛體的波浪非線性繞射問(wèn)題進(jìn)行了試驗(yàn)研究，利用淹沒(méi)結(jié)構(gòu)物后兩點(diǎn)處的波面對(duì)產(chǎn)生的高階諧波進(jìn)行分離，得出高階諧波的波幅無(wú)因次量會(huì)隨入射波浪幅增大而增大最終達(dá)到飽和的結(jié)論；鄭永紅等［7］利用一種改進(jìn)的Boussinesq方程，對(duì)淹沒(méi)潛堤上的波浪變形進(jìn)行了數(shù)值模擬；張洪生等［8］建立了一種以新型Boussinesq方程為控制方程的非線性波傳播的數(shù)值模型，能較好地模擬大水深水域和強(qiáng)非線性波的傳播問(wèn)題。陳麗芬等［9］運(yùn)用高階邊界元方法數(shù)值模擬潛堤地形上波浪的傳播變形，得到弱非線性下，潛堤后的基頻波、二階和三階自由波幅值分別與入射波幅成線性、二次和三次函數(shù)關(guān)系。張婷婷［10］采用計(jì)算域內(nèi)造波的方法建立數(shù)值模型，模擬了波浪在潛堤上的傳播演化，并將波面和各階諧波的計(jì)算值與實(shí)驗(yàn)結(jié)果進(jìn)行了比較，得到吻合很好的數(shù)值結(jié)果。Williams［11］對(duì)靠近自由表面的水平淹沒(méi)平板問(wèn)題進(jìn)行數(shù)值研究，得到了直至三階自由諧波的波浪成分。付韻韻等［12］采用多級(jí)子展開(kāi)法對(duì)作用在有限水深中淹沒(méi)水平圓柱上的透反射系數(shù)和波浪力開(kāi)展了理論分析；并在不同淹沒(méi)深度和半徑條件下進(jìn)行了波浪對(duì)水平圓柱作用力的實(shí)驗(yàn)研究。姜?jiǎng)俪龋?3］建立波浪對(duì)淹沒(méi)垂直圓柱繞射問(wèn)題的解析解，考慮了不同淹沒(méi)深度與相對(duì)厚度對(duì)圓柱作用力影響，表明圓柱所受到的最大波浪力并不總是隨著淹沒(méi)深度的增大而減小，而圓柱厚度對(duì)波浪力的影響則較為復(fù)雜，在不同波數(shù)范圍顯示出不同的特性。

本文采用基于域內(nèi)源造波技術(shù)的時(shí)域高階邊界元方法建立波浪與淹沒(méi)水平圓柱相互作用的完全非線性數(shù)值水槽模型，可以實(shí)現(xiàn)在較小計(jì)算域內(nèi)長(zhǎng)歷時(shí)模擬強(qiáng)非線性波浪對(duì)淹沒(méi)圓柱作用過(guò)程。并通過(guò)兩點(diǎn)分離法得到圓柱后的高階自由諧波和鎖定波，研究及其隨入射波浪非線性、圓柱尺寸、淹沒(méi)水深等因素的變化規(guī)律。

1 數(shù)學(xué)模型

1.1 控制方程和初始邊界條件

考慮波浪和一固定水平淹沒(méi)圓柱相互作用（圖1），建立一笛卡爾坐標(biāo)系Oxz，使得z＝0位于靜水面上，且z軸向上為正，x軸向右為正方向，x＝（x，z）代表任一點(diǎn)的坐標(biāo)。圖中，Ω代表整個(gè)計(jì)算域，自由水面用Гf表示，水底為Гd，圓柱邊界為Гb，入射邊界為ГI，靜水面至圓柱最上端的距離為h，靜水深為H。

圖1 水槽中波浪與一淹沒(méi)水平圓柱作用示意圖Fig.1 A sketch map showing the wave action on a horizontal submerged cylinder in a flume

在理想流體假定下，流域內(nèi)可用速度勢(shì)來(lái)描述。速度勢(shì)通常滿(mǎn)足拉普拉斯方程，為了防止從結(jié)構(gòu)物反射回來(lái)的波浪在入射邊界形成二次反射，本文在計(jì)算域內(nèi)引入內(nèi)嵌造波源，故速度勢(shì)滿(mǎn)足泊松（Possion）方程：

在固定邊界Гd和Гb上滿(mǎn)足速度勢(shì)法向?qū)?shù)為0。在自由水面Гf上滿(mǎn)足完全非線性動(dòng)力學(xué)和運(yùn)動(dòng)學(xué)邊界條件，采用半混合歐拉－拉格朗日方法更新自由水面，同時(shí)采用人工層布置在水槽兩側(cè)，來(lái)消除出流波浪，進(jìn)而自由水面邊界條件可寫(xiě)成如下形式［14］：

式中，η表示自由表面；g表示重力加速度；阻尼項(xiàng)μ（x）表示為如下形式

1.2 邊界積分方程

在整個(gè)流域內(nèi)應(yīng)用格林第二函數(shù)，上述邊值問(wèn)題可轉(zhuǎn)化為如下的邊界積分方程［15］：

式中，ω為波浪圓頻率；L為阻尼層寬度；x1（2）分別為左右阻尼層起始位置；α為阻尼系數(shù)，；本文取L為1.5倍入射波長(zhǎng)，系數(shù)α＝1。

由于是進(jìn)行時(shí)域數(shù)值模擬，需要給定初始邊界條件，這里給定靜初始邊界條件，即：

式中，p＝（x0，z0）為源點(diǎn)；q＝（x，z）為場(chǎng)點(diǎn)；α（p）為固角系數(shù)；Γ為流域邊界，包括自由水面邊界和固體邊界。G為簡(jiǎn)單格林函數(shù)，考慮水底的鏡像作用，可以寫(xiě)成如下形式［16］：

對(duì)式（7）采用三節(jié)點(diǎn)的二次邊界元進(jìn)行離散，對(duì)每個(gè)單元通過(guò)數(shù)學(xué)變換，將其轉(zhuǎn)換成參數(shù)坐標(biāo)（ξ）下的等參單元，采用二次形狀函數(shù)插值方法保證單元內(nèi)物理量分布的連續(xù)性。積分方程經(jīng)離散之后，可以建立關(guān)于速度勢(shì)和速度勢(shì)導(dǎo)數(shù)的線性方程組進(jìn)行求解。由于所有節(jié)點(diǎn)和網(wǎng)格在每一時(shí)間步都要更新，所以方程（7）也要在每一時(shí)間步建立和求解一次，并采用四階Runge－Kutta法計(jì)算下一時(shí)刻的速度勢(shì)和波面。

1.3 淹沒(méi)潛體后的波浪組成

當(dāng)波浪傳播至水平淹沒(méi)圓柱上方時(shí)，水深變淺，非線性作用增強(qiáng)，產(chǎn)生與基頻同速度傳播的高頻諧波，通過(guò)圓柱后，由于水深增大，波浪的非線性作用相對(duì)減弱，部分高頻諧波由鎖定波釋放為自由波。在水平淹沒(méi)圓柱后既包括鎖定波也包含高頻的自由波，因此，水平圓柱后x點(diǎn)處的波面可以寫(xiě)成如下的形式：

對(duì)于式（9）中的未知量，可以采用兩點(diǎn)法［15，17］把潛體后兩點(diǎn)波面時(shí)間歷程帶入到式（9）中，并利用三角變換及三角函數(shù)正交性來(lái)求解得到。

2 數(shù)值結(jié)果和分析

作為算例，考慮一靜水深H＝0.45m的水槽中角頻率ω＝6.594rad/s或ω＝7.66rad/s的規(guī)則波與一淹沒(méi)水平圓柱相互作用，水槽長(zhǎng)度為10倍波長(zhǎng)，淹沒(méi)水平圓柱中心橫坐標(biāo)為5.6m，圖1中4個(gè)測(cè)點(diǎn)分別布置在x1＝5.2m，x2＝5.4m、x3＝6.2m和x4＝6.4m。通過(guò)開(kāi)展數(shù)值收斂性試驗(yàn)，x和z兩方向的網(wǎng)格數(shù)分別300和10，圓柱結(jié)構(gòu)網(wǎng)格數(shù)為40個(gè)。

圖2是入射波幅A＝16mm，靜水面距離圓柱最上端距離h＝50mm，圓柱半徑和波浪角頻率分別為（R＝100mm，ω＝7.66rad/s）和（R＝190mm，ω＝6.594rad/s）兩種種工況下，在t＝23T和t＝25T沿水槽x方向的波面分布，圖中豎向虛線為圓柱中心所在位置（x＝5.6m）。從圖中可以看出，在結(jié)構(gòu)前的波面都是規(guī)則的，而在結(jié)構(gòu)之后的波面是很不規(guī)則的，由此判斷出潛堤前后的波浪組成及其相位角是有很大區(qū)別的，這也與堤前主要以鎖定波為主體、堤后以自由波為主體的說(shuō)法相一致。同時(shí)可以看出兩個(gè)時(shí)刻下的波面歷程重合很好，水槽兩端阻尼層處出流波浪得到很好的吸收，未有明顯的波浪反射現(xiàn)象發(fā)生，證明了本文數(shù)值模型的穩(wěn)定性。

圖2 波面在t＝23T和t＝25T沿水槽x方向的分布Fig.2 Wave surface distributions at t＝23Tand t＝25Talong the xdirection of the flume

圖3是不同工況下堤后線性和二階自由諧波幅值隨入射波幅A的變化關(guān)系，圖中還給出了本文數(shù)值模擬結(jié)果與試驗(yàn)結(jié)果［6］的對(duì)比，自由諧波幅值通過(guò)除以入射波幅值無(wú)量綱化。從圖中可以看出，隨著入射波幅A的增大，堤后線性自由諧波幅值衰減，而二階自由諧波幅值先是增大到一極限值，然后隨著入射波幅的進(jìn)一步增大保持飽和狀態(tài)（圖3a）或遞減狀態(tài)（圖3b和c），在圖3b和c工況中，二階自由諧波幅值最大可接近入射波幅的0.4倍，而線性自由諧波幅值最小可接近入射波幅的0.7倍。通過(guò)圖a和b對(duì)比可知，隨著波浪角頻率的增大，波浪非線性增強(qiáng)，堤后自由諧波幅值變化斜率也更大；通過(guò)圖a和c對(duì)比可知，隨著淹沒(méi)圓柱尺寸的增大，堤后自由諧波幅值變化斜率越大，堤后波浪非線性效果更明顯。同時(shí)可以看出，本文數(shù)值模擬結(jié)果與試驗(yàn)結(jié)果吻合的很好，進(jìn)一步驗(yàn)證本文模型的準(zhǔn)確性。

圖3 堤后線性和二階自由諧波幅值隨入射波幅變化關(guān)系Fig.3 Changes of the first and the second－order free harmonics amplitudes at the lee side of the submerged body with the incident wave amplitude

圖4是對(duì)應(yīng)圖3中非線性最強(qiáng)工況（b、c）的淹沒(méi)圓柱后二階鎖定波幅值隨著入射波幅A的變化規(guī)律，圖中一并給出了與解析解［17］的對(duì)比情況。從圖中可以看出，鎖定波的量級(jí)明顯小于同一工況下的自由波，進(jìn)一步說(shuō)明由于潛體的作用，大部分鎖定波被轉(zhuǎn)化為自由波形態(tài)存在于堤后。本文數(shù)值結(jié)果與理論結(jié)果吻合的很好，進(jìn)一步說(shuō)明本文模型的可靠性。

圖4 堤后二階約束諧波幅值隨入射波幅A變化關(guān)系Fig.4 Changes of the second－order bounded harmonics amplitudes at the lee side of the submerged body with the incident wave amplitude A

圖5是靜水深H＝0.45m，入射波幅A＝5.8mm，圓柱半徑R＝100mm，角頻率ω＝2π×1.05rad/s情況下，基頻波幅值、二階自由波幅值、三階自由波幅值和二階鎖定波幅值隨淹沒(méi)水深h的變化關(guān)系。由圖5可知，隨著淹沒(méi)水深的增加，基頻波波幅增大并趨于入射波幅A，二階自由諧波波幅、三階自由諧波波幅呈高次函數(shù)關(guān)系快速減小，并逐漸趨于N（0）量級(jí)，二階鎖波波幅快速減小至N（0）量級(jí)，且在同等條件下量級(jí)遠(yuǎn)遠(yuǎn)小于二階自由波幅值，這與圖4得到的結(jié)論是一致的。圖5進(jìn)一步說(shuō)明了水深對(duì)非線性的影響，及高階自由諧波的產(chǎn)生與波浪非線性的密切關(guān)系。

圖5 自由波和鎖定波幅值隨淹沒(méi)水深的變化關(guān)系Fig.5 Changes of the free and bounded harmonics amplitudes with the submerged water depth

圖6是靜水深H＝0.45m，入射波幅A＝6.0mm，淹沒(méi)水深h＝50mm，ω＝2π×1.05rad/s情況下，基頻波幅值、二階自由波幅值和三階自由波幅值隨淹沒(méi)圓柱半徑R的變化關(guān)系。從圖中可以看出，隨著圓柱半徑R的增大，基頻波幅值a1/A逐漸遞減，而二階自由諧波波幅值a2/A和三階自由諧波幅值a3/A逐漸增大，但最終都趨于飽和，即達(dá)到一常量，而三階自由諧波最先趨于飽和。這與Grue［7］在試驗(yàn)研究中發(fā)現(xiàn)在特定條件下高頻自由諧波隨著入射波幅值增大而趨于飽和的結(jié)論是相似的。

圖6 堤后自由諧波幅值與淹沒(méi)圓柱半徑的關(guān)系Fig.6 Relationship between the amplitude of free harmonics at the lee side of the submerged body and the radius Rof the submerged cylinder

3 結(jié)語(yǔ)

本文采用基于域內(nèi)源造波技術(shù)的時(shí)域高階邊界元方法建立了波浪與淹沒(méi)潛體相互作用的完全非線性數(shù)值水槽模型，并通過(guò)兩點(diǎn)法分離得到淹沒(méi)潛體后的高階自由波和鎖定波。通過(guò)與試驗(yàn)數(shù)據(jù)和解析解的對(duì)比驗(yàn)證了本文模型的準(zhǔn)確性。通過(guò)數(shù)值試驗(yàn)分析可以得到以下結(jié)論：隨著入射波浪非線性增強(qiáng)，基頻波幅值逐漸減小，高階自由諧波幅值逐漸增大至一最高值呈飽和狀態(tài)或減小變化，而鎖定波波幅盡管也隨著入射波幅增大而增大，但在量級(jí)上遠(yuǎn)小于對(duì)應(yīng)的自由波；潛體淹沒(méi)水深越小，高階自由諧波貢獻(xiàn)越大，基頻波衰減越大。

（References）：

［1］LEE C M.The second order theory of heaving cylinders in a free surface［J］.J.Ship Res.，1968，12：313－327.

［2］VADA T.A numerical solution of the second－order wave diffraction problem for a submerged cylinder arbitrary shape［J］.J.Fluid Mech.，1987，174：23－37.

［3］FRITS A，GRUE J，PALM E.Mathematical Approaches in Hydro－dynamics［M］.Philadelphia：SIAM，1991.

［4］COINTE R.Nonlinear simulation of transient free surface flows［C］//MORI K H.5th Intl.Conf.on Numerical Ship Hydrodynamics.Washington D C：the national academics press，1989：239－250.

［5］GRUE J，GRANLUND K.Impact of nonlinearity upon waves traveling over a submerged cylinder［C］//KORSMEYER F T.3rd Intl Workshop on Water Waves and Floating Bodies，Cambridge：MIT，1988：57－60.

［6］GRUE J.Nonlinear water waves at a submerged obstacle or bottom topography［J］.J.Fluid Mech.，1992，244：455－476.

［7］ZHENG Y H，SHEN Y M，WU X G，et al.Simulation of wave propagation over a submerged bar using a modified Boussinesq equation［J］.Engineering Science，2004，6（4）：34－40.鄭永紅，沈永明，吳修廣，等.用改進(jìn)的Boussinesq方程模擬潛堤上的波浪變形［J］.中國(guó)工程科學(xué)，2004，6（4）：34－40.

［8］ZHANG H S，F(xiàn)ENG W J，WANG Y L，et al.A new approach to numerical simulation of nonlinear wave propagation and its experimental verification［J］.Acta Oceanologica Sinica，2007，29（4）：137－147.張洪生，馮文靜，王亞玲，等.非線性波傳播的新型數(shù)值模擬模型及其實(shí)驗(yàn)驗(yàn)證［J］.海洋學(xué)報(bào)，2007，29（4）：137－147.

［9］CHEN L F，NING D Z，TENG B，et al.Study on the higher free harmonic waves on the lee side of a submerged bar［J］.Acta Oceanologica Sinica，2011，33（6）：1－12.陳麗芬，寧德志，滕斌，等.潛堤后高階自由諧波的研究［J］.海洋學(xué)報(bào)，2011，33（6）：1－12.

［10］ZHANG T T.Study on the evolution of the wave propagation over the Submerged Sill［D］.Dalian：Dalian University of Technology，2011.張婷婷.波浪在潛堤上的傳播演化研究［D］.大連：大連理工大學(xué)，2011.

［11］WILLIAMS J A.A nonlinear problem in surface water waves［D］.Berkeley：University of California，1964.

［12］FU Y Y.Analytical solution for wave action on submerged horizontal circular cylinder［D］.Dalian：Dalian University of Technology，2012.付韻韻.波浪對(duì)水平圓柱作用力的理論分析［D］.大連：大連理工大學(xué)，2012.

［13］JIANG S C，TENG B，NING D Z，et al.An analytical solution of wave diffraction by a submerged vertical cylinder［J］.The Ocean Engineering，2010，28（3）：69－75.姜?jiǎng)俪?，滕斌，寧德志，?波浪對(duì)淹沒(méi)垂直圓柱繞射解析解［J］.海洋工程，2010，28（3）：69－75.

［14］NING D Z.Application of the fast multipole boundary element method in fully nonlinear water wave problems［D］.Dalian：Dalian University of Technology，2005.寧德志.快速多極子邊界元方法在完全非線性水波問(wèn)題中的應(yīng)用［D］.大連：大連理工大學(xué)，2005.

［15］BRORSEN M，LARSEN J.Source generation of nonlinear gravity waves with the boundary integral equation method［J］.Coastal Engineering，1987，11（2）：93－113.

［16］NEWMAN J N.Approximation of free－surface Green functions［M］.Cambridge：Cambridge University Press，1992：107－142.

［17］ZOU Z L.Water wave theories and their applications［M］.Beijing：Science Press，2005.鄒志利.水波理論及其應(yīng)用［M］.北京：科學(xué)出版社，2005.

［18］RYU S，KIM M H，LYNETT P J.Fully nonlinear wave－current interactions and kinematics by a BEM based numerical wave tank［J］.Computational Mechanics，2003，32（4－6）：336－346.

［19］LI Y C，TENG B.Wave Action on Marine Structures：the second edition［M］.Beijing：Ocean Press，2002.李玉成，滕斌.波浪對(duì)海上建筑物的作用：第2版［M］.北京：海洋出版社，2002.

［20］BREBBIA C A，WALKER S.Boundary Element techniques in Engineering［M］.London：Newnes－Butterworths，1980.