張文海，賈會(huì)杰，王炎

（中交一航局第一工程有限公司，天津300456）

深水到淺水域非線性波傳播的數(shù)值模型

張文海，賈會(huì)杰，王炎

（中交一航局第一工程有限公司，天津300456）

統(tǒng)一方程是在Stokes波理論與Boussinesq型方程相結(jié)合的基礎(chǔ)上推導(dǎo)出的，適用于深水及淺水域波浪的傳播。文中首先分析了統(tǒng)一方程的頻散性及其適用性。其次，采用ADI法對(duì)控制方程進(jìn)行離散，并對(duì)控制方程中的非線性項(xiàng)進(jìn)行線性化近似處理，用改進(jìn)的Patankar半隱格式方法求解動(dòng)量方程。直接給定入射邊界條件，出流邊界條件采用Sommerfeld邊界條件和消波層相結(jié)合的方法，從而建立起從深水到淺水域都有效的數(shù)值模型。最后，利用平底與圓形暗礁組合地形上波浪傳播的經(jīng)典物理模型實(shí)驗(yàn)來驗(yàn)證數(shù)值模型的精確性。將實(shí)驗(yàn)結(jié)果與數(shù)值解相比較，兩者吻合較好，說明本文建立的數(shù)學(xué)模型能有效地模擬水深復(fù)雜變化水域波浪傳播，具有較高的適用性。

統(tǒng)一方程；數(shù)值模擬；消波；非線性波；物理模型實(shí)驗(yàn)

0 引言

在引起海岸變形和岸灘演變的眾多因素中，波浪為最重要的因素之一。在探索各種適用方法將原屬于三維波動(dòng)問題的豎向坐標(biāo)變量分離化為二維水平控制方程的過程中，產(chǎn)生了Boussinesq型方程（非線性長(zhǎng)波方程）和Stokes波型方程（緩坡方程）這兩種用于解決水波問題的控制方程[1]。

Boussinesq型方程屬于非線性波的范疇，能較好地描述非線性作用。但只能適用于淺水域。緩坡方程具有完全頻散性的特征，屬于線性波的范疇，適用于淺水域和深水域，但只適合底坡緩慢變化或者弱的非線性波。雖然國內(nèi)外眾多學(xué)者進(jìn)行了許多改進(jìn)以擴(kuò)大其適用范圍，但是相伴而來的是方程求解的復(fù)雜化。考慮到緩坡方程和Boussinesq型方程兩者的優(yōu)缺點(diǎn)，Li[2]推導(dǎo)了一個(gè)不增加額外導(dǎo)數(shù)項(xiàng)，適合深水到淺水域波浪傳播的統(tǒng)一方程。在有限水深情況下統(tǒng)一方程可化為改進(jìn)的Boussinesq方程。在深水情況下線性化的統(tǒng)一方程可化為經(jīng)典波動(dòng)方程。

1 水深任意變化水域非線性波傳播的數(shù)值模型

1.1 控制方程及其分析

文獻(xiàn)[2]中推導(dǎo)了不增加額外導(dǎo)數(shù)項(xiàng)，適用于深水到淺水域非線性波傳播的統(tǒng)一方程，假設(shè)x、y軸位于靜水面，z軸豎直向上，其控制方程為：

式中：η為水面高程；k為波數(shù)；g為重力加速度；h為水深；矢量u為水平速度（u，v）；uα為zα處的水平速度矢量，對(duì)于規(guī)則波傳播zα=βh= -0.66h[2]。相比于經(jīng)典Boussinesq型方程，統(tǒng)一方程式（1）和式（2）沒有增加額外的高階導(dǎo)數(shù)項(xiàng)，且非線性項(xiàng)與經(jīng)典Boussinesq型方程一致。

在淺水域情況下，利用雙曲函數(shù)的泰勒級(jí)數(shù)，并忽略高階導(dǎo)數(shù)項(xiàng)O（μ4），函數(shù)F1、F2、F3、zαF4可以表示為：

將式（7）、式（8）和式（10）代入到式（1）中，式（9）和式（10）代入到式（2）中，統(tǒng)一方程可以變?yōu)椋?/p>

由上述變換所得式（13）和式（14）和改進(jìn)的Boussinesq型方程一致（Nwogu 1993[3]；Liu 1995[4]）。因此，在淺水區(qū)域中如果忽略高階項(xiàng)O（μ4），統(tǒng)一方程可轉(zhuǎn)化為改進(jìn)的Boussinesq型方程。

