房克照,鄒志利

(大連理工大學(xué)海岸和近海工程國(guó)家重點(diǎn)實(shí)驗(yàn)室,遼寧 大連 116024)

基于四階完全非線性Boussinesq水波方程二維波浪傳播數(shù)值模型

房克照,鄒志利

(大連理工大學(xué)海岸和近海工程國(guó)家重點(diǎn)實(shí)驗(yàn)室,遼寧 大連 116024)

建立基于四階完全非線性Boussinesq水波方程的二維波浪傳播數(shù)值模型。采用 Kennedy等提出的渦粘方法模擬波浪破碎。在矩形網(wǎng)格上對(duì)控制方程進(jìn)行離散,采用高精度的數(shù)值格式對(duì)離散方程進(jìn)行數(shù)值求解。對(duì)規(guī)則波在具有三維特征地形上的傳播過(guò)程進(jìn)行了數(shù)值模擬,通過(guò)數(shù)值模擬結(jié)果與實(shí)驗(yàn)結(jié)果的對(duì)比,對(duì)所建立的波浪傳播模型進(jìn)行了驗(yàn)證。同時(shí),為了考察非線性對(duì)波浪傳播的影響,給出和上述模型具有同階色散性、變淺作用性能但僅具有二階完全非線性特征的波浪模型的數(shù)值結(jié)果。通過(guò)對(duì)比兩個(gè)模型的數(shù)值結(jié)果以及實(shí)驗(yàn)數(shù)據(jù),討論非線性在波浪傳播過(guò)程中的作用。研究結(jié)果表明,所建立的Boussinesq水波方程在深水范圍內(nèi)不但具有較精確的色散性和變淺作用性能,而且具有四階完全非線性特征,適合模擬波浪在近岸水域的非線性運(yùn)動(dòng)。

Boussinesq水波方程;波浪傳播模型;波浪破碎;渦粘

Abstract:The present study develops a 2Dwave breakingmodel based on the fourth-order full nonlinear Boussinesq equations.This setof equations possesses higher accuracy in dispersion,shoaling and full nonlinearity characteristics,which enable it applicable in describingwater waves in nearshore zone.The discretizationis conducted on a rectangle grid system,higherorder finite difference formulasare used to approximate the spatial derivatives and high accuracy scheme is used to integrate equations in time.Kennedy et al′s eddy viscosity breakingmechanism is incorporated to mimicwave breaking.Themathematicalmodel is run to simulate waves transformation over typical three dimensional bathymetries.The numerical results are verified against available experimental data.The satisfactory agreement between numerical resultsand experimental data demonstrate the validation and efficiency of the presentwave breakingmodel.The effectof full nonlinearityonwave propagation is discussed by comparing the numerical resultswith those from Boussinesqmodelwith only the second-order full nonlinearity characteristics.

Key words:Boussinesqwave equations;wave propagationmodel;wave breaking;eddy viscosity

