有限長正弦地形上的波浪傳播問題研究

滕 斌，崔 杰

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

有限水深下波浪傳播過程會受到海底地形的影響。周期性的沙紋和沙壩是常見的海岸地形，這些周期地形在很大程度上會影響水波在海岸上的傳播。因此，探究周期地形對水波的作用是非常有意義的。

Davies[1]應用攝動展開法研究了有限長正弦周期地形對水波的反射，發(fā)現當水波波長是地形波長的兩倍時波浪反射率會顯著增加。Heathershaw[2]通過物理模型試驗證實了這一現象。因其類似于固體物理中的布拉格共振現象，故將該現象稱為水波的布拉格共振。布拉格共振現象有助于解釋沙紋和沙壩的生成機理[3-4]，一些文獻[5-6]還提出了運用人工潛堤來防護海灘的建議。

為探究水波的布拉格共振現象，許多方法被用來研究周期地形對水波的反射問題。對于規(guī)則的周期地形，Davies和Heathershaw[7]利用關于地形的攝動展開方法研究了小幅度起伏周期地形對水波的反射和透射問題。該理論可以預測布拉格共振頻率，但其計算的反射率在共振點處大于模型試驗值，同時存在波能流不守恒的缺陷。Mei[3]利用多尺度展開法研究了沙壩誘發(fā)的布拉格共振現象，計算的反射率與模型試驗結果吻合良好。Davies等[8]還利用馬修方程研究了長波在周期地形上的傳播問題，發(fā)現當海床波數近似等于水面波數的兩倍時，水波與海床發(fā)生強烈的相互作用，這時反射系數最大。

對于任意變化的地形，可采用數值方法進行研究。O’Hare和Davies[9]、Cho和Lee[10]將連續(xù)的曲面海底地形近似為高度變化的臺階地形，然后將臺階上的速度勢按特征函數進行展開，運用傳輸矩陣和旋轉矩陣聯立相鄰臺階的波浪場，以此求解了水波在周期地形上的傳播問題，計算了周期地形對水波的反射。Kirby[11]推導出一種擴展的緩坡方程，以此研究了帶有坡度的海底起伏地形上的水波現象。Hsu等[12]利用一種新的緩坡方程方法(EEMSE)，研究了人工系列潛堤引起的水波布拉格共振現象。Liu等[13]采用修正緩坡方程方法(MMSE)研究了正弦地形上的波浪反射問題。Dalrymple等[14]采用邊界元方法計算了Davies和Heathershaw[1, 7]試驗中周期地形對水波的反射，研究了正向和斜向波浪入射下周期地形上的反射系數，該方法不受海底地形坡度和起伏高度的限制，計算結果與Mei[3]的多尺度展開法和試驗結果吻合良好。Liu等[15]采用高階譜方法研究了布拉格共振下的波—波間的相互作用，給出了布拉格共振發(fā)生的三種波浪與海底間的非線性相互作用條件。

下文應用基于勢流理論的頻域高階邊界元方法[16]，計算了正弦周期地形對水波的作用。該方法以拉普拉斯方程為控制方程，對地形不做任何近似簡化，準確描述海底地形變化，可應用于任意的海底坡度和起伏高度。應用該方法系統(tǒng)地研究了正弦周期地形個數、地形高度對反射系數的影響，以及對應于高反射率和低反射率波浪頻率下波面高度的沿程分布，發(fā)現了一些新的波浪特征和現象。

1 數學模型與方法

考慮水波在多個正弦周期地形上的傳播問題，如圖1所示。建立笛卡爾坐標系，x軸位于靜水面上，z軸垂直向上，x坐標零點位于周期地形的起始位置。地形的前后兩側為水平海底，水深分別為d1和d2，海底方程為

式中：局部地形幅度為a，周期地形波長為λ，總的地形長度為XR。為描述方便，圖1中只繪出了兩個周期的地形。假設流體不可壓縮且無黏無旋，有速度勢函數Φ(x,z,t)存在。在小波幅近似下，速度勢滿足的控制方程和線性邊界條件：

?2Φ(x,z,t)=0，f(x)

Φ=Re[φ(x,z)e-iωt](5)