緩坡方程（Berkhoff 1972[5]）具有完全頻散特性，適用于淺水域和深水域，易推廣于底坡緩慢變化和弱非線性，但由于是從線性波理論推導(dǎo)而來的，求解強(qiáng)非線性問題相比Boussinesq型方程較差。要求水底坡度較緩和邊界波浪的非線性是緩坡方程主要的局限性。而對(duì)于統(tǒng)一方程沒有這些限制，若以同樣的方法建立數(shù)值模型，在淺水域中統(tǒng)一方程可以考慮到非線性波的影響，而緩坡方程卻不能。

1.2 控制方程的差分格式

小學(xué)數(shù)學(xué)教學(xué)情境的創(chuàng)設(shè)，要符合不同年齡段兒童的心理特點(diǎn)和認(rèn)知規(guī)律，要根據(jù)不同的教學(xué)內(nèi)容有所變化。多年來，通過參加課題研究和教學(xué)實(shí)踐，從在新課引入時(shí)、在新知的探究中、在知識(shí)鞏固上、在整個(gè)課堂教學(xué)中創(chuàng)設(shè)不同的情境，我探索了一些方法，并取得了一些收獲：

離散控制方程采用有限差分法，式（1）和式（2）的數(shù)值求解是基于空間交錯(cuò)網(wǎng)格系統(tǒng)。對(duì)于一個(gè)時(shí)間步求解統(tǒng)一方程的步驟如下：

1）采用雙向掃描法求解動(dòng)量方程；

2）返回1）進(jìn)行迭代；

3）用交替方向隱格式法（ADI）求解連續(xù)性方程。

控制方程離散格式詳見文獻(xiàn)[2]。

1.3 波浪傳播的邊界條件

1.3.1 入射邊界條件

本文采用直接給出入射波波面和速度分布作為入射條件，對(duì)于不同的地形條件，入射波也不盡相同，本文采用一階Stokes波理論作為入射邊界條件，入射波均為正向入射。其表達(dá)式為：

用ηc*代替ηc，用uc*代替uc，再進(jìn)行下一步計(jì)算。

1.3.3 出流邊界條件

在計(jì)算域的下游邊界采用Sommerfeld[8]（1949）邊界條件。本文數(shù)值模擬中應(yīng)用于外向波沿向岸邊的邊界條件可以寫為：

式中：H為入射波波高；h為水深；k為波數(shù)；0.34kh是由zα=-0.66h[2]得出的。

1.3.2 入射邊界上反射波的消波

直接給定波面和特征流速確定入射邊界條件的方法簡(jiǎn)便直觀，但計(jì)算域內(nèi)存在反射或偽反射波時(shí)，會(huì)導(dǎo)致計(jì)算程序的不穩(wěn)定。因此本文采用張洪生等[6]提出的吸收入射邊界上反射波的方法。在離入射邊界一定長(zhǎng)度區(qū)域內(nèi)設(shè)置海綿層，在每一時(shí)間步后，變量η和u都要除以指數(shù)函數(shù)β（x），作者在文獻(xiàn)[7]中進(jìn)行了進(jìn)一步深化及論證，其表達(dá)式為：

式中：Q表示uα、vα和η；c為波速；γ為邊界法線方向和出射波的方向之間的夾角。通過假設(shè)γ=0，式（20）有不同的有限差分格式。本文采用的方法為：

式中：IB代表邊界網(wǎng)格點(diǎn)。

2 均勻水深水域波浪傳播的數(shù)值模擬

為了驗(yàn)證統(tǒng)一方程的適用性，就均勻水深水域中波浪的傳播進(jìn)行了系統(tǒng)的數(shù)值模擬。設(shè)置長(zhǎng)度為85 m的水槽，在水槽左端采用二階Stokes波作為入射條件，右端為自由出流的開邊界。為了進(jìn)行系統(tǒng)的數(shù)值模擬，設(shè)置固定波周期T=2.0 s，分別模擬了kh=1、kH=0.2，kh=2、kH=0.3，kh=3、kH=0.4，kh=100、kH=0.5四種組合情況。在數(shù)值計(jì)算過程中，空間步長(zhǎng)取為Δx= 0.05 m，時(shí)間步長(zhǎng)取為Δt=0.025 s，數(shù)值模擬結(jié)果如圖1所示。