近幾十年來(lái),基于Boussinesq水波方程的波浪模型得到了長(zhǎng)足的發(fā)展并被廣泛用于模擬近岸區(qū)域各種波浪運(yùn)動(dòng)現(xiàn)象,見(jiàn)綜述性文獻(xiàn)[1]和[2]。建立這類(lèi)波浪模型時(shí),Boussiensq方程的性能是首要考慮的因素,方程的性能包括線性性能(色散性和變淺作用性能)和非線性性能(各階非線性傳遞函數(shù)以及波幅離散性能)。Boussinesq水波方程的發(fā)展過(guò)程很大程度上是提高上述各種性能的過(guò)程,起初諸多研究著重于提高方程的線性性能,代表性的工作有Madsen和Sorensen[3]、Nwogu[4]等,而后來(lái)的研究表明非線性在Boussinesq模型中有非常重要的作用。目前,具有強(qiáng)非線性特征的波浪模型不斷涌現(xiàn),其中Madsen等[5]給出的方程以及Lynett和Liu[6]提出的多層模型頗具代表性。前者最大應(yīng)用水深達(dá)到kh=40,但方程的數(shù)值求解過(guò)程涉及到大型稀疏矩陣的求解[7-8],即便對(duì)于一維問(wèn)題控制方程多達(dá)6個(gè),形成的稀疏矩陣帶寬至少為25[9];后者通過(guò)適當(dāng)增加模型層數(shù),也可以得到性能較好的方程,但隨著層數(shù)增加,控制方程增多,加大了數(shù)值求解工作的難度。相對(duì)而言,具有較低色散性精度但具有強(qiáng)非線性特征的方程數(shù)值求解簡(jiǎn)便,仍不失為進(jìn)行波浪數(shù)值模擬的有力工具。Wei等[10]給出了具有二階完全非線性特征的Boussinesq方程,基于此方程發(fā)展的二維模型FUNWAVE[11]得到了廣泛的應(yīng)用。Madsen和Schaffer[12]給出了具有四階非線性特征的Boussinesq方程,其非線性達(dá)O(εμ4)階(這里ε為表征方程非線性的參數(shù),定義為ε=A/h,A和h分別為特征波幅、特征水深;μ=kh為色散性參數(shù),k為特征波數(shù));Gobbi和Kirby[13]給出了具有四階完全非線性特征的Boussinesq方程FN4;Zou和Fang[14]也推導(dǎo)了具有四階完全非線性的Boussinesq方程BouN4D4。上述模型的相關(guān)研究成果均表明,具有強(qiáng)非線性特征的方程較僅具有弱非線性特征的模型有優(yōu)越性,更適合模擬波浪的非線性運(yùn)動(dòng)。因此,建立具有強(qiáng)非線性特征的波浪模型非常必要,然而上述模型中,FUNWAVE雖然具有二階完全非線性特征,但方程色散性?xún)H為精確色散關(guān)系的Pade[2,2]階近似,限制了模型在深水范圍的應(yīng)用;Madsen和Schaffer[12]給出的方程僅具有弱非線性特征同時(shí)并未建立相應(yīng)的數(shù)值模型;上述兩個(gè)具有完全非線性特征的方程(FN4和BouN4D4)雖然具有強(qiáng)非線性特征,但所建立的模型僅限于一維問(wèn)題,也并未考慮波浪的破碎問(wèn)題,而具有四階完全非線性特征的二維波浪模型更加鮮見(jiàn)。建立基于BouN4D4的二維波浪模型,該方程除具有四階完全非線性特征外,其色散性為精確色散關(guān)系的Pade[4,4]階近似,適用于深水波浪問(wèn)題,同時(shí)具有良好的變淺作用性能。在矩形網(wǎng)格上對(duì)控制方程進(jìn)行離散,建立高精度數(shù)值格式,考慮波浪的破碎問(wèn)題。針對(duì)幾個(gè)典型三維地形上波浪進(jìn)行數(shù)值模擬,通過(guò)對(duì)比分析數(shù)值結(jié)果和實(shí)驗(yàn)數(shù)據(jù),考察非線性特征在波浪傳播過(guò)程中的作用,同時(shí)驗(yàn)證具有四階完全非線性特征波浪傳播模型的有效性。

1 數(shù)值模型

1.1 Boussinesq方程

Zou和Fang[14]推導(dǎo)了具有四階完全非線性特征的Boussinesq水波方程,該組方程以波面η和水深平均速度ˉu表達(dá),經(jīng)擴(kuò)展后(考慮波浪破碎、水低摩擦等效應(yīng))方程的二維無(wú)因次形式為

式中:ˉu為水深平均速度,η為波面升高,h為靜水深,d=h+εη為當(dāng)?shù)厮?g為重力加速度,▽=(?x,?y)為二維偏微分算子 ,α1=1/9,β1=1/945,α2=0.146 488,β2=0.001 996 為常數(shù) ,其中α1、β1確定了方程的色散性,α2、β2決定了方程的變淺作用性能。上式中R定義為R=Rf+Rb+Rs+Rm,其中Rf代表底摩擦阻力,Rb、Rs、Rm分別為波浪破碎引起的耗散項(xiàng)、海綿吸收層和側(cè)混項(xiàng),其表達(dá)式同F(xiàn)UNWAVE中一致,fs為內(nèi)域波浪生成時(shí)引入的源函數(shù),同Gobbi和 Kirby[13]給出的一致。

