宋向榮， 周校軍， 高俊亮,2， 王崗

(1.江蘇科技大學 船舶與海洋工程學院，江蘇 鎮(zhèn)江 212003; 2.江蘇科技大學 江蘇省船舶先進設計制造技術重點實驗室，江蘇 鎮(zhèn)江 212003; 3.河海大學 港口海岸與近海工程學院，江蘇 南京 210098)

港灣共振是指當傳入港內的低頻波浪的頻率與港灣的本征頻率相同時，港內水體會出現(xiàn)大幅波動的現(xiàn)象[1]。港灣共振可以誘發(fā)港內泊船的劇烈運動，顯著影響系泊安全和貨物的裝卸[2]。甚至可能造成人員傷亡和重大的經(jīng)濟損失。夏季由于受太平洋上的熱帶風暴影響，我國東南沿海城市的很多港口也曾受到港灣共振的影響[3]。誘發(fā)港灣共振的原因有很多，其中主要能夠誘發(fā)港灣共振的有氣象波、次重力波、地震波、海嘯、海洋內波和剪切流等[4]。

為了對海嘯波誘發(fā)的港灣共振進行預報和研究，Sobey[5]提出了正交模態(tài)分解法，港灣的本征頻率、共振模態(tài)形狀及各模態(tài)的響應幅值可以由這一方法求得。但該方法在全反射邊界條件的處理上不夠精確。Gao等[6]提出鏡像法對其進行了改進，顯著提高了預測結果的精度。隨后Gao等[7]使用改進的正交模態(tài)分解法并結合Boussinesq數(shù)值模型系統(tǒng)研究了在不同入射波幅和港底剖面形狀條件下孤立波對港內波幅演化、共振波幅及港內總波能和相對波能分布的影響。Gao等[8-9]進一步使用改進的正交模態(tài)分解法和Boussinesq數(shù)值模型研究了連續(xù)孤立波和N波的波形參數(shù)對港內最大爬高及港內相對波能分布的影響，上述研究中港內均設置為平底地形。隨后，鄭子波等[10]引入折線型地形，系統(tǒng)研究了該地形條件對孤立波誘發(fā)的港內最大爬高及波能的影響。Gao等[11]利用Boussinesq數(shù)值模型研究了高模態(tài)下雙色短波群在近海岸礁地形附近引起的狹長港內低頻振蕩。

雖然絕大多數(shù)學者采用孤立波來研究海嘯波[12-13]，但人們從許多實際觀察中發(fā)現(xiàn)大多數(shù)傳至海岸的海嘯波波形像“N”，可以看作由一個大波峰和一個大波谷組成，這類海嘯波被稱為N波[14-15]。然而，目前對由N波誘發(fā)的瞬變港灣共振現(xiàn)象的研究還極為有限。因此，本文進一步系統(tǒng)研究不同港內折線型地形對N波入射誘發(fā)港灣共振條件下的港內最大爬高及波能分布的影響。另外，不同N波入射波幅對于港內波況的影響也被研究。

1 數(shù)值模型和分析方法

1.1 FUNWAVE-TVD數(shù)值模型

FUNWAVE-TVD 是Shi等[16]開發(fā)的完全非線性Boussinesq波浪模型，該模型使用混合有限元有限差分格式求解了Chen等[17]提出的完全非線性Boussinesq方程，并且結合了Kennedy等[18]提出的動參考面理論。該模型在頻散中保留到O(μ2)，在非線性a/h保留了所有階(這里μ=h/l，h為局部靜水深，l為波長，a為波幅)。其控制方程表達式為：

ηt+·M=0

(1)

V2-V3-R}+gηh

(2)

式中：g、η和H=h+η分別表示重力加速度、局部表面位移和總局部水深；=((?/?x),(?/?y))為水平梯度運算符；M為水平體積通量：

(3)

(4)

式中：A=·(huα)，B=·uα。V1和V2為色散項。

(6)

V3為豎直渦量場的O(μ2)分布：

(7)

其中：

ω0=(×uα)·iz=vα,x-uα,y

(8)

zα,x(Ay+zαBy)-zα,y(Ax+zαBx)

(9)