復速度勢φ(x,z)滿足的控制方程和邊界條件為

?2φ(x,z)=0，f(x)

采用二維邊界元模型[16]求解該問題。取圖1中所示的流域為研究對象，其中SL、SR、SF、SB分別表示上游輻射邊界、下游輻射邊界、水面邊界和海底邊界。上、下游輻射邊界SL、SR離開海底地形的距離為LH。

圖1 水底地形示意Fig. 1 Sketch of two periods of underwater terrain

在左端立面SL處分解速度勢為入射勢φI與反射勢φR：

式中：A為入射波波幅，k1為水深d1處的波數，即色散關系ω2=gk1tanh(k1d1)的實根；μ1j(j=1,2,……)為色散方程ω2=-gμ1jtan(μ1jd1)的實根；Rj(j=0,1,2,……)為待定系數，波浪反射系數為KR=|R0|。

取SL遠離海底起伏段，可忽略掉反射速度勢的非傳播模態(tài)，而將其寫為：

φR(x,z)=R0φ1(x,z)(11)

由式(9)和(11)，可得SL處的邊界條件：

同樣，取SR遠離海底起伏段，忽略非傳播模態(tài)，可將透射速度勢近似表示：

φT(x,z)=T0φ2(x,z)(15)

式中：μ2j為色散方程ω2=-gμ2jtan(μ2jd2)的實根，Tj(j=1,2,……)為待定系數。由式(15)，可得SR處的邊界條件：

選取Rankine源為格林函數：

式中：x=(x,z)和x0=(x,z)分別為格林函數的場點和源點。對格林函數和速度勢應用格林第二定理，得到邊界積分方程：

式中：α為固角系數。將邊界條件代入式(19)后，可得：

我很吃驚。現在很多家庭都只有一個孩子，只要求孩子讀好書，其他的事情基本都由家長包辦了。慢慢地，孩子們覺得父母的付出都是理所當然的，這樣的心理造成孩子沒有感恩父母之心。

將邊界離散為三節(jié)點單元。使用多項式函數插值三節(jié)點上的速度勢，然后運用高斯積分法求式(20)中的積分。未知量為：水面和水底的N個節(jié)點處的速度勢，另外加上反射系數R0和透射系數T0。將源點分別布置在水面和水底的各個節(jié)點上，以及左、右邊界的任一節(jié)點上，可建立(N+2)個線性方程：

由此可求得水面、海底邊界節(jié)點處的速度勢，以及反射系數R0和透射系數T0。

2 數值研究

2.1 數值模型的驗證

首先以圖1所示的地形為例，檢驗邊界單元數量和計算域尺寸對計算結果的影響。水底地形由三部分組成：中間為正弦周期地形，前后兩側為等深的水平地形。正弦地形的平均水深和前后水平地形的水深均為d=0.312 5 m，正弦地形波幅為a=0.2 m，正弦地形波長為λ=1 m，正弦地形總長度為XR=10 m(波形個數為n=10)。

表1是單元尺度對反射系數KR影響的研究。計算中采用了四種網格，網格單元的密度分別為：網格一，30單元/m；網格二，60單元/m；網格三，120單元/m；網格四，240單元/m。正弦地形前后水平地形的長度與水深比均為LH/d=16?？梢园l(fā)現，當kd≤2.5時由網格二和網格三計算的結果與網格四的結果直到小數點后第四數字都是相同的。而當kd>2.5時，計算結果的相對誤差則明顯增加。誤差增加的主要原因是：隨kd增加，水波波長減小，在網格尺寸不變的情況下計算精度下降。文中，kλ/π的最大取值為2.5，水深為d=0.312 5 m，正弦地形波長為λ=1 m，此時kd<2.5。所以下文的算例均使用網格二。

表1 各網格下的波浪反射率 (a/d=0.64, n=10)

表2是計算域大小對反射系數影響的研究。計算中分別選取了水平地形長度與水深比LH/d=1.0,2.0,5.0和10.0的四種情況。網格密度采用上述網格二的劃分方式，即60單元/m?？梢园l(fā)現，當kd≤2.5時，LH/d=5下的計算結果與LH/d=10的結果直到小數點后四位數字都是相同的；當kd>2.5時，LH/d=5的計算結果與LH/d=10的結果相比直到小數點后三位數字都是相同的。在本文的計算中，考慮到要分析水平地形上的水波，故將前后水平地形的長度設置為LH=5 m(LH/d=16.0)，此時水平地形長度已經可以滿足計算精度的需求。