后文中將上述模型稱(chēng)為BouN4D4,N4代表4階完全非線性,D4代表色散關(guān)系為精確色散關(guān)系的Pade[4,4]階近似。該方程保留了精確到O(μ4)階的所有非線性項(xiàng),最高階達(dá)O(ε5μ4),使得方程具有四階完全非線性特征,更適合描述具有非線性特征的波浪運(yùn)動(dòng)現(xiàn)象,方程的推導(dǎo)以及性能分析過(guò)程見(jiàn)文獻(xiàn)[14]。忽略上述方程中所有O(μ4)階非線性,則得到僅具有二階完全非線性特征的方程,稱(chēng)之為BouN2D4(N2代表2階完全非線性,D4代表色散關(guān)系為精確色散關(guān)系的Pade[4,4]階近似),其將用于和BouN4D4進(jìn)行對(duì)比,以考察非線性在波浪傳播過(guò)程中的作用。

1.2 數(shù)值求解

有限差分方法由于具有直觀性和容易程序?qū)崿F(xiàn)的優(yōu)點(diǎn),被廣泛用來(lái)數(shù)值求解Boussinesq方程,發(fā)展得比較完善,如FUNWAVE[11]以及 Gobbi和 Kirby[13]。因此,這里仍然采用有限差分方法對(duì)上述方程進(jìn)行數(shù)值求解,但值得注意的是若采用傳統(tǒng)的手工編程方法對(duì)方程(1)～(7))進(jìn)行程序?qū)崿F(xiàn)非常繁瑣,而通過(guò)充分利用MAPLE強(qiáng)大的數(shù)學(xué)公式推導(dǎo)功能及其向Fortran語(yǔ)言的轉(zhuǎn)化功能,可實(shí)現(xiàn)快速、準(zhǔn)確地撰寫(xiě)程序代碼的目的。文中采用的數(shù)值求解格式同上述文獻(xiàn)基本類(lèi)似,這里不再一一詳述,詳細(xì)過(guò)程可參考上述相關(guān)文獻(xiàn),下面僅給出主要步驟。

1.2.1 方程整理

將方程中含時(shí)間導(dǎo)數(shù)的線性項(xiàng)移在方程左側(cè),其它項(xiàng)移至方程右側(cè):

注意,上式中ξ,U,V中包含了所有含ηt,,的線性項(xiàng)。E,Et,F,Ft以及G,Gt分別為質(zhì)量方程,動(dòng)量方程中剩余的所有空間項(xiàng)和時(shí)間項(xiàng),r項(xiàng)為式(2)中R在x,y方向的分量。

1.2.2 方程離散及差分公式

在圖1所示的矩形網(wǎng)格上對(duì)控制方程進(jìn)行空間離散,x,y方向網(wǎng)格標(biāo)記分別為i(i=1,mi)和j(j=1,mj)。對(duì)于方程右端空間導(dǎo)數(shù)采用七點(diǎn)差分格式,內(nèi)點(diǎn)采用中心格式,邊界點(diǎn)采用偏心格式,對(duì)于單方向的差分統(tǒng)一用如下公式表示

式中:δfl表示函數(shù)在f網(wǎng)格點(diǎn)l處的空間導(dǎo)數(shù)(i-3≤l≤i+3)。當(dāng)l=i時(shí)為中心差分,當(dāng)l≠i時(shí)為偏心差分,γ,αm為差分格式的系數(shù),其具體取值以及各階差分格式精度參見(jiàn)附錄A.1。對(duì)于水平二維問(wèn)題,還涉及到交叉導(dǎo)數(shù)項(xiàng),其具體實(shí)施辦法在附錄A.2中列出?？梢钥闯?對(duì)于任意一個(gè)內(nèi)點(diǎn)(i,j),其差分近似涉及到周?chē)?5個(gè)點(diǎn),見(jiàn)圖1。對(duì)于方程左端線性時(shí)間項(xiàng)的系數(shù)采用五點(diǎn)中心差分格式,則按上述網(wǎng)格離散后,形成一個(gè)系數(shù)矩陣帶寬為五的線性方程組。對(duì)于方程右端對(duì)時(shí)間的一階導(dǎo)數(shù)項(xiàng),采用 Gobbi和 Kirby[13]文中公式進(jìn)行計(jì)算,這里不再給出。