式中：iz為豎直方向的單位向量；R為包含底部摩擦和次網(wǎng)格橫向紊流混合的耗散和擴散項。

FUNWAVE-TVD通過一個高階激波捕捉總變差下降(total variation diminishing，TVD)格式求解，這使得該模型不依賴經(jīng)驗公式就可以模擬波浪破碎并結合波能耗散[15]。此外，該模型采用保持強穩(wěn)定性的三階Runge-Kutta格式來自動調整時間步長，通過帶有非阻塞通信的消息傳遞接口(MPI)實現(xiàn)并行運算。通過這些改進，該模型在預測近岸波浪過程時變得更加穩(wěn)定，這些波浪過程包括波浪淺水、折射、衍射、破碎及在平臺和自然沙灘上的爬高和下沖[15]。由于Gao等[9]對FUNWAVE-TVD模型在模擬N波誘發(fā)的瞬變港灣振蕩的適用性進行了驗證，故本文沒有再對該模型進行驗證。

1.2 分析方法

本文各實驗均使用正交模態(tài)分解法來對港內波況進行分析處理。該方法由Sobey[5]提出，用來預測港灣的本征頻率和與之對應的本征模態(tài)以及確定海嘯誘發(fā)的港灣內各共振模態(tài)的響應幅值。該方法包括2個計算步驟。第1步為計算港口的本征頻率和模態(tài)形狀。第2步通過一個多維優(yōu)化問題，對由海嘯波所誘發(fā)的港灣共振過程中的各模態(tài)響應幅值進行預測。本征頻率和模態(tài)形狀將作為已知量參與到第2步的預測中。Gao等[6]隨后對該方法進行了改進以提高本征頻率和模態(tài)形狀的計算精度。由于正交模態(tài)分解法是基于線性理論得到的，只有在港內波浪非線性參數(shù)較小的情況下才能精確地預測港內各模態(tài)的響應幅值[13]。故本文只針對入射N波波幅較小的實驗進行了分離，然后再研究各實驗的港內波能分布。

2 數(shù)值實驗

圖1呈現(xiàn)了所有數(shù)值實驗所使用的狹長型矩形港口俯視圖、港內各地形的剖面圖以及坐標系的布置。所有港口的長度為1 500 m，寬度2b=30 m；坐標系的x、y軸布置在靜水面上，x軸在沿港口長度方向的中心線上，y軸在港口口門處，z軸豎直向上。本文所研究的港內地形為折線型，分別為凸型底、凹型底、平臺型底和斜直線型底。對于凸型底與凹型底地形而言，有且只有一個折點，且折點位于沿港口長度的垂直平分線上，各組中折點處的水深有所區(qū)別。對于平臺型地形來說，可以將其看作是具有2個折點的折線型地形，2個折點處的水深一致并且關于港口長度垂直平分線對稱。圖中h1、h0、hp和l分別代表港口口門處水深、港口底墻處水深、平臺處水深和平臺的長度。在本文中，港口口門處的水深和外海水深均為12 m，底墻處的水深為3 m，平臺處的水深為7.5 m；各組實驗的底部剖面形狀各不相同。港內的平均水深為：

(10)

圖1 不同底部剖面形狀的狹長型矩形港口概圖Fig.1 Elongated harbors with different bottom profile shapes

圖2 數(shù)值波浪水槽示意Fig.2 Schematic diagram of numerical wave flume

3 結果和分析

3.1 波面時間序列

