彭延建，劉應(yīng)中，時鐘

（上海交通大學(xué) 船舶海洋與建筑工程學(xué)院 港口與海岸工程系，海洋工程國家重點(diǎn)實(shí)驗(yàn)室，上海 200030）

波浪在不平坦海岸底床上傳播的平面二維數(shù)學(xué)模型

彭延建，劉應(yīng)中，時鐘*

（上海交通大學(xué) 船舶海洋與建筑工程學(xué)院 港口與海岸工程系，海洋工程國家重點(diǎn)實(shí)驗(yàn)室，上海 200030）

基 于線性波 浪 理論，利用攝動 理 論 中 的 多 重 尺度 法 ， 以坡度 ε 為小 參 數(shù) ， Liu 和 Shi[1]推導(dǎo)出 波 浪 在不 平 坦底床傳播的波幅方程，以及一階勢函數(shù)的解析解。為驗(yàn)證該模型是否可用于平面二維波浪傳播問題，利用有限差分法，分別對波幅方程和邊界條件進(jìn)行離散和求解，得到因折射、繞射引起的波幅和波向的變化。通過引入繞射因子，使該 模型的應(yīng)用范圍更為廣泛，并且計算精度得到提高。將計算得到 的波向和波 高與 Berkhoff[2]、Berkhoff 等[3]橢圓形淺灘試驗(yàn)值進(jìn)行比較，表明該模型可以應(yīng)用于平面二維問題。

多重尺度法；波浪傳播；不平坦底床

0 引言

波浪在不平坦海岸底床上的傳播是波浪研究的重要課題之一。波浪在從深水區(qū)向近岸傳播過程中，因?yàn)楹０兜匦巫兓戎T多因素，會發(fā)生淺化、折射與繞射等現(xiàn)象，波況也隨之變化，比如：波長變短、波速變慢、波高先減小后增大、波向也會趨于與岸線垂直。本文擬對線性波浪在不平坦海岸底床上的傳播進(jìn)行初步研究，對海岸工程學(xué)科具有科學(xué)意義。

對波浪在不平底床傳播問題的研究[4-9]，大體有兩種研究方法：

第一種是基于緩坡方程和 Boussinesq 方程等已經(jīng)較為成熟的模型，尋求數(shù)學(xué)模型和數(shù)值求解方法的改進(jìn)，比如在緩坡方程中引入底摩擦、波浪破碎、非線性、波流相互作用等效應(yīng)，以及對緩坡方程的簡化和近似，例如雙曲形緩坡方程、拋物形緩坡方程等，國內(nèi)外眾多學(xué)者在這方面的努力使基于緩坡方程和 Boussinesq 方程的模型在計算精度、計算速度、穩(wěn)定性等方面有了較大的提高，已經(jīng)能夠較為成功地應(yīng)用于某些工程案例。這類模型也存在一些缺點(diǎn)，比如有的只適合大面積計算區(qū)域，有的只適合小面積計算區(qū)域，求解短周期波浪場需要較高的網(wǎng)格密度，需要較大的計算機(jī)內(nèi)存，求解過程也較為耗時等。

第二種是理論解析的方法，例如通過直接求解底部邊界條件、線性自由表面條件來推導(dǎo)波浪勢函數(shù)的解，或通過攝動展開求波浪勢函數(shù)的近似解。這類模型存在的問題是，給出的解析解形式一般是局限在某些特定條件下，例如二維的斜坡床底等特定水深地形。

根據(jù)以上論述，對波浪在不平底床傳播問題的研究，多數(shù)研究還是基于緩坡方程和Boussinesq 方程，進(jìn)行模型的改進(jìn)和數(shù)值方法完善，而在進(jìn)一步發(fā)展波浪在不平岸灘傳播的數(shù)學(xué)模型方面的研究較少。所以，進(jìn)一步發(fā)展波浪在不平底床傳播的數(shù)學(xué)模型，對波浪在不平底床傳播的研究是有意義的。