圖1 空間離散示意以及差分模板Fig.1 The sketch of spatial discretization and difference stencil

1.2.3 時(shí)間積分、邊界處理等

按照上述網(wǎng)格以及差分公式對(duì)方程離散后,采用五階Adams-Bashforth預(yù)報(bào),六階Adams-Moulton校正格式進(jìn)行時(shí)間步進(jìn)求解。對(duì)于波浪入射邊界、固壁邊界的處理同用 Gobbi和 Kirby[13]文中完全一致(注意文獻(xiàn)中為一維),不再給出。波浪傳播過(guò)程中,由于非線性或者地形突變?nèi)菀滓鸶哳l數(shù)值震蕩,計(jì)算中采用SG光滑技術(shù)[7]對(duì)數(shù)值結(jié)果進(jìn)行光滑處理。

2 模型驗(yàn)證和對(duì)比

將對(duì)上面建立的數(shù)值模型進(jìn)行驗(yàn)證,針對(duì)幾個(gè)典型三維地形進(jìn)行波浪傳播的數(shù)值模擬。通過(guò)對(duì)比數(shù)值模擬結(jié)果和實(shí)驗(yàn)結(jié)果對(duì)模型進(jìn)行驗(yàn)證;對(duì)比模型BouN4D4和BouN2D4的數(shù)值結(jié)果,討論非線性作用在波浪傳播過(guò)程中的作用。

2.1 圓型淺灘上波浪的傳播

Chawla和Kirby[15]進(jìn)行了波浪跨通過(guò)圓型淺灘傳播的實(shí)驗(yàn)(見(jiàn)圖2)。該實(shí)驗(yàn)在水池中進(jìn)行(長(zhǎng)18 m,寬18.2m),水池中央放置圓型淺灘,其半徑為2.57m,淺灘中心位于(x=8.0 m,y=8.98 m)處,在七個(gè)不同斷面上進(jìn)行波高測(cè)量。圖2為計(jì)算區(qū)域以及測(cè)量斷面示意圖,取計(jì)算域長(zhǎng)23m,寬18.2 m,在計(jì)算域兩端分別設(shè)置2m、3m的海綿層,兩側(cè)面為固體邊界,空間網(wǎng)格步長(zhǎng)為0.10 m,時(shí)間步長(zhǎng)0.02 s,每0.1個(gè)周期對(duì)波面和速度應(yīng)用一次SG光滑器。這里僅選取其中兩組規(guī)則波實(shí)驗(yàn)進(jìn)行模擬。其中,第一組波浪沒(méi)有發(fā)生破碎,其對(duì)應(yīng)的波浪要素為水深0.45 m(圓型淺灘頂部水深0.08m),波高0.011 8m,周期1 s;第二組存在波浪破碎現(xiàn)象,對(duì)應(yīng)的波浪要素為水深0.395m,波高0.02m,周期1 s。該實(shí)驗(yàn)中除淺灘外均為常水深,這將導(dǎo)致波浪在淺灘后產(chǎn)生較強(qiáng)的波浪幅聚現(xiàn)象。

對(duì)波浪不破碎情況,模型運(yùn)行40 s,取最后6 s的數(shù)據(jù)進(jìn)行波幅計(jì)算(30 s以后波浪場(chǎng)基本達(dá)到穩(wěn)定狀態(tài))并將結(jié)果在圖3中給出。可以看出,BouN4D4和BouN2D4的數(shù)值結(jié)果基本重合,同時(shí)和實(shí)驗(yàn)數(shù)據(jù)吻合良好,這表明這一算例對(duì)高階非線性的要求較低。沿?cái)嗝鍭-A,模型較準(zhǔn)確地模擬了由于波浪淺化引起的波高增大現(xiàn)象以及淺灘后波高減小現(xiàn)象。還可以看出,淺灘處波高幾乎高達(dá)入射波高的2.7倍,這進(jìn)一步表明該處有強(qiáng)烈的波浪匯聚現(xiàn)象。對(duì)于其它的測(cè)量斷面,模型也給出了較好的模擬結(jié)果,這表明其較準(zhǔn)確地模擬了波浪的折射、繞射過(guò)程。從圖中還可以看出,由于淺灘的位置并未處于計(jì)算域中心位置,計(jì)算結(jié)果和實(shí)驗(yàn)數(shù)據(jù)關(guān)于淺灘中心軸線并不完全對(duì)稱(chēng)。