圖1計(jì)算值與二階Stokes波理論值的比較Fig.1Comparison of numerical results with Stokes second-order wave

圖1 為各計(jì)算組合情況下相對(duì)高程沿程變化的過程線，可以看出在不同水深、不同參數(shù)的情況下數(shù)值解與理論解都吻合良好。表明均勻水深情況下模型對(duì)深水域波在較長(zhǎng)的時(shí)間和傳播距離內(nèi)可進(jìn)行穩(wěn)定準(zhǔn)確地?cái)?shù)值模擬。因此統(tǒng)一方程可適用于模擬淺水到深水幾乎整個(gè)水域弱非線性波的傳播。

3 平底與圓形暗礁組合地形上波浪傳播的數(shù)值模擬

Ito 1972[9]對(duì)平底與圓形暗礁組合的地形上的波浪傳播進(jìn)行了物模實(shí)驗(yàn)，這一物理模型試驗(yàn)是檢驗(yàn)數(shù)學(xué)模型精度的經(jīng)典試驗(yàn)之一，試驗(yàn)地形見圖2。

圖2 Ito物理模型試驗(yàn)示意圖（單位：m）Fig.2Sketch of physical experiment by Ito（m）

入射波周期T為6.3 s，波高為1.0 m，上邊界為自由出流的開邊界，側(cè)邊界為全反射的固壁邊界。選取空間步長(zhǎng)為L(zhǎng)0/14，時(shí)間步長(zhǎng)為T/32，數(shù)值模擬結(jié)果如圖3～圖6所示。

圖3 相對(duì)波高的計(jì)算結(jié)果與實(shí)驗(yàn)數(shù)據(jù)比較圖（a斷面）Fig.3Comparison of computational results with experimentaldata of relative wave height(Section a)

圖4 相對(duì)波高的計(jì)算結(jié)果與實(shí)驗(yàn)數(shù)據(jù)比較圖（b斷面）Fig.4Comparison of computational results with experimental data of relative wave height(Section b)

圖5 不同測(cè)點(diǎn)處波面隨時(shí)間變化的數(shù)值解Fig.5Numerical Simulation of wave surface changing with time at various positions

圖6數(shù)值模擬的波高立體圖Fig.6Wave height stereogram of numerical simulation

圖3 、圖4分別是a、b斷面處的相對(duì)波高，由圖可見數(shù)值解與實(shí)驗(yàn)值吻合良好。圖5為A、B、C、D四個(gè)點(diǎn)的波面隨時(shí)間變化過程曲線圖，由圖可以看出波動(dòng)過程隨時(shí)間變化十分穩(wěn)定，說明計(jì)算結(jié)果可靠。深水點(diǎn)A處波形呈余弦波，波谷和波峰對(duì)稱。由于橫斷面上中心凸起，在由深水到淺水（B向C）的傳播過程中，波高沿程逐漸減小，出了圓形暗礁后的平坦水域，水槽中心線的波高沿程逐漸增加。圖6為數(shù)值模擬的波高立體圖。通過數(shù)值計(jì)算結(jié)果與物模試驗(yàn)值的比較表明，此數(shù)值模型能夠有效地模擬復(fù)雜變化地形上的波浪傳播。

4 結(jié)語

相比于Boussinesq型方程和緩坡方程這兩種基本方程，統(tǒng)一方程顯著的改進(jìn)是其能準(zhǔn)確的滿足規(guī)則波在深水和淺水域的頻散關(guān)系，其淺化梯度系數(shù)與適合深水域的Stokes波理論值吻合良好。若忽略高階項(xiàng)O（μ4），在淺水情況下統(tǒng)一方程可化為改進(jìn)的Boussinesq型方程。本文以統(tǒng)一方程為控制方程，基于空間交錯(cuò)網(wǎng)格系統(tǒng)對(duì)其進(jìn)行數(shù)值求解，采用交替方向隱格式法（ADI）求解連續(xù)性方程，用改進(jìn)的Patankar（1980）半隱格式方法求解動(dòng)量方程，從而建立了新的數(shù)值模型。在給定初始條件和入射邊界條件后，采用Sommerfeld型邊界條件和消波層相結(jié)合的方法處理出流邊界。采用兩者相結(jié)合的邊界條件，能夠有效地吸收波浪的能量，消除或者減少波浪能量的反射，從而保證數(shù)值模擬的精度。利用Ito對(duì)平底與圓形暗礁組合的地形上波浪傳播的物模實(shí)驗(yàn)對(duì)數(shù)值模型的精度及適用性進(jìn)行驗(yàn)證。兩者吻合良好，表明本文所建立的數(shù)值模型可以較為有效地模擬水深復(fù)雜變化地形上非線性波的傳播與變形，具有較高的適用性。本模型未考慮波浪水流等對(duì)波浪變形的影響，這有待于進(jìn)一步的探討和改善，其對(duì)近岸工程、遠(yuǎn)海生產(chǎn)作業(yè)、近岸環(huán)境保護(hù)等方面具有十分重要的意義。