圖3顯示了A組、C組、E組實驗中入射N波波幅A0=0.04 m時測點G01、G76和G151處的波面時間序列，為了更直觀地觀察波面時間序列所以通過A0對其進行了無因次化處理。Ani和Anr分別表示的是入射和反射波波幅，其中n代表的是測點的位置。由于N波入射波幅均相同，港口口門外各水體條件也均一致，所以當波浪傳播到測點G01處時各組的波幅基本均于t1=44.28 s達到峰值。隨后入射N波逐漸向港內傳播，因為港內地形不同，各組的波面時間序列也會有不同的變化。A～E組的平均水深逐漸減小，港內波浪的平均傳播速度也逐漸減小。故從t3開始，A～E組的波幅達到峰值的時間逐漸增大，且下沖的極值受平均水深的影響比爬高要更大一些。由于底墻不可穿透，因此t5和t6時刻對應的是港內最大的爬高和下沖。例如，由于淺水效應，t5從A組的209.48 s逐漸增加到E組的249.16 s，增幅為39.68 s。對于波峰而言，入射N波的放大因子Am/A0逐漸從A組的4.171增加到E組的4.274；對于波谷而言，放大因子Am/A0的變化規(guī)律不是很明顯。在波峰到達底墻之前，不管是入射還是反射，越靠近底墻處波峰與波谷的放大因子Am/A0之差越大。通過觀察G01和G76測點處的波面時間序列，不難發(fā)現(xiàn)經(jīng)過底墻反射后波浪的波峰和波谷的幅值均小于入射時波浪的波峰和波谷的幅值。這是因為隨著波浪向港內傳播，港內水深均逐漸減小，波浪越往港內傳播,其發(fā)生的反射越多，波能在不斷的被消耗，故反射回來的波幅較小。

圖3 A、C、E組中A0=0.04 m時測點G01、G76和G151處的波面時間序列Fig.3 Time histories of the wave elevations at gauges G01, G76 and G151 for A0=0.04 m in Group A, C, E

3.2 港內最大波面

圖4呈現(xiàn)了各組底墻處的最大爬高隨入射波波幅的變化曲線。圖中Am表示入射波在港內的最大爬高。為了更好地區(qū)分各組實驗的最大爬高，通過A0對其進行了無因次化處理。由圖4可以看出，入射N波波幅、港內地形的平均水深以及折點的數(shù)量都對最大爬高有一定的影響。對于各組而言，隨著入射波波幅的增加最大爬高，均呈現(xiàn)先增大后緩慢減小的趨勢。從縱向來看，可以得出港內平均水深越大，則最大爬高越大，C組與F組的平均水深相等，但C組的放大因子要略小于F組，所以折點的存在將會略微減小最大爬高。

3.3 港內波能

由于正交模態(tài)分解法只有在港內波浪非線性較弱的時候才能精確分離各共振模態(tài)的響應幅值[12]。故本文所分析的入射波幅為0.01～0.07 m，間隔為0.01 m。圖5為波幅A0=0.02 m時B組和D組通過FUNWAVE-TVD模型模擬和使用正交模態(tài)分解法擬合得到的自由波面的對比。本文只研究了最低的40個共振模態(tài)。由圖5不難發(fā)現(xiàn)兩者吻合良好。B組實驗中模擬和擬合的自由波面達到最大的時刻分別為235.6 s與235.1 s，最大爬高分別為0.083 3 m和0.083 5 m，最大下沖均為-0.093 9 m。D組實驗中模擬和擬合的自由波面達到最大波面的時刻分別為255.6 s和254.7 s，最大爬高分別為0.085 1 m與0.084 8 m，最大下沖均為-0.093 7 m?？梢钥闯觯M和擬合的自由波面存在一定的誤差。考慮到這些最大爬高和下沖的特殊性，將正交模態(tài)分解法的數(shù)值擬合誤差NFE重新定義為模擬和擬合的最大爬高和下沖的最大相對誤差，即：

NFE=

(11)

式中：Au指模擬的最大爬高;(Au)fitted指擬合的最大爬高；Ad指模擬的最大下沖;(Ad)fitted指擬合的最大下沖。NFE反映了不同共振模態(tài)響應幅值的精確性。B、D各組的NEF分別為0.24%、0.35%。圖6呈現(xiàn)了所有實驗中正交模態(tài)分解法的數(shù)值擬合誤差。如圖6所示，各組實驗的NEF均隨著入射波波幅的增大而增大，這是由于波幅的增加使港內共振的非線性變強，從而使數(shù)值擬合誤差增大。這與Gao等[12]的研究結果一致。為保證不同共振模態(tài)響應波幅模擬結果的精確性，將各組實驗的NFE控制在5%以內。

圖4 各實驗組次的港內最大爬高Fig.4 Maximum runups inside the harbor for all groups

圖5 A0=0.02 m時模擬與擬合的自由波面的對比Fig.5 Comparison of simulated and fitted free wave surfaces for A0=0.02 m