對(duì)于波浪破碎情況,模型運(yùn)行50 s,取最后10 s數(shù)據(jù)進(jìn)行均方根波高計(jì)算,圖4給出了數(shù)值結(jié)果同實(shí)驗(yàn)數(shù)據(jù)的對(duì)比(圖中H0為入射波高)?？梢钥闯?BouN4D4和BouN2D4給出的數(shù)值結(jié)果相近,但較波浪不破碎情況差異加大,這可能是由于入射波浪非線性增強(qiáng)所致。在斷面E-E之前,BouN4D4在y=10m區(qū)域附近給出的幅值略高于BouN2D4并且更接近實(shí)驗(yàn)數(shù)據(jù);在淺灘頂部(斷面F-F),兩模型模擬的波高均高于實(shí)驗(yàn)值但BouN4D4較BouN2D4和實(shí)驗(yàn)數(shù)據(jù)吻合程度稍好?？傮w而言,兩個(gè)模型的數(shù)值表現(xiàn)良好,但較不破碎情況差。這可能是由于以下幾個(gè)原因:1)用于模擬波浪破碎的渦粘公式是經(jīng)驗(yàn)性的,不能很好地描述波浪破碎這一復(fù)雜物理現(xiàn)象,同時(shí)由于實(shí)驗(yàn)中波浪破碎時(shí)具有三維特征,而采用的渦粘方法在數(shù)值實(shí)現(xiàn)時(shí)僅能沿某幾個(gè)固定方向?qū)Σɡ似扑檫M(jìn)行捕捉,存在誤差;2)波浪破碎時(shí)波面會(huì)變得非常陡峭,因此為了更好地刻畫(huà)波面,計(jì)算時(shí)所用的空間、時(shí)間步長(zhǎng)要比較小,但出于計(jì)算時(shí)間和對(duì)比考慮,這里參數(shù)的取值同第一個(gè)算例一致。

圖5給出了針對(duì)這兩個(gè)算例t=35 s時(shí)的波面瞬時(shí)圖。可以看出,由于波浪的淺化、折射、繞射以及非線性的相互作用,最終的流場(chǎng)呈現(xiàn)復(fù)雜的三維形態(tài)。其中,波浪不破碎和破碎兩種情況存在以下共同特征:波峰在淺灘上發(fā)生了彎曲,這是由于波浪折射所致,同時(shí)由于波浪淺化作用波高增大;在淺灘末端出現(xiàn)了明顯的波浪匯聚現(xiàn)象;由于淺灘的存在其后部都存在波浪掩護(hù)區(qū);觀察波峰線,都明顯存在波浪繞射現(xiàn)象。不同的是:波浪破碎情況下,淺灘上波面極其陡峭并且最大值出現(xiàn)在淺灘上,而波浪不破碎時(shí)淺灘上波面變化較緩且最大值出現(xiàn)在淺灘末端。

圖2 Chawla和kirby實(shí)驗(yàn)布置以及測(cè)量斷面分布Fig.2 Experiment setup and measurement transects(Chawla and Kirby)

圖3 各測(cè)量斷面波高分布Fig.3 Wave height comparisons along different transects

圖4 各測(cè)量斷面波高分布Fig.4 Wave height comparisons along different transects

圖5 穩(wěn)態(tài)后流場(chǎng)波高分布示意(BouN4D4)Fig.5 Elevation contoursof the steadywave fields(BouN4D4)

圖6給出了流場(chǎng)的時(shí)均流速圖。可以看到,在潛堤頂端,由于波浪破碎引起強(qiáng)烈的水流,幅值明顯大于其它區(qū)域波浪運(yùn)動(dòng)的水流速度。這同一般的波浪破碎理論相吻合,波浪破碎總是伴隨著流的產(chǎn)生,實(shí)驗(yàn)過(guò)程中也觀測(cè)到了這一現(xiàn)象[15]。

2.2 凸型淺灘上波浪的傳播