[1]張洪生.非線性波傳播的數(shù)值模擬[D].南京：河海大學(xué)，2000.

ZHANG Hong-sheng.Numerical modeling of nonlinear wave propagation[D].Nanjing:Hohai University,2000.

[2]LI B.Wave equations for regular and irregular water wave propagation[J].Journal of Waterway,Port,Coastal and Ocean Engineering, 2008,134(2):121-142.

[3]NWOGU O.Alternative form of Boussinesq equations for nearshore wave propagation[J].Journal of Waterway,Port,Coastal and Ocean Engineering,1993,119(6):618-638.

[4]LIU L F,CHO Y S,BRIGGS M J,et al.Runup of solitary waves on a circular island[J].Journal of Fluid Mechanics,1995,302:259-285.

[5]BERKHOFF J C W.Computation of combined refraction-diffraction [C]//Proceedings of 13th conference on coastal engineering,ASCE, 1972:471-490.

[6]張洪生，商輝.對(duì)波浪入射邊界上反射波的消波及其驗(yàn)證[J].上海交通大學(xué)學(xué)報(bào)，2008，42（4）：674-678.

ZHANG Hong-sheng,SHANG Hui.A method of absorbing the reflectedwavesonincidentboundaryand its verification[J].Journal of Shanghai Jiaotong University,2008,42(4):674-678.

[7]ZHANG Hong-sheng,WANG Yan,XU Chun-hui,et al.Tests and applications of an approach to absorbing reflected waves towards incident boundary[J].China Ocean Engineering,2013,27(6):703-718.

[8]SOMMERFELDA.Partialdifferentialequationsinphysics: lectures on theoretical physics,Vol.6[M].Salt Lake City,UT: Academic Press,1949.

[9]ITO Y,TANIMOTO K.A method numerical analysis of wave propagation application to wave diffraction and refraction[C]/Proceedings of 13th conference on coastal engineering,ASCE,1972: 503-522.

Numerical modeling of nonlinear wave propagation from deep water to shallow water

ZHANG Wen-hai,JIA Hui-jie,WANG Yan
（No.1 Engineering Co.,Ltd.of CCCC First Harbor Engineering Co.,Ltd.,Tianjin 300456,China）

The unified equations are derived from the Stokes second-order wave theory and the Boussinesq-type equations.It is suitable for the propagation of waves in deep and shallow seas.We firstly analyzed the dispersion and applicability of the unified equations,then used the ADI method to disperse the governing equations,processed the nonlinear terms of the governing equations by linear approximation,and used the modified Patankar with semi-implicit schemes to solve the momentum equations.Given the boundary conditions,the outflow boundary conditions are combined with the Sommerfeld boundary condition and the wave elimination layer,so as to establish a valid numerical model suitable for wave transformation from deep water to shallow water.At last,the experiment data from physical model of wave propagation and deformation in the complicated water is used to verify the accuracy of present numerical model.The experimental results are in good agreement with those of numerical solution.It indicates that the numerical model can effectively simulate the wave propagation in the water with varying topography,and has a high applicability.

unified equations;numerical simulation;absorbing waves;nonlinear wave;physical experiment

U651.3；P731.2

A

2095-7874（2017）06-0036-05

10.7640/zggwjs201706008

2017-01-25

2017-04-14

張文海（1976—），男，天津市人，高級(jí)工程師，副總經(jīng)理，港航工程專業(yè)。E-mail：33162727@qq.com