本文的目的為通過引入繞射因子，進(jìn)一步擴(kuò)展 Liu 和 Shi[1]提出的模型，并對其進(jìn)行驗(yàn)證。

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

本文采用的是 Liu 和 Shi[1]中的數(shù)學(xué)模型，為了全文的完整性，本章將數(shù)學(xué)模型的推導(dǎo)過程在此做簡要介紹。本文數(shù)學(xué)模型的推導(dǎo)基于經(jīng)典的線性波理論。以笛卡爾坐標(biāo)系為參考坐標(biāo)，水平x軸位于靜水面，與岸線的交點(diǎn)為坐標(biāo)原點(diǎn)，指向外海。z軸垂直于靜水面，向上為正。靜水水深為（x，y）的函數(shù)，即 h=h（x，y）。假設(shè)外海入射波為單色波，無窮遠(yuǎn)處水深不變，比如 z=-h0，h0是常數(shù)。其速度勢函數(shù)如下：

式中：φ 為速度勢；g為重力加速度；a0為深海波幅；i為虛數(shù)單位，ω 為波浪頻率；k0為深海波數(shù)。

ζ為波面函數(shù)，代表自由水位高程，表達(dá)式為：

波浪的頻率 ω 與波數(shù) k0滿足色散關(guān)系式：

波浪在向近岸傳播過程中會因地形的變化發(fā)生變形效應(yīng)，引起折射和繞射，如果地形變化劇烈還可能發(fā)生反射。

假設(shè)海底地形為緩變海底，（|▽h|/（kh））為一小量ε，即在一個波長范圍內(nèi)，h的變化很小。引入兩種尺度，長尺度 （x?，y?）= ε （x，y） 和 短尺度（x，y），前者對應(yīng)于地形變化的尺度，后者對應(yīng)于波浪傳播的尺度。利用多重尺度法，可以得到：

z方向只有一種尺度。假設(shè)解具有以下形式：

其中：φ0為零階近似解；φ1為一階近似解。將這些關(guān)系代入通常的線性波浪問題的控制方程和邊界條件，可以得出：φ0滿足拉普拉斯方程和自由表面邊界條件及底部邊界條件：

φ1滿足：

無窮遠(yuǎn)處：

式 （6d） 說明，在無窮遠(yuǎn)處只存在零階近似解。此處的一階近似解是在線性波理論框架下考慮緩坡的高階近似。

從式 （6a） ~ （6d） 可解得一階近似解。對 x求偏導(dǎo)：

對 y，y?的推導(dǎo)和 x，x?類似。線性色散關(guān)系兩邊分別對 x?和y?求偏導(dǎo)，得到：

代入式（6a）~（6d），一階近似解的方程和條件變?yōu)椋?/p>

無窮遠(yuǎn)處，φ1=0。

其中：

以上推導(dǎo)中的參數(shù)均為關(guān)于水平尺度（x，y）的函數(shù)。假設(shè)一階解的形式為：

將式 （10） 代入式 （9a），比較等式兩邊的系數(shù)，得到：

需要指出的是，系數(shù) αi，（i=1，2，3，4，5）均為實(shí)數(shù)，零階近似解的波幅函數(shù)也在實(shí)數(shù)域內(nèi)，所以一階勢函數(shù) φ1= φa1meiS，幅值函數(shù) φa1m也是實(shí)數(shù)。

結(jié)合色散關(guān)系式，將式 （10） 代入式 （9b），自由表面邊界條件為：

將 αi（i=1，2，5）代入式 （12）：

將 A ，C1，C2的表達(dá)式代入式 （13），得到波幅方程：

這里的波幅函數(shù)是零階近似解的波幅。

2 計算方法

在每 個網(wǎng) 格點(diǎn) ，k（x?，y?，t）與當(dāng)?shù)?深度 h（x?，y?，t）的關(guān)系是如同 h 為常數(shù)一樣的色散關(guān)系。根據(jù)已知每個網(wǎng)格點(diǎn)的水深，及波浪在此傳播過程中周期不變，可求得每個網(wǎng)格點(diǎn)的波數(shù) k。波數(shù)無旋性條件還可寫成如下形式：