Whalin[16]在波浪水槽中進(jìn)行了一系列波浪在淺灘上的非線性折射和繞射物理實(shí)驗(yàn),實(shí)驗(yàn)水槽尺度為25.603m×6.096m。實(shí)驗(yàn)地形由以下公式給出式中:所有變量單位為米。圖7為計(jì)算域以及地形示意圖。在圖示地形上,波浪會(huì)發(fā)生變淺幅聚作用,由于非線性作用高次諧波變得尤為突出。因此,該實(shí)驗(yàn)被廣泛地用來(lái)驗(yàn)證各類(lèi)波浪模型[7-8,18]。實(shí)驗(yàn)可按照波浪周期不同分為三組,各組波浪要素及計(jì)算參數(shù)取值見(jiàn)表1,從表中可以看出,對(duì)應(yīng)于某一周期,入射波浪的非線性依次增強(qiáng),因此通過(guò)對(duì)這些工況的模擬,可以考察模型非線性精度對(duì)數(shù)值模擬結(jié)果的影響。下面,將模型BouN4D4和BouN2D4分別用于模擬Whalin的實(shí)驗(yàn)并進(jìn)行對(duì)比分析。

圖6 流場(chǎng)時(shí)均流速分布Fig.6 Phase-averaged velocity field

圖7 計(jì)算域及地形示意Fig.7 The computation domain and bathymetry sketch

表1 Whalin實(shí)驗(yàn)波浪要素及數(shù)值模擬網(wǎng)格尺寸Tab.1 Wave conditions for Whalin experiment and grid size in simulations

圖8給出了T=1 s時(shí),模型BouN4D4和BouN2D4的計(jì)算結(jié)果同實(shí)驗(yàn)數(shù)據(jù)的對(duì)比?？梢钥闯?數(shù)值結(jié)果同實(shí)驗(yàn)數(shù)據(jù)吻合良好。其中,當(dāng)H=0.019 4 m時(shí),兩個(gè)模型的數(shù)值結(jié)果幾乎沒(méi)有差別,都較好地模擬了各次諧波沿空間的變化,但隨著波浪非線性的增強(qiáng),兩個(gè)模型的差異變大,即當(dāng)H=0.039 m時(shí),BouN2D4給出的數(shù)值結(jié)果在x大于15 m時(shí)高于BouN4D4的數(shù)值結(jié)果,而B(niǎo)ouN4D4和實(shí)驗(yàn)數(shù)據(jù)更為接近。由于兩模型的色散性相同因此這一差別反映了非線性對(duì)數(shù)值結(jié)果的影響。圖9給出了T=2 s時(shí),模型BouN4D4和BouN2D4的計(jì)算結(jié)果同實(shí)驗(yàn)數(shù)據(jù)的對(duì)比。可以看出,隨著波浪非線性增強(qiáng)模型表現(xiàn)出的差異增大,這與上一算例一致,不同的是,兩個(gè)模型的數(shù)值結(jié)果差別較大。具體的,BouN2D4模擬的一、二次諧波幅值均高于實(shí)驗(yàn)數(shù)據(jù)和BouN4D4的結(jié)果而后者同實(shí)驗(yàn)數(shù)據(jù)更為吻合,其中一次諧波差別最大。對(duì)于三次諧波,BouN2D4給出的結(jié)果同實(shí)驗(yàn)數(shù)據(jù)吻合良好而B(niǎo)ouN4D4模擬的幅值略低于實(shí)驗(yàn)數(shù)據(jù),但考慮到三次諧波的幅值非常小(毫米量級(jí)),這一不足可以忽略。由于兩個(gè)模型的色散性以及變淺性能相同,因此這一差別反映了非線性對(duì)數(shù)值結(jié)果的影響。圖9也給出了T=3 s時(shí),模型BouN4D4和BouN2D4的計(jì)算結(jié)果同實(shí)驗(yàn)數(shù)據(jù)的對(duì)比。可以看出,模型結(jié)果同實(shí)驗(yàn)數(shù)據(jù)差別較大,目前所見(jiàn)各類(lèi)模型普遍存在這一現(xiàn)象[7-8,18]。不同于上述兩個(gè)算例,兩模型的數(shù)值結(jié)果幾乎沒(méi)有差別,這是由于該算例值最小,兩模型的非線性性能曲線在此處差別較小所致。因此,對(duì)于Whalin的實(shí)驗(yàn),當(dāng)周期較小,非線性較大時(shí),模型的非線性性能對(duì)數(shù)值結(jié)果的影響較大,總體而言,BouN4D4模擬的數(shù)值結(jié)果較BouN2D4更接近實(shí)驗(yàn)數(shù)據(jù)。

