張時斌， 寧德志， 陳麗芬， 張崇偉， 滕斌

(大連理工大學 海岸及近海工程國家重點實驗室，遼寧 大連 116024)

在近海區(qū)域，人們常能發(fā)現(xiàn)周期性規(guī)則沙壩[1-4]。這些沙壩的壩高最高可達1.5 m，而間距從數(shù)米到數(shù)百米不等[5]。Davies[6]通過攝動方法研究波浪與周期性沙壩相互作用時，發(fā)現(xiàn)當入射波浪波長為沙壩波長的兩倍時，連續(xù)沙壩產(chǎn)生的同相位弱反射波將在地形迎浪側疊加形成強反射波，產(chǎn)生高度共振的自由表面運動，即布拉格共振現(xiàn)象。該現(xiàn)象被Heathershaw[7]的物理實驗所證實。Davies等[6,8]的攝動方法在共振反射區(qū)，即當2k1/kb=1時(k1、kb分別為表面波和地形波的波數(shù))，求得的反射系數(shù)超過1，與物理現(xiàn)實不符。Mei[5]利用多尺度攝動法解決了Davies的這一理論缺陷。Yu等[9]基于線性淺水方程研究發(fā)現(xiàn)不僅存在2k1/kb=1的共振態(tài)，更有2k1/kb=2, 3, 4等高階共振模態(tài)的存在。Chamberlain等[10]及O′Hare等[11]通過數(shù)值計算，同樣得到了2k1/kb=2的高階結果，更是在雙周期沙壩地形(沙壩波數(shù)分別為kb1和kb2，且kb2=2kb1)條件下發(fā)現(xiàn)了2k1/kb=3,4，5時的共振現(xiàn)象。Liu等[12]采用多極展開結合極坐標變換的方法得到的解析結果顯示，斜向入射波條件下也存在高階共振模態(tài)。Zhang等[13]發(fā)現(xiàn)由于傳統(tǒng)緩坡方程忽略了高階項，局限了其在陡坡地形條件下的應用。

波浪與水下結構物相互作用時，由于淺水作用而產(chǎn)生的高階諧波[14-18]，有可能同時誘發(fā)以上多種布拉格共振模態(tài)。并且布拉格共振的發(fā)生可能會反過來改變高階諧波的特性，從而對周期性結構物周圍的波浪場造成不可忽略的影響。

鑒于布拉格共振和高階諧波都具有重要的工程意義，并且目前缺少對這2種現(xiàn)象的聯(lián)合研究，本文將在前人的研究結果基礎上，利用時域高階邊界元方法，在布拉格共振條件下，進行高階諧波的特性研究。

1 布拉格共振的數(shù)值計算方法

在笛卡爾直角坐標系下考慮波浪與水下正弦形周期性地形(沙壩)相互作用的問題。笛卡爾直角坐標系Oxz的z軸向上為正，x軸向右為正(與入射波浪傳播方向同向)。設置靜水面為z=0平面，平底上方靜水深為h。水槽底部布置波幅為b，波長為λb的n周期(下文稱n為沙壩相對長度)正弦形沙壩結構ζ(x)：

(1)

式中：x1表示沙壩結構的起點坐標；kb為正弦形沙壩結構的波數(shù)。

假定流體為無粘、不可壓縮且運動無旋的理想流體，引入勢流理論來描述水槽內的流體運動。在數(shù)值模型中，本文采用內嵌源造波技術產(chǎn)生入射波浪，并且在水槽兩端分別布置長為Lb的數(shù)值阻尼層以消除來自出流邊界的波浪反射及造波邊界的波浪二次反射，因此數(shù)值模型的控制方程和初邊值條件為：

(2)

(3)

(4)

(5)

φ(x,z,t)=0，η(x,t)=0，t≤0

(6)

(7)

式中：ν為阻尼系數(shù)(=1)；Lb為阻尼層寬度，取為2倍波長；x01和x02分別為水槽兩端阻尼層的起始位置。