式 （15） 為一階的非線性偏微分方程，通過有限差分可以求得其數(shù)值解。

差分格式采用以下形式：

可用迭代法求解。w為權(quán)重因子，α為松弛因子，用以加速數(shù)值迭代之速度。當(dāng) w=1.0 時，方程為全隱格式。因?yàn)槭?（16） 中，x 方向?yàn)橄蚯安罘?，y方向?yàn)橹行牟罘?，引進(jìn)權(quán)重因子 w，考慮第 i行 y方向中心差分的影響。松弛因子 α的作用是考慮 j-1 列和 j+1 列的影響。波浪傳播方向的初值用 Snell定律求解，然后用迭代的方法求解式 （16），當(dāng)兩次求解的波向角相差小于某個小參數(shù)δ時，迭代停止。

下面為構(gòu)造方程式 （14） 的 Crank-Nicolson[10]格式離散方程：

k1（i，j）代表 k1在 （i，j）處的值，k2（i，j）代表 k2在 （i，j）處的值。式 （17） 在波數(shù) k 的離散和波幅的 y 方向做的 Crank-Nicolson[10]離散，x 方向是從外海邊界到海岸方向做步進(jìn)求解 （a0為波幅邊界條件），計算格式實(shí)現(xiàn)了穩(wěn)定，編程的工作量也相對減輕。

y方向透射邊界條件：

3 計算實(shí)例

荷蘭 Berkhoff[2]、Berkhoff等[3]的橢圓形淺灘試驗(yàn)是緩坡方程被引用最多的試驗(yàn)，試驗(yàn)地形和測量斷面如圖1所示。為了驗(yàn)證本文模型，根據(jù)上節(jié)中的方法，將本文模型應(yīng)用于上述算例。采用等間距網(wǎng)格，Δx 和 Δy 均取為 0.25m。

圖1 Berkhoff[2]、Berkhoff 等[3]橢圓形淺灘試驗(yàn)水深等值線 （單位：m）及 8個波幅測量斷面分布Fig.1 Contoursofwater depth（m）and distribution of eight sections for wave am p litudemeasurements of Berkhoff[2]、Berkhoff etal.[3]elliptic shoal test

根據(jù)上節(jié)中波浪傳播方向的求解方法，計算得到的越過 Berkhoff[2]、Berkhoff等[3]試驗(yàn)橢圓形淺灘的波向角計算結(jié)果如圖 2 所示。因?yàn)閷?Berkhoff[2]、Berkhoff等[3]淺灘試驗(yàn)而言，波浪傳播方向的變化主要發(fā)生在淺灘之后，淺灘之前和靠近側(cè)邊界的地方波向基本沒有變化，所以，在圖2中，只選取了部分區(qū)域，比如在圖 2 （a） 中，取 5m

圖2 本模型計算得到的波浪越過 Berkhoff[2]、Berkhoff等[3]試驗(yàn)中橢圓形淺灘的波向量分布Fig.2 Distributionsofwave vectorscalculated by the presentm odelover the elliptic shoalof Berkhoff[2]or Berkhoff etal.[3]’s test

因 為 波 浪 傳 播 方 向 很 難 測 量 ， Berkhoff[2]、Berkhoff 等[3]試驗(yàn)中 沒有提供波向角的試驗(yàn)數(shù)據(jù)。本文模型計算得到的波向角結(jié)果，以及其他模型的波向角結(jié)果，均不能通過試驗(yàn)數(shù)據(jù)驗(yàn)證精確性。鑒于以上原因，本文通過與其他模型的波向角計算結(jié)果比較，分析計算結(jié)果的合理性。Maa 等[11]針對應(yīng)用較為廣泛的幾種模型，如 REF/DIF-1、RCPWAVE 等，做了歸類和比較，給出了各種不同模型關(guān)于波向角的計算結(jié)果。總結(jié)各個模型波向角計算結(jié)果的特點(diǎn)，波向的改變主要發(fā)生在淺灘后，并且向淺灘后的中心線相聚，這與本文的結(jié)果圖 2 中 （b） 區(qū)域是符合的。