圖8T=1.0 s時(shí)諧波幅值對(duì)比Fig.8 Comparisonsof harmonic amplitude forwave periodT=1.0 s

圖9 不同時(shí)刻諧波幅值對(duì)比Fig.9 Comparisonsof harmonic amplitude forwave periodT=2.0 s(left)andT=3.0 s(right)

3 結(jié) 語(yǔ)

建立基于四階完全非線性Boussinesq水波方程的二維波浪破碎模型。選用的控制方程BouN4D4具有強(qiáng)非線性特征同時(shí)具有較好的色散性和變淺作用性能,更適合于描述波浪的非線性運(yùn)動(dòng)。在矩形網(wǎng)格上對(duì)上述方程進(jìn)行空間離散,采用高階導(dǎo)數(shù)近似方程中的時(shí)、空間導(dǎo)數(shù),用高精度ABM預(yù)報(bào)校正格式進(jìn)行時(shí)間積分,建立了高精度的數(shù)值求解格式,同時(shí)通過(guò)采用Kennedy等的渦粘破碎方法處理波浪破碎問(wèn)題。

將所建立的具有四階完全非線性特征的波浪傳播模型BouN4D4和僅具有二階完全非線性特征的模型BouN2D4用于模擬規(guī)則波在具有三維特征地形上的傳播過(guò)程。通過(guò)對(duì)比兩個(gè)模型的數(shù)值結(jié)果以及相關(guān)的實(shí)驗(yàn)數(shù)據(jù),驗(yàn)證了所建立模型的有效性,討論分析了非線性在波浪傳播過(guò)程中的作用。結(jié)果表明:1)實(shí)驗(yàn)數(shù)據(jù)同數(shù)值結(jié)果吻合良好,這表明本研究所建立的模型可以較準(zhǔn)確的模擬波浪在三維地形上的傳播、破碎過(guò)程;2)入射波浪非線性較強(qiáng)時(shí),模型的強(qiáng)非線性特征發(fā)揮較大作用,具有抑制波幅過(guò)度增長(zhǎng)的作用。

[1] Kirby J T.Boussinesqmodels and application to nearshorewave propagation,surfzone processes and wave-induced currents[M]∥advances in Coastal Engineering.Lakhan V C(ed),Elsevier,2003:1-41.

[2] 李孟國(guó),王正林,蔣德才.近岸波浪傳播變形數(shù)學(xué)模型的研究與進(jìn)展[J].海洋工程,2002,20(4):43-57.

[3] Madsen PA,Sφrensen O R.A new form of the Boussinesq equationswith improved linear dispersion characteristics.Part 2:A slowlyvarying bathymetry[J].Coast.Engrg.,1992,18:183-204.

[4] Nwogu O.An alternative form of the Boussinesq equations for nearshorewave propagation[J].J.Wtrwy,Port,Coast.and Oc.Engrg.,1993,119:618-638.

[5] Madsen PA,Bingham H B,Schaffer H A.Boussinesq-type formulations for fully nonlinear and extremely dispersive water waves:derivation and analysis[J].Proc.R.Soc.London.A,2003,459(2033):1075-1104.

[6] Lynett P,Liu PL F.A two layer approach to wavemodelling[J].Proc.R.Soc.London.A,2004,460:2637-2669.

[7] Fuhramn D.Numerical solutionsof Boussinesq equations for fully nonlinear and extremely dispersivewaterwaves[D].Denmark:Technical Univ.of Denmark,2006.

[8] 王本龍.基于高階Boussinesq方程的海岸破波帶數(shù)學(xué)模型研究[D].上海:上海交通大學(xué),2005.

[9] 張洪生,馮文靜,王亞玲,等.非線性波傳播的新型數(shù)值模擬模型及其實(shí)驗(yàn)驗(yàn)證[J].海洋學(xué)報(bào),2007,29(4):137-147.

