直墻前多向不規(guī)則波波高平面分布的數(shù)值模擬研究

羅天翔，朱良生*，徐潤剛，程趙德

（1.華南理工大學土木與交通學院，廣東 廣州 510640；2.中交第四航務工程勘察設計院有限公司，廣東 廣州 510230）

0 引言

由于離岸建筑物的布置在港灣工程的重要性逐漸提高，如今工程設計不僅關(guān)注波要素對直墻的作用，對于直墻前沿水域的波高平面分布也要有全面的了解，以便能夠合理地確定直墻前沿的泊位、浮碼頭等工程布置位置。實際波浪往往是不規(guī)則且方向不單一的，多向不規(guī)則波更符合天然波浪的特性。因此研究正向、斜向入射后的多向不規(guī)則波經(jīng)歷反射、疊加作用后波浪在直墻建筑物前沿水域的平面空間分布具有重要意義。

陳漢寶[1]使用數(shù)值模擬和物理模型探討了斜向波與直立式防波堤的相互作用，波高最大值與斜角的關(guān)系，給出了對應的波態(tài)圖，定性分析了波浪特性；馮衛(wèi)兵等[2]利用淺水定常莆田經(jīng)驗譜分解結(jié)果等模擬非線性隨機波在直墻前的運動和動力過程，形成了波高、周期等一系列極值序列、分布曲線，但二階非線性理論解不能完全表現(xiàn)近岸波浪的非線性，對波高等要素的統(tǒng)計分布停留在直墻上的概率分布；俞聿修等[3]探討了波浪斜向性和多向性對不規(guī)則波作用于直墻堤上波浪荷載的影響，并提出了實用的計算方法；駱俊彬等[4]采用物理模型試驗證明不規(guī)則波作用下的直墻前的波谷壓力概率分布可以用韋伯分布來描述，并與相對波高有關(guān)；姚國權(quán)等[5]在從堤頭多角度入射的波況下研究了規(guī)則波與不規(guī)則波反射問題，推導了規(guī)則波沿直立堤外側(cè)波高分布的計算式，沿直立堤方向橫截面的不規(guī)則波和規(guī)則波的波高分布，這里未強調(diào)波浪的非線性；周枝榮等[6]結(jié)合物理試驗和相關(guān)研究成果，提出斜向波作用在矩形重力墩結(jié)構(gòu)上最大總水平波浪力簡化計算公式；張文海等[7]分析了統(tǒng)一方程，并以此為基礎建立了從深水到淺水域都有效的數(shù)值模型。劉思等[8]首先根據(jù)給定的波要素、群性參數(shù)和方向分布參數(shù)，給出了多向不規(guī)則波群的數(shù)值模擬方法，在此基礎上利用有限元法求解改進的Boussinesq 方程的數(shù)值計算模型，建立了數(shù)值水池在指定位置模擬所要求群性的多向不規(guī)則波入射波浪邊界條件的計算方法。欒英妮[9]采用Mike21 BW 數(shù)學模型，并選用文獻[10]中的計算模型進行計算并與其提供的試驗結(jié)果比較，分析結(jié)果表明，Mike21 BW 數(shù)學模型與前人試驗結(jié)果相近，適用于非線性較強、淺水以及有限水深情況下的計算；張娜等[11]使用Mike21 BW 數(shù)學模型對波浪在緩坡上的破碎現(xiàn)象及爬坡過程進行了精確模擬，并將此方法成功應用于實際工程中。

關(guān)于波浪傳播入射后直墻前波高分布的研究，國內(nèi)外學者已對其進行了大量的理論推算、物理模型以及數(shù)值模擬，并得到有效的計算方法和經(jīng)驗公式。但前人的研究多考慮于直墻上波浪力、反射系數(shù)，對波高分布的研究討論集中在直墻與波浪交匯處以及沿堤方向的橫切面，相比前人的研究，本文完全考慮不規(guī)則波的多向性、非線性，且對直墻前整體水域平面的波高分布以及波長和腹點、節(jié)點位置、間距的關(guān)系進行了較系統(tǒng)研究。本研究采用基于Boussinesq 方程的Mike21 BW 模型，對正向、斜向入射的多向不規(guī)則波在直墻前的反射、疊加過程進行數(shù)值模擬，探討入射波與直墻前有效波高的平面空間分布關(guān)系，不同波況的有效波高分布以及得到腹點、節(jié)點間距和位置與入射波波長的關(guān)系。

1 試驗研究方法

1.1 Mike21 BW 數(shù)學模型簡介