位于淺灘后的8個斷面的計算值和試驗(yàn)值的對比如圖 3。對于 Berkhoff[2]、Berkhoff等[3]的特定地形，本文計算值給出了一個較為合理的結(jié)果，波高分布與試驗(yàn)值較為相似。斷面SⅠ的計算值與試驗(yàn)值非常吻合；斷面SⅡ的計算值與試驗(yàn)值在淺灘兩側(cè)吻合的較好，在淺灘上計算值要高于試驗(yàn)值；斷面 SⅢ、斷面 SⅣ和斷面 SⅤ位于淺灘后，由于折射作用，會發(fā)生波能的相聚，波高增大，某些區(qū)域達(dá)到了入射波高的兩倍多，本文計算值給出了一個在波高和變化趨勢上非常接近的結(jié)果，但是沒有反映出兩側(cè)波高“減小—增大—減小”波動的過程。在向岸方向的3個斷面中，斷面SⅥ計算值與試驗(yàn)值吻合程度不夠理想，斷面 SⅦ和斷面 SⅧ計算值與試驗(yàn)值吻合較好。

圖3 本模型計算值與 Berkhoff[2]、Berkhoff 等[3]橢圓形淺灘試驗(yàn) 8 個斷面（SⅠ、SⅡ、SⅢ、SⅣ、SⅤ、SⅥ、SⅦ和 SⅧ） 的相對波高試驗(yàn)值比較Fig.3 Com parison between calculated relativewave heightsby the presentmodeland experimental data for eight sections（SⅠ，SⅡ，SⅢ，SⅣ，SⅤ，SⅥ，SⅦ，and SⅧ） of Berkhoff[2]、Berkhoff etal.[3]'selliptic shoal test

本文數(shù)值結(jié)果的不足及其原因：

1） 如圖 3 所示，在剛過淺灘的地方，斷面SⅡ與斷面 SⅢ之間，計算值略高于試驗(yàn)值，這可能是由于本文模型為線性模型，沒有考慮非線性作用。如果考慮非線性作用，能夠得到更好的結(jié)果。在本模型以后的研究中，可以引入非線性色散關(guān)系來考慮非線性效應(yīng)。

2）與試驗(yàn)數(shù)據(jù)相比，在淺灘的后面，即入射邊界 10m 之后，斷面 SⅣ和斷面 SⅤ缺少因繞射而引起的波高在淺灘兩側(cè)的波動，顯得平滑。原因是因?yàn)橐肓斯獬毯瘮?shù)方程。

本文的目的是對 Liu 和 Shi[1]提出的 數(shù)學(xué)模式進(jìn)行數(shù)值計算驗(yàn)證和擴(kuò)展，使其可以應(yīng)用于平面二維波浪變形問題，本文重點(diǎn)在于驗(yàn)證 Liu 和Shi[1]模型應(yīng)用與平面二維問題的有效性。模擬的波高結(jié)果在部分?jǐn)嗝嫔系挠嬎憔壬行杼岣?，今后，可以將計算結(jié)果與其它模型進(jìn)行定量對比。

4 結(jié)語

本文基于 Liu 和 Shi[1]的數(shù)學(xué)模型，引入繞射因子，利用有限差分法，對數(shù)學(xué)模型的波幅方程和邊界條件進(jìn)行數(shù)值離散，空間步進(jìn)求解波浪由外海向近岸傳播時，因?yàn)檎凵?、繞射和淺化等波浪變形現(xiàn)象，引起的波高、波向等波浪要素的變化。將本文的數(shù)值模 型應(yīng)用于經(jīng)典的 Berkhoff[2]、Berkhoff等[3]橢圓形淺灘算例，并與其試驗(yàn)結(jié)果做比較驗(yàn)證，證實(shí)本文模型求解二維地形下的波浪折射、繞射是可行的，具有方程求解簡單、可適用于大面積海域的優(yōu)點(diǎn)。并且，光程函數(shù)的引入，使 Liu 和 Shi[1]的模型應(yīng)用范圍更廣泛，在求解平面二維問題時精度得到提高。