[10] Wei G,Kirby J T,Grilli ST.A fully nonlinear Boussinesqmodel for surfacewaves.Part I.Highly nonlinear unsteadywaves[J].J.Fluid Mech.,1995,294:71-92.

[11] Kirby J T,Wei G,Chen Q,et al.FUNWAVE 1.0,Fully nonlinear Boussinesq wavemodel documentation and user′smanual[R].CACR search report,University of Delaware,1998:CACR-98-06.

[12] Madsen PA,Sch?ffer H A.Higher-order Boussinesq-type formulation for surface gravity waves:Derivation and analysis[J].Phil.Trans.R.Soc.Lond.A,1998,356:3123-3184.

[13] Gobbi F,Kirby J T.Wave evolution over submerged sills:testsof a high-order Boussinesqmodel[J].Coast.Engrg.,1999,37(1):57-96.

[14] Zou Z L,Fang K Z.Alternative forms of the higher-order Boussinesq equations:Derivations and validations[J].Coast.Engrg.,2008,55(6):506-521.

[15] Chawla A,Kirby J T.Wave transformationover a submerged shoal[R].CACR research report,University of Delaware,1996:CACR-96-03.

[16] Whalin RW.the limitof applicability of linearwave refraction theory on a convergence zone[R].Res.Rep.,Vicksburg,MS:U.S.Army Corpsof Engineers,Waterways Expt.Station.1971:Res.Rep.H-71.

[17] Kennedy A B,Chen Q,Kirby J T,et al.Boussinesqmodeling of wave transformation,breaking,and runup[J].I:1D.J.Wtrwy,Port,Coast.and Oc.Engrg.,2000,47:39-47.

[18] Chan Y,Liu PL F.Modified Boussinesq equations and associated parabolicmodels forwaterwave propagation[J].J.Fluid Mech.,1995,288:351-381.

附錄:

1 單方向差分格式

導(dǎo)數(shù)階數(shù)/m 所在點(diǎn) γ αi?3 αi?2 αi?1 αi αi±1 αi±2 αi±3i90 1i?1i?2i?3±60-12-98-10-147-24-77 360-45-35 150-450 80-100 400 45-30 50-225-15 72 1-12-10 2 3 4ii?1i?2i?3i i?1i?2i?3i i?1i?2i?3 180±8 6 2-12 137 812 1-1-15-49-151 7 35-27 228-147-3 132-8-85 6 232 12-18-84-186 270-420-255 5 265 13 35-83-461-39 21 171 411-490 200 470-5 080 0-48 64 496 56-4-184-484 270 15-285 2 970-13 29-29-307-39-9 111 321-27-12 93-972 8-88 104 12 6-36-114 22-13 137-11-1-15-1-151 7i4-505 5i?1i?2i?3±2-1-3-5-7 16 28 40-35-65-95 40 80 120-25-55-85-482 0 32 1-1-3-5

2 混合導(dǎo)數(shù)差分格式

混合導(dǎo)數(shù)差分格式可以表示為

式中:m和n分別表示沿x,y方向?qū)?shù)的階數(shù);Fm,n為系數(shù)矩陣,分別為

對(duì)于F1,2,F2,3,F1,3,F1,4則可通過(guò)對(duì)上面公式的簡(jiǎn)單變換得到。

A 2D numericalwavemodel based on fourth order fully nonlinear Boussinesq equations

FANG Ke-zhao,ZOU Zhi-li
(The State Key Laboratory of Coastal and Offshore Engineering,Dalian University of Technology,Dalian 116024,China)

TV148

A

1005-9865(2011)01-0032-09

2010-04-07

中央高校基本科研業(yè)務(wù)費(fèi)專(zhuān)項(xiàng)資金資助項(xiàng)目;大連理工大學(xué)海岸和近海工程國(guó)家重點(diǎn)實(shí)驗(yàn)室、河海大學(xué)海岸災(zāi)害及防護(hù)教育部重點(diǎn)實(shí)驗(yàn)室開(kāi)放基金資助項(xiàng)目;國(guó)家自然科學(xué)基金資助項(xiàng)目(51009018)

房克照(1980-),男,山東日照人,講師,博士,從事海岸工程研究。E-mail:kfang@dlut.edu.cn