Davies等[19]用水槽底部流體粒子水平位移幅度2Ab與地形波長λb的比值2Ab/λb(=2ga0k1/(λbω2cosh(k1h)，其中g為重力加速度，h為左端水平地形上的靜水深，a0/k1/ω1分別為入射波的波幅/波數(shù)/角頻率)來判別正弦波地形上流體的流動狀態(tài)。當該參數(shù)～O(10-1)時，認為此時的流分離現(xiàn)象不嚴重，即流體粘性作用可以忽略。本文考慮的工況中，2Ab/λb最大為0.3～O(10-1)，因此，基于勢流理論建立數(shù)值模型用以模擬研究波浪與周期性地形相互作用是合理的。同樣，本文主要關注非破碎陡波在周期性地形上的傳播變形；波浪破碎等會對波浪非線性、布拉格共振及高階諧波特性產(chǎn)生一定的影響，但不在本文考慮的范圍內，是本文下一步需要深入研究的內容。

本文采用高階邊界元方法求解上述方程(2)～(6)，具體求解方法詳見文獻[20]。數(shù)值計算中，通過收斂性分析，水槽自由表面網(wǎng)格尺寸和沙壩結構網(wǎng)格尺寸分別選為Δx=λ/30和Δx=λb/10，時間步長為Δt=T/60 (λ和T分別為表面波的波長和周期)。

2 布拉格共振對高階諧波特性的影響

為了證明本文結果的準確性，將本文模型結果和已發(fā)表的實驗結果及理論解進行對比。本文依照文獻[7]及文獻[19]的實驗參數(shù)進行如下設置：沙壩波幅b=0.05 m，沙壩波幅與水深之比b/h=0.16，沙壩波長λb=1 m，沙壩相對長度n=10，入射波幅a0=0.01 m。

利用四點處(堤前1.6 m，間隔0.1波長布置)的波高信息進行入反射分離得到堤前的反射系數(shù)(具體的分離方法參見文獻[21])，其隨入射波數(shù)的變化規(guī)律見圖1?？梢园l(fā)現(xiàn)，本文模型結果與理論解吻合良好，但在布拉格共振條件附近與實驗結果相比，高估了共振反射系數(shù)。

圖1 堤前反射系數(shù)隨入射波數(shù)的變化Fig.1 Reflection coefficients as a function of the wavenumber

圖2給出了布拉格共振條件下，反射系數(shù)的空間分布(沿水槽中線布置一系列浪高儀，同樣利用四點法進行入反射波分離)。同樣發(fā)現(xiàn)，在沙壩的迎浪側，本文模型結果相對于實驗結果高估了反射系數(shù)。反射系數(shù)的高估可能是實驗結果受到來自造波板二次反射的影響，且這也可能是造成地形左端水平地形上反射系數(shù)沿程變化的原因。而在沙壩的背浪測，實驗結果則可能受到水槽右側消波區(qū)反射波(實驗消波區(qū)未能吸收所有透射波)的影響而產(chǎn)生波動。

圖2 布拉格共振下反射系數(shù)的空間變化Fig.2 The spatial distribution of the reflection coefficient under the Bragg resonant condition

將數(shù)值模型的源造波法轉換為用速度勢入口法進行造波(速度勢同樣根據(jù)二階斯托克斯方法給定)，則此時行進波在遇到結構物(周期性地形)后產(chǎn)生的反射波會在入口(造波板)處產(chǎn)生二次反射。圖1、2中也給出了考慮二次反射的計算結果(黑虛線)，可以發(fā)現(xiàn)，新的結果和實驗結果十分吻合，由此證明了上述猜想以及本文數(shù)值模型模擬布拉格共振現(xiàn)象的可靠性。

2.1 沙壩波幅的影響

所建立/驗證的模型現(xiàn)用于模擬研究沙壩地形參數(shù)對布拉格共振和高階自由諧波特性的影響規(guī)律。

圖3給出了3組不同沙壩波幅條件下，反射系數(shù)隨入射波數(shù)的變化情況。圖中虛線指示出布拉格共振的理論位置(2k1/kb=1)。可以發(fā)現(xiàn)反射系數(shù)在2k1/kb=1附近取得最大值，并且其峰值和帶寬均隨沙壩波幅增加而增加(峰值最大值為1)。另外，共振峰值頻率相對于理論共振頻率(2k1/kb=1)出現(xiàn)低頻偏移現(xiàn)象，且偏移程度隨波幅增加而增加。該趨勢與Liu等[22]據(jù)攝動理論計算給出的低頻偏移現(xiàn)象一致。另外在高階布拉格共振條件2k1/kb=2處，發(fā)現(xiàn)隨著沙壩波幅b/h增加，反射系數(shù)逐漸呈現(xiàn)局部峰值。這是由于該共振對應的波頻較高，波浪能量主要分布于水體表層，只有當沙壩波幅較大時表面波才能受到沙壩地形調制從而產(chǎn)生共振。

圖3 不同沙壩波幅下反射系數(shù)隨入射波數(shù)變化Fig.3 Reflection coefficients as a function of the wavenumber ratio for the three bar amplitudes

為了研究布拉格共振條件(2k1/kb=1)下，沙壩波幅對各階諧波的影響，本文分別采用兩點法[23]和四點法[21]分離沙壩背浪側和迎浪側的波浪組分。圖4給出了沙壩迎/背浪側的各階波成分透反射系數(shù)隨相對沙壩波幅b/h的變化。圖例中上標r和t分別表示反射波和透射波；下標f和l分別表示自由波和鎖定波；下標(1)和(2)分別表示一階波浪和二階波浪。

從圖4可以看出，隨著沙壩波幅的增加，基頻波的共振反射系數(shù)逐漸增強并最終收斂至1。這一趨勢與圖3呈現(xiàn)的規(guī)律一致。在圖4(b)和(c)中，二階鎖定波的透射系數(shù)和反射系數(shù)具有相反的變化趨勢，即隨著沙壩波幅增加，反射的二階鎖定波能量增加而透射的二階鎖定波能量減少。為了便于分析，圖4(c)中用虛線標示了平底條件下測得的二階鎖定波無量綱波幅0.045(可視為入射波中的二階鎖定波成分的無量綱波幅)?？梢园l(fā)現(xiàn)，圖中二階鎖定波的反射系數(shù)最終收斂至略高于此虛線處。因為二階鎖定波的波數(shù)為2k1，其同樣滿足圖3中發(fā)現(xiàn)的高階布拉格共振條件(2k1/kb=2)，所以可以推測布拉格共振條件下，二階鎖定波也將發(fā)生共振反射。而由于淺水效應誘導的基頻波能量轉移，其最終反射系數(shù)略高于入射二階鎖定波的波幅。此外二階自由波透射系數(shù)和反射系數(shù)均隨著沙壩波幅的增加而增加，前者在b/h=0.224時達到約0.05的飽和值，然后緩慢降低，后者則持續(xù)增加。

圖4 布拉格共振條件下各階諧波透/反射系數(shù)隨相對沙壩波幅變化Fig.4 Variation of the reflection coefficient of the first two harmonics with the relative bar amplitude under the Bragg resonant condition

圖5給出了布拉格共振條件(2k1/kb=1)下，不同沙壩波幅條件(b/h=0，0.16，0.32，0.48)下各階諧波幅值的沿程變化(圖中虛線指示沙壩左右邊界)。水槽內每點的各階諧波幅值是通過對空間內各點波面的時間歷程進行快速傅里葉變換(FFT)得到?？梢栽趫D5(a)中清楚地發(fā)現(xiàn)在沙壩迎浪側和上方，一階諧波幅值出現(xiàn)了以1 m為空間周期的振蕩。這是由一階入射波和反射波線性疊加形成駐波所致(駐波理論節(jié)間距也為π/k1=1 m)。布拉格共振反射使沙壩上方的一階透射波和反射波波幅沿程線性減小，進而導致一階諧波幅值的振蕩幅度沿沙壩長度方向線性減小并在沙壩末端收斂為一定值(即沙壩背浪側的透射波波幅)。

圖5 布拉格共振條件下各階諧波幅值沿程變化Fig.5 Spatial distribution of the wave amplitude of the first two harmonics under the Bragg resonant condition

有趣的是二階諧波波幅的空間變化，首先在圖5(a)中b/h=0的工況下，一階諧波波幅沿水槽基本恒定。這說明在該工況下沒有發(fā)生3波相互作用[24](3波相互作用將導致基頻波和二階波浪之間發(fā)生能量交換)。其次，可以在圖5(b)中發(fā)現(xiàn)沙壩背浪側的二階諧波幅值出現(xiàn)了被稱作“拍”的周期性振蕩。譬如，在b/h=0.16時，實際測量所得的拍長為1.975 m，與理論拍長公式所得的結果1.939 m (=2π/(k2-2k1))較為接近。對比不同沙壩波幅條件下沙壩背浪側的“拍”現(xiàn)象，可以發(fā)現(xiàn)沙壩波幅大小不改變拍長,但是會導致拍在空間上出現(xiàn)相位偏差。

在沙壩迎浪側，二階諧波幅值呈現(xiàn)出更復雜的空間變化(體現(xiàn)為在拍長1.975 m的主拍基礎上出現(xiàn)數(shù)個寄生拍)，并且幅值變化的劇烈程度隨沙壩波幅增加而增加。

對于“拍”的形成，Hansen等[25]給出了數(shù)學解釋，即在平底條件下，水槽中二階波面(忽略初始相位)通過三角函數(shù)變換可以寫為：

(8)

式中第2項屬于波群效應，其幅值隨空間振蕩但不隨時間變化。波群成分的空間振蕩周期(最大波幅值或最小波幅值的間距)，就是上文提到的拍長L，其滿足公式(k2-2k1)L/2=π,即拍長L=2π/(k2-2k1)。

從上述推導中可以發(fā)現(xiàn)，“拍”是二階鎖定波和二階自由波之間相位不匹配導致的空間波幅周期性變化，因此如果缺少其中一種波成分，二階諧波的“拍”現(xiàn)象將會消失而變?yōu)橐欢ㄖ?。例如b/h=0.48的工況下，在圖4(b)可以發(fā)現(xiàn)沙壩背浪側二階鎖定波波幅較小，導致在圖5(b)中相應的幅值變化曲線近似于水平直線。