表2 水平段長度對波浪反射系數的影響 (a/d=0.64, n=10)

最后，將本模型結果與Heathershaw[2]的物理模型試驗結果做了對比，以驗證本模型的正確性。正弦地形的平均水深和前后水平地形的水深為d=0.312 5 m，正弦地形波幅為a=0.05 m，正弦地形波長為λ=1 m，正弦地形總長度為XR=10 m(波形個數為n=10)。計算網格采用上述網格二的方式劃分，地形前后的水平段長度均為LH=5 m。圖2為計算與試驗結果的對比。可見，本文邊界元法計算的結果與物理模型試驗吻合良好，可用于研究周期地形上的波浪傳播問題。

圖2 與Heathershaw物理模型試驗結果的對比 (a/d=0.16, n=10)Fig. 2 Comparison with Heathershaw’s model test (a/d=0.16 and n=10)

2.2 正弦地形幅度對波浪反射的影響

圖3是正弦地形幅度對波浪傳播影響的對比研究。正弦地形的平均水深和前后水平地形的水深為d=0.312 5 m，正弦地形波長為λ=1m，正弦地形總長度為XR=20 m(波形個數為n=20)，正弦地形波幅分別為a=0.05 m,0.10 m和0.20 m (a/d=0.16,0.32,0.64)。

圖3 地形幅度對波浪反射的影響 (n=20)Fig. 3 Influence of sea bed amplitude on the wave reflection (n=20)

圖3中同時繪制了按照Mei[3]的理論解計算得到的結果。從圖中可以看到，當地形波幅較小a/d=0.16時(圖3(a))，BEM的結果與Mei的理論解吻合良好；而當波幅增大為a/d=0.32和0.64時(圖3(b)、(c))，Mei的理論解無法計算出布拉格共振帶下移、以及第二共振區(qū)的結果。由此可以發(fā)現，對于地形波幅較大的情況，Mei的理論解有一定限制。而本文的BEM法可以用于計算任意的水底地形，對地形波幅的大小沒有限制。

2.3 正弦地形個數對波浪反射的影響

圖4是正弦地形個數對水波反射率的影響。正弦地形的平均水深和前后水平地形的水深均仍為d=0.312 5 m，正弦地形波長為λ=1 m，正弦地形波幅為a=0.2 m(a/d=0.64)。圖4分別為波浪在長度分別為5 m、10 m和20 m正弦地形(波形個數為n=5、10和20)上反射系數隨相對波數kλ/π的變化關系。從三種情況的對比中可以看出，在第一個布拉格共振帶附近，反射系數的外包絡基本一致，第一個布拉格共振帶的寬度和反射系數大小基本上不隨海床地形長度而變化，但反射系數波動周期隨海底地形長度而變化。海底地形較長(個數較多)的情形，反射系數的波動周期較短，反射系數的波動周期與海底周期地形長度成倒數關系(∝1/XR)。第二個布拉格共振帶的寬度隨海床地形長度的減少而變寬，反射系數隨海床地形長度的減少而減小。

圖4 正弦地形個數對波浪反射的影響Fig. 4 Influence of sea bed terrain number on the wave reflection

2.4 波幅的沿程分布

對波浪在正弦地形上傳播過程中的波幅分布進行分析。波浪從左側入射，正弦地形的平均水深和前后水平地形的水深均為d=0.312 5 m，正弦地形波幅為a=0.2 m(a/d=0.64)，正弦地形波長為λ=1 m，正弦地形總長度為XR=20 m(波形個數為n=20)。設入射波的波幅為A0，沿程的波動幅值為A(x)。

2.4.1 布拉格共振區(qū)域內的波浪傳播