計算采用丹麥DHI 的Mike21 中的BW 波浪模型，該模型通過求解沿垂向積分的Boussinesq方程獲得沿水深的平均流速、水位變化以及波高等物理量。模型中考慮了波浪的折射、繞射、反射、底部損耗等因素，可以很好地滿足計算內(nèi)容的要求。對方程采用ADI（交替方向法）算法，式中各項保證二階精度。新形式的Boussinesq 方程改進了色散關(guān)系式，使BW 模型適合模擬從深水到淺水的傳播。

1.2 基本方程

式中：x、y 為水平坐標，m；t 為時間，s；ξ 為高出平均水位的水面高度，m；p、q 分別為x、y 方向的流量密度；h 為水深，m；D 為平均水深，m；c為謝才阻力系數(shù)，c=M·h1/6，m0.5/s；M 為曼寧系數(shù)，m1/3/s；E 為紊動“渦黏”系數(shù)；g 為重力加速度，g=9.81 m/s2。

1.3 方向譜

模型中，本文關(guān)于方向譜的設置包括對主波向的定義、方向譜的形式以及具體傳播參數(shù)的設置，波浪頻譜選擇Jonswap 譜：

式中：f 為波頻；f0為譜峰頻率，s-1；γ 為譜峰尖度因子；H1/3為有效波高，m；T0為平均周期，s；β為峰型參數(shù)；α 為γ 的函數(shù)。方向分布類型選擇頻率相關(guān)分布，方向分布函數(shù)G（θ）具體表現(xiàn)形式為：

式中：s 為頻率相關(guān)形狀參數(shù)；θmain為主波向，（°）。具體傳播參數(shù)設置為主波向θmain，最大偏離角度Δθmax兩項，為了簡化試驗，本文中Δθmax固定為常用值30°。

1.4 計算實例

1.4.1 試驗概況及計算條件

計算水域尺度為1 000 m×7 000 m，水深8 m。造波線距離左側(cè)邊界200 m，左、上、下側(cè)邊界分別設置海綿層，右側(cè)邊界設置為直立實體結(jié)構(gòu)，具體布置情況如圖1 所示。本文數(shù)值計算多向不規(guī)則波經(jīng)直立實體結(jié)構(gòu)全反射后全水域波浪的水平分布，各計算工況條件見表1，其中主波向角度θmain指波浪的主波向與造波線法向的夾角。

圖1 試驗地形示意圖Fig.1 Schematic diagram of experimental topography

表1 波浪初始入射條件Table 1 Initial wave incident conditions

1.4.2 模型設置

采用矩形網(wǎng)格劃分，空間步長Δx=Δy=1 m，時間步長Δt=0.1 s。近岸采用直立實體結(jié)構(gòu)，反射系數(shù)為1，模型四周在開邊界上設置海綿層，吸收超出模型范圍的波浪，造波線后面產(chǎn)生的波浪也被置于該處的海綿層吸收。

試驗造波方式選擇多向不規(guī)則波，Jonswap 譜參數(shù)的譜峰尖度因子γ=3.3，σa=0.07，σb=0.09。取穩(wěn)定后200 個左右波的實時波面數(shù)據(jù)進行統(tǒng)計分析，測點間距1 m，采取的測點與直墻的最遠距離為入射波的2 倍深水波長，造波線設置在距離左、上以及下邊界200 m 處，依據(jù)模型的計算要求，最小周期（截斷周期）設置為0.7 s。

消波系數(shù)按式（6）計算：

式中：α 和r 為給定的常數(shù)；Nsponge為消波線的數(shù)量。海綿層的寬度一般為波長的1～2 倍。

2 數(shù)學模型驗證

選取Demirbilek 和Nwogu 的物模試驗工況對數(shù)學模型進行對比、驗證。模型布置見圖2。試驗從造波端至右側(cè)斜坡邊界一共布置9 個波高儀，造波端左側(cè)布置的海綿層厚度為5 m，右側(cè)斜坡邊界的坡度為1∶12。波浪在深水區(qū)經(jīng)歷一段等水深傳播至由3 段斜坡組成的復合斜坡上，深水區(qū)的水深為0.531 m，復合斜坡平均坡度為1∶10.2。經(jīng)過斜坡后波浪在等水深的淺水段傳播，淺水區(qū)的水深為0.031 m。模擬波浪為不規(guī)則波，有效波高HS=0.075 m，譜峰周期TP=1.5 s。

圖2 Demirbilek 和Nwogu 的試驗布置示意圖Fig.2 Schematic diagram of Demirbilek and Nwogu experimental layout