[1]LIU Y Z，SHIJZ.A theoretical formulation for wave propagations overuneven bottoms[J].Ocean Engineering，2008，35（3-4）：426-432.

[2]BERKHOFF J C W.Computation of combined refractiondiffraction[C]//Proceedingsof the 13th international conference on coastalengineering.Vancouver，ASCE，1972：471-490.

[3]BERKHOFF J C W，BOOIJ N，RADDER A C.Verification of numerical wave propagation models for simple harmonic linear waterwaves[J].Coastal Engineering，1982（6）：255-279.

[4]ECKARTC.The propagation of gravitywaves from deep to shallow water[M]//Gravitywaves.Washington DC：National Bureau of StandardsCircular 521，1952：165-173.

[5]BIESEL F.Study ofwave progression in waterofgradually varying depth[M]//Gravity waves.Washington D C：National Bureau of Standards Circular521，1952：243-253.

[6]BATTJESJA.Refraction ofwaterwaves[J].Journalof theWaterways，Harbors and Coastal Engineering Division ASCE，1968，WW4：437-451.

[7]MEIC C，LEMEHAUTE B.Note on the equations of long waves over an uneven bottom[J].Journal of Geophysical Research，1968，71：393-400.

[8]MEIC C，TLAPA G A，EAGLESON PS.An asymptotic theory for water waves on beaches of mild slope[J].JournalofGeophysical Research，1968，73：4 555-4 560.

[9]LIU P L F，DINGEMANS M W.Derivation of the third-order evolution equations for weakly nonlinear water waves propagating over uneven bottoms[J].WaveMotion，1989，11：41-64.

[10]CRANK J，NICOLSON P.A practical method for numerical evaluation ofsolutions ofpartialdifferentialequationsof the heatconduction type[J].Proceedings of the Cambridge Philosophical Society，1947，43：50-67.

[11]MAA JP Y，HSU TW，TSAIC H，JUANG W J.Comparison of wave refraction and diffractionmodels[J].Journal ofCoastal Research，2000，16：1 073-1 082.

A two-dimensional horizontal numericalmodel of wave propagation over an uneven coastal bottom

PENGYan-jian，LIUYing-zhong， SHIJohn Z.
（SchoolofNavalArchitecture，Ocean and CivilEngineering，Shanghai Jiao Tong University，Shanghai200030，China）

Based on the linearwave theory，by using themultiple-scalemethod and a small slope parameter ε，Liu and Shi[1]obtain the first-order approximate solutions of the wave potential function，aswell aswave amplitude equation， for the problem of wave propagation over an uneven coastal bottom.This study introduces the diffraction factor，which expands Liu and Shi[1]＇s application range and improves itsaccuracy.Both wave amplitude equation and boundary conditions are discretized and solved by the finite differencemethod.Variations in wave amp litude and direction ofwave propagation caused by the refraction and diffraction are obtained.The validation study hasbeenmade useing of theexperimental results of Berkhoff[2]or Berkhoff et al[3]＇selliptic shoal test.It is suggested that the presentmodelbe applied to two-dimensionalwave propagation problems.

multiple-scalemethod；wave propagation；uneven bottom

TV139.2

A

1003-3688（2014）01-0008-06

10.7640/zggw js201401002

2013-04-12

2013-09-02

海洋工程國家重點(diǎn)實(shí)驗(yàn)室自主研究課題 （GKZD010056-9）

彭延建 （1982 — ），男，山東德州市人，碩士，工程師，港口

海岸及近海工程專業(yè)。

* 通訊作者：時鐘，E-mail：zshi@sjtu.edu.cn