圖5為kλ/π=0.860位于第一共振區(qū)域內，波幅A(x)/A0的沿程變化。如圖3(c)所示，該頻率下反射率為1。因反射系數為1，在上游的水平地形上形成了由入射波與反射波疊加而成的駐波，其波節(jié)的間距為入射波波長的一半；在下游地形上，波幅為零，沒有波浪通過。在正弦地形上，波幅A(x)/A0也呈立波狀，存在多個波幅為零的節(jié)點，波幅A(x)/A0的包絡沿程從1衰減到0。同時應該注意到，波幅在周期地形上并非呈線性衰減。

圖5 kλ/π=0.860下波面幅值的沿程分布Fig. 5 Distribution of wave amplitude along the x coordinate (kλ/π=0.860)

2.4.2 布拉格共振區(qū)域外高反射區(qū)的波浪傳播

圖6是布拉格共振區(qū)域外高反射區(qū)kλ/π=0.696,0.721,1.007,1.035處，A(x)/A0的沿程變化。kλ/π=0.696和0.721為第一布拉格共振區(qū)左側的第二和第一個反射率為峰值的位置；kλ/π=1.035和1.007為第一布拉格共振區(qū)右側的第一和第二個反射率為峰值的位置，如圖3(c)所示。

圖6 波面幅值沿水平空間x的分布 (kλ/π=0.696, 0.721, 1.007, 1.035)Fig. 6 Distribution of wave amplitude along the x coordinate (kλ/π=0.696, 0.721, 1.007, 1.035)

可以發(fā)現，由于反射系數小于1，左側水平地形上的波浪為入射波浪與不完全反射波浪的疊加，右側水平地形上為波幅減小了的穩(wěn)定透射波浪。在正弦地形上，波幅A(x)/A0除了自身波動外，還存在波動的包絡。在左、右第一個反射率峰值處，存在1個包絡；而在左、右第二個反射率峰值處，存在2個包絡。

2.4.3 布拉格共振區(qū)域外低反射區(qū)的波浪傳播

圖7是布拉格共振區(qū)外低反射區(qū)kλ/π=0.708,0.729,0.999,1.021處，A(x)/A0的沿程變化。kλ/π=0.708和0.729為第一布拉格共振區(qū)左端第二和第一個反射率為零的位置；kλ/π=0.999和 1.021為第一布拉格共振區(qū)右側第一和二個反射率為零的位置(如圖3(c)所示)。

圖7 波面幅值沿水平空間x的分布 (kλ/π=0.708, 0.729, 0.999, 1.021)Fig. 7 Distribution of wave amplitude along the x coordinate (kλ/π=0.708, 0.729, 0.999, 1.021)

在低反射區(qū)，反射系數接近于零。左側和右側水平地形處的波幅近似為1，正弦地形上(0≤x≤20)波幅A(x)/A0呈現出包絡形狀，并在包絡線內振蕩。在第一個布拉格共振區(qū)左、右側第一反射率為零的頻率處，出現1個包絡，最大波幅為入射波幅的6倍左右。而在第一個布拉格共振區(qū)左、右側第二反射率為零的頻率處，最大波幅為入射波幅的3倍左右。

當入射波波幅不是很小時，上述的周期地形上的波幅增大現象可能導致強烈的水面波動，以及波浪的破碎，實際計算中需加以注意。

3 結 語

利用基于勢流理論的頻域邊界元方法，對正弦地形上的波浪傳播問題做了數值計算研究，得到如下結論：

1) 邊界元方法可以被用來研究任意周期地形上的水波傳播問題，對水底地形的波動幅度沒有限制。當水底地形幅度較小時，邊界元方法得到的結果與Mei的理論解較為吻合。隨海底地形波動幅度的增加，邊界元方法得到的布拉格共振區(qū)域明顯向低頻移動，這是Mei的理論所無法準確預測的。

2) 海底地形長度(或個數)對反射系數的外包絡影響不大，但對反射系數的振蕩周期有顯著的影響，反射系數的振蕩周期與海底地形長度成倒數關系。

3) 在非布拉格共振區(qū)域，周期地形上的波浪高度沿程空間分布呈現包絡形狀。當頻率位于布拉格共振區(qū)域邊緣且反射率為零時，包絡中的最大波幅遠大于入射波波幅。在文中的算例計算中，最大波幅可達到入射波幅的6倍左右。