本文將相對幅值定義為各共振模態(tài)響應幅值與相應的入射波波幅之比。圖7顯示了在入射波波幅分別為0.1、0.3、0.5 m條件下各組實驗最低的40個共振模態(tài)的相對幅值分布。通過對圖7的分析可以得出，當入射波幅相同時，港內地形的平均水深對港內波能的分布有一定的影響。總的來說，港內波能分布會隨著平均水深的逐漸增加而向高模態(tài)偏移，港內波能分布更加均勻。以A0=0.05 m為例，A～E組港內水深逐漸增加。A組的港內的波能主要分布在前8個共振模態(tài)上，占據(jù)最高波能的為第3模態(tài)。E組的港內的波能主要分布在前10個共振模態(tài)上，占據(jù)最高波能的為第4模態(tài)。F組為直線型地形，其平均水深與C組相同，且與其他組相比不存在折點，不難發(fā)現(xiàn)F組中各模態(tài)的相對幅值分布較其他組要略小。

港內波能的分布不僅與港內平均水深有關，與入射N波波幅的大小也有一定的關系。圖8呈現(xiàn)的是在A組、C組和E組中不同入射波幅對于港內響應幅值分布的影響，如圖所示，當入射波波幅較小時，港內的波能主要分布在最低的幾個共振模態(tài)上，只有很小的一部分波能分布在更高的共振模態(tài)上，且在較高的共振模態(tài)上，其波能的分布更加均勻。隨著入射波幅的不斷增加，波能的分布開始向較高的模態(tài)上偏移。即入射波幅越大，港內相對波能分布越均勻。以C組為例，當入射波幅較小時，港內的波能主要分布在前6個共振模態(tài)上，占據(jù)最高波能的共振模態(tài)為第2模態(tài)。

隨著入射波幅的增加，港內的波能主要分布在前21個共振模態(tài)上，占據(jù)較高波能的共振模態(tài)為第4模態(tài)和第12模態(tài)。

圖6 各組數(shù)值實驗的擬合誤差(NEF)Fig.6 Numerical fitting errors (NEF) for all cases

圖7 A～F組中港內共振模態(tài)的相對幅值分布的比較Fig.7 Comparison of relative amplitude distributions of resonance modes in group A～F

圖8 在A組、C組、E組中不同入射波幅對于港內響應幅值分布的影響Fig.8 The influence of different incident wave heights on the distribution of response amplitude in pot in group A, group C, group E

由于將港內自由波面視為各共振模態(tài)線性疊加，故港內總波能為：

(12)

式中Ai(i=1,2,…,40)表示的是第i個共振模態(tài)的響應幅值。圖9顯示了入射波波幅A0為0.01、0.03、0.05 m時各組實驗港內總波能隨著平均水深的變化情況及A～E組港內總波能的線性擬合曲線。為了更直觀觀察圖形，本文對總波能進行了無因次化處理(除以入射波幅的平方)。從圖中可以看出對于入射波幅為0.01、0.03、0.05 m時A～E組的折線型地形而言，港內總波能隨著港內平均水深增加而線性增加。C、F這2組實驗的港內平均水深相等，當入射波波幅較小的時候，平臺型地形的港內總波能要小于直線型，但當入射波幅較大時，前者反而大于后者。

圖9 A～F組實驗的港內總波能和A～E組實驗港內總波能的線性擬合Fig.9 Total wave energy in the harbor of the A～F group and linear fitting of the total wave energy in group A～E

4 結論

1)在本文所研究的地形及入射N波波幅范圍內，最大爬高隨著港內平均水深的增大而減小，隨著入射N波波幅的增大而增大，而折點的存在會略微降低最大爬高。

2)在本文研究的地形以及港內平均水深范圍內來說，港內總波能隨著平均水深的增加而呈線性增加。

3)當入射波為N波時，在港內水深相同的件下平臺型地形與直線型地形的港內總波能的相對大小與入射波幅有關。當入射波波幅較小的時候，平臺型地形的港內總波能要小于直線型，但當入射波幅較大時，前者反而大于后者。

本文的結論僅限于本文所考慮的細長形的港口形狀，港內地形和入射N波的波幅變化范圍。