數(shù)值模型參數(shù)：空間步長Δx=0.05 m，時間步長Δt=0.002 s，底摩擦系數(shù)fw=0.011，換算成謝才系數(shù)（Chezy Number）為55。在出口處布置1~2倍波長的海綿層。模型模擬時長為900 s，模擬100 s后波浪特征趨于穩(wěn)定，穩(wěn)定后的波浪用作統(tǒng)計分析。

沿程波高的數(shù)模值和物模值對比見圖3，試驗結(jié)果表明波高的數(shù)模試驗值與物模試驗值比較一致，說明模型可用于不規(guī)則波在復雜地形上的傳播模擬。

圖3 Demirbilek 和Nwogu 物模試驗值與MIKE21 BW 模型有效波高試驗值對比Fig.3 Comparison of the experimental values of the Demirbilek and Nwogu physical model with the experimental values of the effective wave height of the MIKE21 BW model

3 計算結(jié)果分析

3.1 離岸方向HS 空間分布

3.1.1 正向入射波況

在不規(guī)則波入射直墻后反射的情況中，也將最大波高Hmax值處視為腹點，最小波高Hmin值處視為節(jié)點。本文采用兩種波浪周期求得波長值進行分析，分別為譜峰周期TP對應波長LP和平均周期對應波長L。

圖4 為對3 種主波向為0°的工況結(jié)果進行統(tǒng)計得到的各測點的HS。圖中x 軸為距直墻相對距離s/LP或s/L，s 為距直墻實際距離，m；LP為入射波當前水深譜峰周期TP對應的波長值，m；L為入射波當前水深平均周期對應的波長值；y 軸為該點有效波高HS與入射波有效波高HS0的比波高HS/HS0。

圖4 主波向為0°時直墻前的波高空間分布曲線Fig.4 Spatial distribution curve of wave height in front of a straight wall when the main wave direction is 0°

相對距離分別為s/L、s/LP時，不同工況下的波高空間分布曲線一致，其中相對距離為s/LP時，其曲線的腹點和節(jié)點的分布位置與規(guī)則波入射直墻反射后形成的駐波相同。在遠離直墻的方向上，腹、節(jié)點處的波高值變化很大，振幅隨著離堤距離的增加變小，2 個平均波長距離后腹、節(jié)點處的波高值變化與平均波高之比小于5%。根據(jù)波高與波浪能量存在H2-E 的關(guān)系，波浪入射直墻一定時間后，直墻前波能為2E，HS為入射波的倍，即HS/HS0≈1.4，而試驗工況離堤距離達到一定值，波高分布曲線穩(wěn)定后腹點處的比波高值在1.4 左右，這也說明了試驗的合理性。這表明不同工況的直墻前波高空間分布取決于入射波當前水深下的波長，其中譜峰周期TP對應波長LP與之關(guān)聯(lián)性最強。

3.1.2 斜向入射波況

限于篇幅，在此給出平均周期為12 s 時的不同主波向工況條件下的有效波高分布圖及其腹點包絡線，各入射波向波況下波高空間分布曲線如圖5 如示。

圖5 各入射波向波況下波高空間分布曲線Fig.5 Wave height spatial distribution curve under various incident wave directions and wave conditions

在靠近直墻的部分，主波向為0°的工況波高整體高于斜向入射的工況，且入射角度越大，比波高整體越低。當將波浪分解成橫向、縱向傳遞，入射角度越小，正向傳遞的波浪能量越大，直墻前波浪的反射、疊加作用越大。但在遠離直墻一定距離之后，由于斜向入射條件下的波高空間分布曲線整體相位被拉長，腹、節(jié)點的分布并不均勻，這是因為在遠離直墻的位置，由于海浪的多向性，疊加作用強于主波向正向的工況，主波向角度大的曲線不一定低于主波向角度小的曲線，靠近直墻的腹點不一定比遠離直墻的腹點高。直墻上的比波高值普遍小于2。由圖5 可知，不論何種波況下，直墻前平面波高最低點為離墻最近的節(jié)點處，且該點波高比之后的節(jié)點處波高值小很多。

3.2 反射后腹點、節(jié)點間距Lf 與波長L 的關(guān)系

分別對腹點和節(jié)點的間距統(tǒng)計分析，Lf與L的關(guān)系尚不明確，下面通過彌散方程以及試驗數(shù)據(jù)探討兩者之間的聯(lián)系。首先確定Lf如何統(tǒng)計分析，Lf的測量方式有兩種，一種是腹點的間距Lf1，第二種是節(jié)點的間距Lf2。將各種工況的波高空間分布曲線中，根據(jù)上述對不規(guī)則波入射直墻后的腹點、節(jié)點的判斷方式對腹點、節(jié)點的位置和波高值進行統(tǒng)計，相鄰的腹點、節(jié)點的間距即為Lf的計算方法。將各工況的Lf1和Lf2統(tǒng)計，并引入公式計算其離散度σ。其表達式如下：