與式(8)類似，沙壩迎浪側的二階波面(忽略初始相位)通過三角函數(shù)變換，可以寫為：

(9)

式(9)顯示沙壩迎浪側的波成分可以轉變?yōu)樾羞M波、具有波群特性的波和2個駐波的組合。前兩者與式(8)相同，而后兩者(即駐波)是新增加的成分。新駐波波幅的空間變化頻率顯著高于波群的波幅，這與圖5(b)中的寄生拍的高頻空間振蕩現(xiàn)象一致，因此推斷新駐波成分是寄生拍的成因。

當沙壩背浪側某點的各階波浪處于相同相位時，將形成較大波高的自由表面波浪，進而可能對附近的海岸結構物造成重要影響。因此，對于各階波浪的相位進行研究分析是必要的，因為各階鎖定波處于鎖相狀態(tài)(同相位)，所以只需對各階自由波相位展開研究。

圖6顯示布拉格共振條件(2k1/kb=1)下，沙壩背浪側各階自由波的初始相位隨沙壩相對波幅b/h的變化。其中，初始相位是基于沙壩背浪側G6(堤后1.6波長處)和G7(G6后0.1波長處)號探針處的波面時間歷程，利用兩點法得到。可以發(fā)現(xiàn)，自由波相位隨著沙壩波幅增加而線性增加。通過線性擬合可以得到一、二、三階自由波的初始相位變化斜率分別為1.65、3.60、4.85。其中，后兩者分別約為前者的2倍和3倍。利用線性擬合直線，可以預測不同潛堤波幅條件下布拉格共振時堤后各階自由波的初始相位，進一步可以確定堤后波高最大處，即各階波浪同相位點。

