王培濤,于福江,2
(1.國家海洋環(huán)境預報中心,北京 100081;2.國家海洋局海洋災害預報技術(shù)研究重點實驗室,北京 100081)
基于Boussinesq方程近岸波浪演變的數(shù)值研究
王培濤1,于福江1,2
(1.國家海洋環(huán)境預報中心,北京 100081;2.國家海洋局海洋災害預報技術(shù)研究重點實驗室,北京 100081)
首先對目前描述近岸波浪傳播變形的數(shù)學模型進行了回顧與總結(jié);對不同數(shù)學模型的特點、適用范圍和發(fā)展情況進行了闡述與對比。應用基于Boussinesq方程的Coulwave模式針對幾個經(jīng)典實驗地形進行了數(shù)值實驗,數(shù)值結(jié)果和實驗實測數(shù)據(jù)吻合較好。此外,分別采用不同的近岸波浪模型模擬了某漁港附近波浪的傳播變形,結(jié)果表明:當考慮波浪的折射、繞射、反射聯(lián)合作用時,Coulwave模式計算結(jié)果明顯較緩坡方程及SWAN模型計算結(jié)果更加合理。
Boussinesq型方程;近岸波浪模型;緩坡方程SWAN模型;Coulwave模型
近岸水域是與人類生活發(fā)生聯(lián)系最多的區(qū)域,大氣和海洋每年總交換熱量的近一半分布于近岸水域,如港口航道工程、近海取排水設(shè)施、海水養(yǎng)殖設(shè)施等,都集中分布在近岸海域內(nèi)。在海岸動力場中,海浪一直是施加在近岸建筑物上的最重要的環(huán)境荷載之一。當波浪由外海深水區(qū)傳播到海岸附近的淺水區(qū)時,受到水深、地形、底摩擦、障礙物(如建筑物、島嶼、岬角、沙洲、岸壁)等的影響,波浪的傳播速度、波長、波高、波面狀況、及其產(chǎn)生的流場和壓力場都會發(fā)生較大的變化。對于復雜地形和復雜邊界上波浪傳播變形的研究方法主要包括現(xiàn)場觀測、物理模型、數(shù)值模擬。然而伴隨著計算機軟件、硬件技術(shù)以及計算方法的快速發(fā)展,利用數(shù)值模型來模擬近岸波浪的傳播變形變得簡單有效,已逐漸成為近岸波浪研究的主要手段。
目前研究波浪在近岸區(qū)域傳播變形的模型主要有:基于光學snell定律的射線理論(ray theory)、基于微幅波理論的緩坡方程(mild-slope equation)、基于波能和波作用量守恒方程的譜模型、基于holmheltz方程的波浪繞射模型、基于直接求解Navier-stokes方程模型以及考慮波浪的非線性和色散效應的Boussinesq方程模型[1-7]。然而,上述幾種數(shù)學模型都依賴于特定的物理假設(shè)和數(shù)學近似,所以各自的適用范圍就受到了限制。
射線理論本身為線性理論,適用于緩變地形,且假定無能量跨過波向線,在海底坡度較大,地形較復雜時應用射線理論會出現(xiàn)焦散現(xiàn)象和盲區(qū)現(xiàn)象(波向線不能進入的區(qū)域)[8-9],但射線理論具有簡單省時省力的特點,計算精度為一階精度;拋物型緩坡方程僅考慮沿波峰線方向的繞射作用,而忽略波向線方向的繞射作用,且沒有考慮地形和障礙物的反射作用,在反射作用比較強的區(qū)域不適用。此外,也不適用于對于入射波向偏離主方向太大的情況;橢圓型緩坡方程變量為復變量,方程本身具有不可分離的性質(zhì),對它的離散求解不能應用迭代法,對于大的計算域求解比較困難;雙曲型緩坡方程沒有忽略波浪的反射特性,且計算精度與橢圓形緩坡方程相當,但其為時變型方程,為使模型收斂,需要用小的時間步長;譜模型用源函數(shù)的方式分別考慮了風應力作用、波浪破碎、底摩擦耗散、波-波相互作用等各種物理過程,可用于較大范圍較長時間的計算,但由于其在波浪折射、繞射等方面的不足,對不規(guī)則地形下近岸波浪計算,特別是港池、航道等復雜地形下的波浪傳播變形并不適用;基于holmheltz方程的波浪繞射模型是基于線性簡諧波及等水深水域的波浪繞射方程,可用于任意形狀港灣和防波堤掩護水域的波場計算[10-11];基于不可壓縮完整的三維Navier-stokes方程和連續(xù)方程的波浪數(shù)學模型可以描述全水深、任意坡度的線性和非線性的波浪傳播變形,并可計算速度的垂向分布,但將完整的三維Navier-stokes方程應用于實際還有很多工作要做。
Boussinesq方程用流速和水位描述水體變化[12],假定垂向速度沿水深為線性分布,由此得出水平速度和壓力沿水深為拋物型分布,考慮了波浪的非線性及色散性特征,運動保持質(zhì)量和動量守恒,因此能較好地模擬波浪傳播過程中多種因素的相互作用,其適用性廣泛。隨著計算方法的逐漸成熟,Boussinesq方程已得到越來越廣泛的應用與發(fā)展,特別是在拓展和改進模型的非線性及頻散特征方面取得了較大的進步。本文將重點介紹基于Boussinesq方程的Coulwave模型對近岸波浪傳播變形的模擬。
Boussinesq(1872)[13]假定水平速度上下均勻,垂向速度由底到自由表面是線性分布,得出了一維非線性控制方程,稱為Boussinesq方程。Peregrine(1967)[14]推導出了變水深條件下的二維Boussinesq方程,稱為經(jīng)典Boussinesq方程:
式中:x、y為與靜止水面重合的直角坐標系坐標,u(x,y,t)、v(x,y,t)分別為水深平均的質(zhì)點速度矢量沿x、y方向的分量,ζ為波面到靜止水面的距離。
由于經(jīng)典的Boussinesq方程包含非線性色散性,能夠反映和描述近岸區(qū)域波浪的變形現(xiàn)象,被認為是短波數(shù)值模擬領(lǐng)域的一個重大突破。但Boussinesq方程具有弱色散性、弱非線性以及方程本身沒有考慮底摩擦、波浪破碎和環(huán)境水流的影響,其適用性也受到了限制,適用水深范圍h/L≤0.12。近20年來,Boussinesq方程的理論和應用都得到了快速發(fā)展,特別是在色散性及變淺作用、非線性、波浪破碎作用、環(huán)境流速等方面改進顯著,從而拓展了其水深適用范圍。
Madsen等采用在動量方程中加入含有1個待定系數(shù)常數(shù)Boussinesq三階導數(shù)項的方法得到了1個常水深情況下的用單寬流量(或稱水深積分速度)和波面表達的二維改進型Boussinesq方程,在進一步的工作中,Madsen等將其推廣到適用緩變坡度變水深情況[15]。當B=1/15時,改進型Boussinesq方程的線性色散性能達到Airy波精確解的Pade[2,2]階近似。Madsen等同時提出了線性變淺梯度(linear shoaling gradient)概念作為衡量改進型Boussinesq方程變淺作用性能的指標。Madsen等的改進型Boussinesq方程的變淺作用性能的精度達到O(μ4),能夠適用于h/L≤0.5。在這一工作中,Madsen還詳細討論了方程色散性與方程速度變量選擇的關(guān)系,指出了方程用不同的速度變量表達時,所得的色散關(guān)系是不同的,其中,用水深平均速度給出結(jié)果最好,而用自由表面速度結(jié)果最差。這一思想在Nwogu(1993)的研究中得到了很好的發(fā)展,他取任意水深處速度表達Boussinesq方程[16],并通過與線性色散方程的Pade’,展開擬合,獲得最優(yōu)速度,所得方程色散關(guān)系與Witting的方程和Madsen的方程一致。Liu(1994)[17]和Wei等(1995)[18]擴展Nwogu方法于完全非線性波,導致模型不僅可以適用于中等水深,也適用于模擬強非線性相互作用。Gobbi(2000)[19]使用四階多項式推導的模型,線性色散性精度到kh=6。Madsen等(2002)基于Agnon(1999)的方法提出了一個Boussinesq類模型,準確性達到深水(kh=40)。本方法采用Laplace方程的最優(yōu)展開,通過使用水柱不同層的多級展開,只需要五階空間導數(shù)就可以達到深水。Lynett[20-23]和Philip Liu提出了基于Boussinesq方程的波浪分層數(shù)學模型,假定在每一層中的垂向流速按線性分布,分布在各層對垂向流速進行垂向積分,求其平均流速,用這一平均流速來代替每一層的垂向流速。這樣的處理既提高了方程的精度,又改進了方程的色散性,使方程能應用于更大的水深范圍。圖1是Boussinesq方程各階段主要研究階段色散關(guān)系精度對比情況。
圖1 幾種改進后的Boussinesq方程線性色散關(guān)系與Airy波色散關(guān)系精確解比較
Coulwave模型在垂直方向上將水柱分解成多層,采用深度積分模型來代替對垂直流場的高階多項式近似,這種方法不僅提高了Boussinesq方程的計算精度,而且對方程的線性色散精度也進一步改善,更重要的是方程的物理階數(shù)卻沒有提高,避免了數(shù)值求解高階導數(shù)項帶來的數(shù)值振蕩問題。經(jīng)拓展后的Coulwave模型在源函數(shù)造波、底摩擦效應、波浪破碎、海綿層濾波方面得到了完善,模型可以用來模擬水波淺化、折射、繞射和反射等綜合現(xiàn)象,特別是在刻畫結(jié)構(gòu)物附近的波浪場,有其顯著的優(yōu)勢。
本節(jié)將以兩層模型為例進行介紹及應用。二層的Boussinesq方程模型在一層模型基礎(chǔ)上,對水柱面再進行分層,在每一層中的垂向流速都按線性分布,分別對其進行沿水深方向積分求其平均流速,用此平均流速代替各層的流速(見圖2)。
圖2 二層模型示意圖
對于兩層模型來說,定義其水平流速向量為:
模型的連續(xù)方程及動量方程為:
其中u2可以用u1的函數(shù)關(guān)系式表達如:
底層流速可以通過上層流速計算求解。式(6)—(8)即為二層系統(tǒng)的聯(lián)合控制方程。
為了驗證Coulwave模型的可用性,本文針對2個經(jīng)典試驗地形下波浪演變進行了數(shù)值研究,并將計算結(jié)果與實驗數(shù)據(jù)進行了對比。
3.2.1波浪跨越潛堤的傳播和變形
Dingemans(1994)的實驗模型(見圖3),左端x=0(h0=0.86 m)的位置設(shè)置造波機,取空間步長dx=0.10 m,時間步長Δt=0.0189 s,模擬時間120 s。在模擬中左右兩端分別設(shè)置了1.25 L的海綿層。在圖示的四個靶位放置測波儀測量波面時間序列。x1=26.04 m,x2=33.64 m,x3=37.04 m,x4=41.04 m
圖3 Dingemans實驗模型示意圖
在實驗中,入射波(波幅H=0.02 m,波周期T= 2.86,ε=a0/h0=0.023kh=0.7,L=7.7 m,h=0.86 m)傳播經(jīng)過#1位置后,這時由于水深變淺,波浪發(fā)生淺水效應,使波浪變陡,同時由于非線性作用的影響,波浪的能量不斷從主頻傳向高頻,波浪中的高頻諧波將不斷增加。由于各諧波的相速度不同,以至于在堤后,高次諧波被釋放出來,波浪的形狀及相速度發(fā)生改變。我們分別用一層模型和兩層模型對這次傳播過程進行了數(shù)值模擬。從模擬結(jié)果圖4來看:在坡頂(#1位置處)出現(xiàn)了明顯的非線性特征,表現(xiàn)為波峰變尖,波谷變平。當波浪到達坡頂附近時,波幅達到最大值,出現(xiàn)明顯的次峰,次峰隨著地形的變化而變化。從圖4可以看到對于#1位置,兩種方法模擬結(jié)果與實驗結(jié)果吻合良好。兩層模型與一層模型模擬的效果相當,并沒有顯示出兩層模型的優(yōu)勢所在。這是因為此時的kh<2.0,一層、兩層模式的色散性都能保證合理的精度。當波列行進到坡后,特別是#3,#4位,實驗結(jié)果與模擬結(jié)果有了一些差異,這是由于當高次諧波釋放后形成較短的波,μ=k0h0值較大,在坡后μ=kh≈4,超出了一層模型的預報精度,同時波浪在堤后的變形更加劇烈,這時對模型非線性和色散性精度要求更高,兩層模型模擬的結(jié)果與實驗數(shù)據(jù)更加吻合,說明采用兩層處速度近似可以大大提高方程的色散性能(見圖4)。
圖4 波面高度歷時變化與實驗數(shù)據(jù)對比
3.2.2 波浪在Berkhoff實驗地形上的傳播演變
1982年Berkhoff等設(shè)計和進行了著名的橢圓緩坡地形波浪傳播、變形的物理模型實驗[24]。本實驗可反映出波浪在淺水變形、折射、反射和非線性彌散關(guān)系等諸多因素影響下的變形和傳播情況,已經(jīng)成為數(shù)值模型計算驗證的經(jīng)典實驗之一。為了驗證本文模型在波浪淺水效應、反射、折射和非線性彌散等諸多因素共同作用下的適用性,本小節(jié)將分別用本文模型和拋物型緩坡方程對實驗地形進行了數(shù)值模擬,并與實驗數(shù)據(jù)進行了比較。
具體水下地形和實驗測量斷面分布見圖5,計算區(qū)域為25 m×20 m,整個水域底坡為1∶50,斜坡梯度方向與波浪入射邊界的法向夾角為20°。以x0=10m,y0=10m為中心,有一長半軸為4 m,短半軸為3 m的橢圓淺灘。斜坡旋轉(zhuǎn)后的坐標為:
圖5 Berkhoff地形實驗設(shè)置
圖6 計算所得T=10 s、30 s時刻三維波面
(x0,y0)位于淺灘中心,水深為0.1332 m。淺灘的邊界為:
淺灘外水深如下:
淺灘上水深由下式給出
入射波波高為0.0464 m,波浪周期T=1.0 s,入射方向從x正方向入射。計算空間步長為0.05 m,時間步長dt=0.01s,總時間步數(shù)為3000,即30個波周期。模擬啟動的初始條件為“冷啟動”,在x=1.5 m處設(shè)置波源,左、右兩端設(shè)置l=1.5 m的海綿濾波層,上下兩個側(cè)邊界為全反射邊界(見圖5)。
從模擬的結(jié)果(見圖6)可以看出,當入射波前行遇到斜坡和橢圓形淺灘后開始變形,在斜坡作用下波浪發(fā)生折射,波峰線開始變得不再平行,一側(cè)波長變短波高增加,一側(cè)波長變長波高減小。當入射波傳播至淺灘后方,由于兩側(cè)波浪繞過淺灘后的疊加,波能集中,波高顯著增加。從斷面比波高對比來看(見圖7),緩坡方程及改進的Boussinesq方程模型對于波浪在緩坡地形上的傳播變形都能得到較滿意的結(jié)果;但對于斷面6及斷面8后半段而言,本文的模型跟實驗數(shù)據(jù)吻合更好,究其原因可能為斜坡和淺灘引起波浪復雜的波-波相互作用所致。
上一節(jié)所述的幾個例子都是在邊界相對簡單,地形不太復雜的情況下進行的,實際的海底地形和岸邊界情況往往比這復雜的多,計算區(qū)域也大得多。本節(jié)將利用本文模式和緩坡模式、SWAN模式對某漁港(見圖8)附近區(qū)域波浪傳播進行模擬對比。同時也將本文模型模擬波浪在漁港內(nèi)的演變(見圖9),以及利用上述三種近岸模式模擬外海波浪在特定工程點處的統(tǒng)計波高情況進行了研究,選取經(jīng)常侵襲港口水域E(X=0)方向的入射波進行模擬。計算區(qū)域1500 m×1500 m,波長(L)150 m,波高(H)2.5m,周期(T)12.542s,計算網(wǎng)格步長Δx=Δy=5 m,時間步長Δt=0.2s。側(cè)邊界采用海綿層濾波。
圖7 各個斷面比波高對比
圖8 漁港三維地形圖及控制點分布(單位:m)
圖9 T=30 s、90 s、120 s、200 s波浪演進三維立體圖
表1三種模型工程測點波高計算對比(單位/m)
不同時刻波面變化(見圖9)中可以看到:波浪傳播到口門防波堤前趾時發(fā)生繞射,部分波能向港池內(nèi)擴散。防波堤前方很明顯看到,入射波和反射波疊加,波幅明顯增大。隨著時間推移,反射波的影響范圍越來越大。同時由于水深和地形的影響波向發(fā)生折射,波面變得越來越不規(guī)則,在堤前區(qū)域形成短峰波。由于防波堤的存在,在港池內(nèi)波幅明顯減小。在波浪折射、繞射、反射聯(lián)合作用下,可以看到在槽道內(nèi)出現(xiàn)了明顯的波幅。港池內(nèi),波浪傳播遇到直立堤,同樣發(fā)生反射與進行波疊加,使得波幅增大。
對比表1可以知道,在開闊水域(1、2、3、4)三類模型計算的波高相差不大,在反射作用較強的區(qū)域(5、6、7、8),只有Coulwave模式模擬的結(jié)果比較真實,緩坡方程次之,SWAN模式基本沒有表現(xiàn)出任何反射效應。當波浪繞射作用較強時(9、10、11),Coulwave模式的優(yōu)勢就更加明顯,其他兩類模型未能刻畫出波浪在堤后的繞射效應。而對于能夠反映波浪折射、繞射及反射聯(lián)合作用的控制點(12—18),緩坡模型及swan模型計算結(jié)果在量級上與Boussinesq方程模型偏差較大。
Coulwave模式全面考慮了波浪的折射、繞射、反射效應,而拋物型緩坡方程沒有考慮波浪的反射效應,測點處得到的波高小于Coulwave的結(jié)果就在意料之中了。而Swan模型對波浪的繞射、反射效應模擬能力相對比較差,得到的結(jié)果小于其他兩類模式的結(jié)果,并且Swan模型在模擬規(guī)則波的反射、折射、繞射方面均不如前兩類模式模擬的效果好,特別是波高在一兩個波長水平范圍內(nèi)變化比較大時,模型的計算結(jié)果將出現(xiàn)更大偏差。因此,Swan模式在模擬由近岸海底地形和建筑物引起的波浪繞射和反射效應時會產(chǎn)生較大的誤差,這時一般不采用此類模式。
Coulwave模型采用分層Boussinesq方程方法來改進經(jīng)典的Boussinesq方程的色散性,使得改進后的模型的水深范圍得到了極大拓展。從本文的模擬結(jié)果可以看出:基于Boussinesq方程的Coulwave模型對復雜邊界問題處理能力與緩坡方程模型、Swan模型相比具有明顯優(yōu)勢,同時Coulwave模型很好的模擬了外海波浪傳播到近岸過程中出現(xiàn)的折射、反射、繞射、淺水效應和破碎等波浪在近岸發(fā)生的變形現(xiàn)象,模擬結(jié)果明顯優(yōu)于另外兩種波浪模型,特別是對波浪反射過程的模擬,模型更加適用于復雜地形水域(有人工島、防波堤等建筑物存在的區(qū)域)附近波浪傳播變形的模擬研究。
由于海洋中實際波浪一般呈現(xiàn)不規(guī)則方向譜形式,波浪不僅具有不規(guī)則性,而且具有多向性。本研究僅限于對規(guī)則入射波,且與地形的聯(lián)合作用采用二維兩層Boussinesq模型進行了討論,對不規(guī)則方向譜的模擬未能涉及。研究表明分層Boussinesq波浪模型中層數(shù)分得越多,其色散性和非線性的效果體現(xiàn)越好,計算模型的精度越高,適用水深范圍越廣。因此對不規(guī)則波入射傳播及更高層數(shù)的Boussinesq方程模型的探討將是今后研究的重點。
[1]鄒華志.近岸波浪的Boussinesq型方程模型及其應用[J].廣東水利水電,2009(4):1-5.
[2]李孟國,王正林,蔣德才.近岸波浪傳播變形數(shù)學模型的研究與進展[J].海洋工程,2002,20(4):43-57.
[3]文圣常,余宙文.海浪理論與計算原理[M].北京:科學出版社, 1984.
[4]李孟國,蔣德才.近岸波浪傳播折射變形的數(shù)學模型概述[J].海岸工程,1999,18(4):100-109.
[5]李春穎.基于Boussinesq方程的海岸地區(qū)波浪數(shù)學模型研究[D].天津:天津大學,2003.
[6]王培濤.基于Boussinesq方程的近岸波浪模型的研究及應用[D].北京:國家海洋環(huán)境預報中心,2007.
[7]王海瓏.基于Boussinesq方程波浪分層數(shù)值模擬方法[D].天津:天津大學,2004.
[8]Arthur R S.Refraction of water waves by islands and shoals with circular bottom contours[J].Trans.Amer.Geophys.Un.,1946,27 (2):168-177.
[9]Munk W H,Arthur R S.Wave intensity along a refacted ray.In: Gravity waves[J].Natl.Bur.Standards,1952,Circular521:98-108.
[10]龔崇準,戴功虎.雙突堤波浪繞射的數(shù)學模型[J].海洋學報, 1984,6(1):99-116.
[11]龔崇準,張長寬.任意反射率邊界條件的港灣波浪繞射數(shù)學模型[J].河海大學學報,1987,15(6):72-79.
[12]鄒志利.水波理論及其應用[M].北京:科學出版社,2005
[13]Boussinesq J.Theory of wave and swells propagated in long horizontal rectangular canal and imparting to the liquid contained in this canal[J].Journal de Mathmatiques Pures et Appliquees,1872,17(2):55-108.
[14]Peregrine D H.Long waves on a beach[J].J Fluid Mech.,1967, 27:815-827.
[15]Madsen P A,Sorensen O R.A new form of Boussinesq equations with improved linear dispersion characteristics:PartⅡ.A slowly varying bathymetry[J].Coast Eng,1992,18:183-204.
[16]Nwogu O.Alternative from of Boussinesq equations for nearshore wave propagation[J].Journal of Waterway,Port,Coastal and Ocean Engineering,1993,119(6):618-639.
[17]Liu P L F,Model equations for wave propagation from deep to shallow water[J].Advances in Coastal Engineering,1994,1:125-157.
[18]Wei G,Kirby J T,Grilli S T,et al.A fully nonlinear Boussinesq model for surface waves:Part I.Highly nonlinear unsteady waves[J].J.Fluid Mech.1995,294:71-92.
[19]Gobbi M F,Kirby J T.Wave evolution over submerged sills: tests of a high-order Boussinesq model[J].Coast Eng,1999,37: 57-96.
[20]Lynett P J,Wu T R,Liu P L F,Modelling wave runup with depth-integrated equations[J],CoastalEngineering,2002,46: 89-107.
[21]Lynett P,Liu P L F,Hwung H H,et al.Multi-layer modeling of wave groups from deep to shallow water,Proceedings of OMAE 2003[R].Ocean Engineering Symposium,2003,13:8-13.
[22]Lynett P,Liu P L F.A two-dimensional,depth-integrated model for internal wave propagation[J].Wave Motion in press.2002,36: 221-240.
[23]Lynett P.A multi-layer approach to modeling generation propagation and interaction of water waves[D].USA:CornellUniversity,2002.
[24]Berkhoff J C W,Booy N,Radder A C.Verification of numerical wave propagation models for simple harmonic linear water waves[J].Coastal Engineering,1982,6:255-279.
Numerical study on the nearshore wave evolution based on Boussinesq-type equations
WANG Pei-tao1,YU Fu-jiang1,2
(1.National Marine Environmental Forecasting Center,Beijing 100081,China;2.Key Laboratory of Research on Marine Hazards Forecasting, State Oceanic Administration,Beijing 100081 China)
The paper reviewed several advanced mathematical models that can be used for the description of wave deformation in port and coastal engineering.The different features of mathematical models,the applicable scale and the evaluation were discussed.The Coulwave(Cornell University Long and Intermediate Wave Modeling)model,based on Boussinesq-type equations,is used to simulate the wave propagation and deformation for several classic topography experiments with specific incident waves.The numerical results agree well with the experimental data.In addition,we compared the results of Coulwave with the SWAN and models with mild-slope equation.The results show that:when considering the interaction of the wave refraction,diffraction, reflection,Coulwave model shows high performance than the model based on mild-slope equations and SWAN model.
Boussinesq-type equations;nearshore wave models;mild-slope equations;SWAN model;Coulwave model
P731
A
1003-0239(2012)04-0007-17
2011-12-29
海洋公益性行業(yè)科研專項(200905013)
王培濤(1981-),男,助理研究員,主要從事海嘯、風暴潮理論及預警報技術(shù)研究。E-mail:wpt@nmefc.gov.cn