由表2、表3 可知，腹點和節(jié)點在主波向為0°時分布均勻，它們的間距可以整體被考慮，在主波向角度變大時，分布不再均勻且整體并無明顯規(guī)律。

表2 HS 統(tǒng)計中各波況Lf1 的離散度Table 2 The dispersion degree of each wave condition Lf1 in HS statistics

表3 HS 統(tǒng)計中各波況Lf2 的離散度Table 3 The dispersion degree of each wave condition Lf2 in HS statistics

根據(jù)3.1 節(jié)中對主波向正向入射的波況，波高空間分布曲線與當前水深、譜峰周期對應波長LP，再加上入射角度、疊加作用等因素影響，直墻前Lf、T 和L 的關(guān)系有待進一步確定。受限于篇幅，下面只列出了與Lf值對應關(guān)系最強的LP的計算結(jié)果和結(jié)果分析。

3.2.1 正向入射波況

表4 直墻前第1 個節(jié)點和其相鄰腹點中點比值Table 4 The ratio of the midpoint between the first node in front of the straight wall and its adjacent antinode

表5 正向入射時各波況下的Lf /（1 2Lp）Table 5 LLp）under each wave condition at normal incidence f /（1 2主波向θmain/（°）平均周期/s 8 1.02 1.01 1.01 0 10 12

3.2.2 斜向入射波況

隨著主波向角度的變大，Lf值的離散度增大，波高分布曲線變得越來越不規(guī)則，很難得到準確的與L 的關(guān)系式，難以計算出腹點或節(jié)點的準確位置。

由表4 可知，對于不同入射方向，第1 個節(jié)點與相鄰腹點的中點吻合。為了討論方便，將直墻前第1 段腹點間距定義為L1f。即在L1f中，腹點和節(jié)點的對應關(guān)系不變。因此可以通過入射波況確定L1f對應腹點、節(jié)點的位置，可見對直墻前第1 段腹點間距和波長的討論是有意義的。下面列出L1f與Lp以及不同入射角θmain的對應關(guān)系。

在其他入射波要素相同的情況下，隨著主波向角度θmain的變大，L1f變大。

表6 各波況下的L1f/（1 2cos θmain Lp）Table 6 L1f/（1 2cos θmain Lp）under various wave conditions主波向θmain/（°）平均周期/s 8 10 12 1.02 1.01 1.01 11.25 0.97 1.01 1.00 22.5 0.98 1.00 1.03 33.75 0.99 0.98 1.00 45 0.98 1.00 1.02 0

3.3 直墻前HS 橫向切面的空間分布情況

在波浪正向入射的試驗中，直墻前橫切面上的等波高線大多與直墻平行，由于入射波是多向不規(guī)則波，等波高線不嚴格平行十分正常；在波浪斜向入射的試驗中，直墻前的等波高線與橫切面并不平行，波高分布呈現(xiàn)相位差。

4 結(jié)語

1）多向不規(guī)則波入射后，離岸方向有效波高空間分布情況如下：主波向正向時，不同周期的條件下，堤前相同相對距離s/LP內(nèi)有效波高空間分布曲線一致，腹、節(jié)點分布規(guī)律與規(guī)則波入射直墻的波況一致，在遠離直墻方向腹點對應的比波高值逐漸變小，節(jié)點對應比波高值逐漸變大；斜向入射波況中，隨著主波向偏離法向入射的角度變大，波高分布曲線不規(guī)則程度變大，腹、節(jié)點間距變大且不均勻；各工況在距離直墻一定距離后波高值趨于穩(wěn)定，比波高值在1.4 左右；整體波高分布最小值處為堤前第1 個節(jié)點，且該點波高遠小于其他節(jié)點處波高。

2）多向不規(guī)則波入射直墻后，腹、節(jié)點間距與譜峰周期計算得到波長LP有關(guān)。正向入射時，腹點間距Lf1相較與節(jié)點間距Lf2離散度較小，腹點間距Lf1與LP關(guān)系式為斜向入射時，堤前第1 段L1f與 LP以及主波向角度θmain對應，關(guān)系式為以上關(guān)系式可以依據(jù)入射波確定直墻前兩個腹點的位置和間距以及第1 個節(jié)點的位置，在工程實踐中可以考慮將泊位或浮碼頭的位置設置在第1 個節(jié)點位置及其位置附近。

3）波浪正向入射的試驗中，直墻前橫切面的波高分布基本為等波高線；波浪斜向入射的試驗中，直墻前的等波高線基本不平行于橫切面，波高分布出現(xiàn)相位差。