圖6 布拉格共振條件下沙壩背浪側各階自由波初始相位隨地形幅值的變化Fig.6 Variation of the initial phases of the free waves behind the submerged dikes with the bar amplitude under the Bragg resonant condition

2.2 沙壩長度的影響

為了增強高階諧波的能量，下文結果集中于相對沙壩波幅為0.32的算例。圖7給出了沙壩相對長度n分別為2、5、10時，基頻波反射系數(shù)隨入射波數(shù)變化的結果?？梢园l(fā)現(xiàn)，隨著n變大，布拉格共振峰值和次峰數(shù)量增加，而共振帶寬減小。

圖7 基頻波反射系數(shù)隨入射波數(shù)變化Fig.7 Reflection coefficients of the fundamental wave as a function of the wavenumber ratio

圖8給出了實測共振峰值頻率ωp隨沙壩相對長度n的變化?？梢园l(fā)現(xiàn)隨著n增大，共振峰值頻率向著高頻方向偏移，并且收斂于ωp=4.72附近。這一趨勢與Wen等[26]得到的規(guī)律一致。根據(jù)布拉格共振的物理機制，即相鄰沙壩產(chǎn)生同相位的反射波，則布拉格共振的準確數(shù)學條件可以表達為：

(10)

式中第1/2項表示入/反射波在單位波長沙壩上傳播的相位變化;k1(x)與入射波頻ω1之間滿足線性色散關系ω12=gk1(x)tanhk1(x)[h-ζ(x)]，其中ζ(x)的表達式見式(1)。

圖8 布拉格共振峰值頻率隨沙壩相對長度的變化Fig.8 Variation of the Bragg resonant frequency with the relative bar length

式(10)中取m=1時，即為布拉格共振條件(2k1/kb=1)。經(jīng)計算，相應的理論共振頻率為4.75 rad/s(通過該方法無法考慮沙壩長度的影響，其可被視為沙壩無窮長工況的結果)。該值相比未考慮水底波動地形得到的理論共振值[22]4.82 rad/s，更接近于圖8中共振峰值頻率的近似收斂值(ωp=4.72)。前者相對誤差較于后者減小1.5%，即僅為0.6%。

經(jīng)過分析，式(10)與數(shù)值結果不同的原因有2點：1)圖8中獲得的收斂頻率并非為精確的無窮長地形條件下的共振頻率(后者略大于前者)；2)因為波浪無法隨著水深變化瞬時改變其波數(shù)，而上述推導中默認了波數(shù)的瞬變假設。該結果雖然依舊存在誤差，但是相對于原理論共振頻率[22]，一定程度提高了長沙壩工況下共振峰值頻率的預測準確度。

(11)

式(10)可以改寫為：

(12)

接著在前文布置的基礎上，研究布拉格共振條件(2k1/kb=1)下，沙壩相對長度n對各階諧波的影響。圖9分別給出了基頻波和二階諧波的透反射系數(shù)隨相對長度n的變化。

圖9 布拉格共振條件下各階諧波透/反射系數(shù)隨沙壩相對長度變化Fig.9 Variation of the reflection and transmission coefficients of the first two harmonics with the relative length under the Bragg resonant condition

可以發(fā)現(xiàn)，圖9中基頻波與二階鎖定波透反射系數(shù)的整體趨勢與隨沙壩波幅變化的相應結果(圖8)基本一致，這意味著沙壩相對長度n和沙壩相對波幅b/h對這兩種波成分在布拉格共振條件下的影響類似(此處不再贅述)；二階自由波透反射系數(shù)均隨n增大而增大，并最終各自穩(wěn)定于一定值附近振蕩。值得注意的是隨n變化的透反射率振蕩現(xiàn)象在多個成分波浪中出現(xiàn)，其中二階透射自由波尤為顯著，這可能與圖10(b)中的二階自由諧波初始相位震蕩現(xiàn)象有關。

圖10給出了布拉格共振條件(2k1/kb=1)下堤后各階自由波初始相位隨沙壩相對長度n的變化?？梢园l(fā)現(xiàn)二階自由波初始相位φ2以約為π的幅值，整體呈線性緩慢減小的趨勢。該相位的振蕩變化可能與圖9中二階諧波透射系數(shù)的振蕩有關。此外基頻波初始相位φ1初始線性增長；當沙壩長度繼續(xù)增加，φ1出現(xiàn)小幅振蕩并且增長率變緩。

圖10 布拉格共振條件下沙壩背浪側各階自由波的初始相位隨地形長度的變化Fig.10 Variation of the initial phases of the free waves behind the submerged dikes with the bar length under the Bragg resonant condition

3 結論

1)隨著沙壩波幅增加，布拉格共振峰值頻率向低頻偏移，共振帶寬變寬，共振反射系數(shù)增加并收斂至1。

2)隨著沙壩長度增加，布拉格共振峰值頻率向高頻偏移，共振帶寬變窄，共振反射系數(shù)增加并收斂至1。

3)在布拉格共振條件下，二階鎖定波發(fā)生共振反射，并且反射系數(shù)隨著沙壩波幅或沙壩長度而增加并最終收斂于入射二階鎖定波無量綱波幅附近。

4)隨著沙壩波幅增加，二階自由波透射系數(shù)先增大達到約0.05的飽和值，然后緩慢減小。二階自由波反射系數(shù)則持續(xù)增大。

5)隨著沙壩長度增加，二階自由波透反射系數(shù)均增大，并最終各自穩(wěn)定于一定值附近振蕩。

6)布拉格共振條件下，在沙壩背浪側的二階諧波波幅會呈現(xiàn)“拍”的現(xiàn)象，而在沙壩的迎浪側，除了“拍”以外，還存在附著于“拍”上的“寄生拍”，增加了極端大波出現(xiàn)的概率。