劉必勁，張振偉，劉忠波，傅丹娟，陳小云

(1.廈門理工學院 土木工程與建筑學院，福建 廈門 361024；2.大連海事大學 交通運輸工程學院，遼寧 大連 116026)

1 引言

畸形波是一種非線性很強的波浪，其波高很大，通常是有效波高的2.2 倍以上，這種波浪對海上建筑物或過往的船只帶來嚴重威脅。為了深入掌握畸形波的形成機制，國內外學者通過理論分析、數(shù)值模擬或實驗模擬等方式展開了數(shù)十年的研究。2003 年，Kharif 和Pelinovsky[1]較為全面地分析了其形成機制，認為在相位聚焦、波?流相互作用、地形變化、表面風、波浪不穩(wěn)定性等影響下都可能產生畸形波。

為了研究畸形波，實驗室中一般都采用相位聚焦的方式產生較大的波浪，如Baldock 等[2]在二維水槽中采用線性相位聚焦的方式產生了聚焦波；柳淑學和洪起庸[3]在三維水池中開展了三維極限波浪的實驗研究。在理論分析方面，裴玉國等[4]采用線性模型模擬了定點生成畸形波的方法；趙西增等[5]則借助于聚焦波和不同形式波浪的組合模擬畸形波。在數(shù)值模擬方面，國內外學者則通過利用N-S 方程模型、非靜壓模型、Boussinesq 方程數(shù)值模型、勢流方程模型、高階譜模型等不同的數(shù)學模型模擬聚焦波。如趙西增等[6]利用高階譜方法建立了模型，并應用于畸形波的模擬；Li 和Liu[7]建立了非周期邊界的高階譜聚焦波模型，并研究設定聚焦波幅、中心頻率、頻率寬度、波幅分布以及方向角等對聚焦波的影響，研究表明，寬的方向譜減弱了非線性聚焦波幅；Ai 等[8]建立了不同網(wǎng)格下的非靜壓模型，并模擬了水槽和水池中的聚焦波；Ning 等[9]基于勢流理論建立了高階邊界元數(shù)值模型，模擬了深水聚焦波；Li 等[10]在CIP 模型（利用CIP 求解對流項）中引入內部造波法，模擬了潛堤上聚焦波的傳播變形問題。

關于Boussineq 模型在聚焦波模擬的研究文獻相對較少，能成功模擬聚焦波速度分布特征的僅有少數(shù)Boussinesq 水波方程，這是因為多數(shù)Boussinesq 方程本身在色散性、非線性，特別是速度分布特征方面存在嚴重不足。分析相關數(shù)值模擬聚焦波的研究可知，準確模擬波面?zhèn)鞑プ冃涡枰獢?shù)值模型具有良好的色散性與非線性特征，準確模擬速度分布必須要求模型具備良好的速度分布特征。當前能具備這些性能的Boussineq 數(shù)值模型并不多，在分析大量Boussineq方程文獻基礎上，本文選用了一組單層Boussinesq 水波方程，即Liu 和Fang[11]給出的雙層Boussineq 模型的簡化版。事實上，該簡化版本的模型在平底情況下與Madsen 等[12]給出的Boussineq 水波方程完全一致。本文將針對單層Boussinesq 水波方程建立數(shù)值模型，并考察其是否能適用于模擬強非線性波浪傳播和聚焦波的演化問題。

2 Boussinesq 數(shù)值模型

2.1 控制方程

假設水體無旋、無黏，在保留垂向速度下，Liu 和Fang[11]經嚴格數(shù)學推導給出了適合緩變地形的雙層Boussinesq 方程，將其簡化成為單層方程。立面二維模型由7 個方程組成，含有7 個待求變量：波面位移η，自由表面處水平和垂向速度uη和wη，靜止水位處的水平和垂向速度u0和w0以及計算水平和垂向速度u*和w*，方程的具體表達形式如下：

（1）自由面 η處的運動學邊界條件和動力學邊界條件為

（2）自由面速度（uη，wη）與靜止水位處速度（u0，w0）的關系為

（3）靜止水位處速度（u0，w0）與計算速度（u*，w*）的關系為

（4）水底處滿足的運動學邊界條件為

式中，g是重力加速度；下標x表示對x求導；h表示水深；其中(α,β12,β13)=(0.5,0.95,0.55)。

整個流場的表達式分為兩部分，從靜止水位到波面處的速度表達式為

式中，z表示原點位于靜水面時向上為正的垂向位置。從水底到靜止水位的速度表達式為

方程（1）至方程（7）構成單層Boussinesq 水波方程，圖1 給出了變淺梯度與解析解的對比以及無因次相速度。由圖可見，在2%誤差范圍內，相速度可適用于最大水深k0h=10（k0為深水波數(shù)，h為水深）[12]；速度分布特征在1%誤差范圍內適用的最大水深為k0h=3.5；變淺性能可適用0＜k0h＜10 范圍①林鵬程，劉忠波，劉勇.基于Boussinesq 數(shù)值模型的波浪速度垂向分布模擬研究.海洋湖沼通報，已采用.。

2.2 數(shù)值模型的求解

以上垂向二維Boussinesq 方程中的最高空間導數(shù)是3 階，模型本身與開源代碼FUNWAVE 的最高空間導數(shù)是一致的[13]，為此時空差分格式均采用與該開源代碼一樣的格式。在時間差分格式上采用混合4 階Adams-Bashforth-Moulton 格式，預報時利用3 階Adams-Bashforth 格式求解方程（1）和（2），可以得到波面位移以及相應的波面處水平速度的預報值；分別求解方程（3）和（5）可得到2 個水平速度的預報值；求解方程（7）、（6）和（4），可依次得到3 個垂向速度的預報值。校正時利用4 階Adams-Moulton 時間步進格式，求解得到所有變量的校正值。當7 個變量的校正值與預報值之間的誤差在0.000 1 內，當前迭代結束，否則用校正值和預報值的加權重新賦值給變量的預報值，本文算例中校正加權值R為0.1～0.2，而預報加權為1-R，并多次重復校正過程，直到滿足條件為止。以上過程可以參見文獻[13-14]。

圖1 方程的變淺梯度（a）和無因次相速度（b）Fig.1 The shoaling gradient (a) and non-dimensional (b)phase celerity of the present model

3 數(shù)值模型的驗證

為了驗證數(shù)值模型，我們采用流函數(shù)波浪的數(shù)值解作為近似的解析解來考察模型的適用情況。數(shù)值計算長度設置為10L（L為流函數(shù)波浪的波長），在入射邊界2 個波長范圍內用流函數(shù)波浪的解析解作為模擬的驅動條件，并在該范圍內采用了Fuhrman[15]提出的邊界松弛造波法，這相當于將近似解析解逐漸轉化為方程數(shù)值解的一種方式。在出口邊界，設置了2L的海綿邊界層進行消波。在固定波浪周期為6 s情況下，針對水深為20 m、30 m 和40 m 情況下，對波高為6 m 的流函數(shù)波浪進行了數(shù)值模擬研究。3 種情況對應的流函數(shù)波長分別為60.469 m、61.524 m 和61.679 m，無因次水深kh分別為2.078、3.064 和4.074；對應的波陡H/L=0.099、0.098 和0.097；這幾組工況的波高約是極限波高的0.7、0.65 和0.64 倍左右，因此屬于強非線性波浪算例。

圖2 水深為20 m 的計算結果與解析結果的比較Fig.2 Comparisons of analytical solution and numerical simulation for water depth 20 m

圖3 水深為20 m 波峰面下的水平速度分布與解析解的比較（x=4L）Fig.3 Comparisons of velocity profile between numerical simulation and analytical solution for water depth 20 m (x=4L)

每一個工況計算時間為300 s，空間計算步長取為L/32，時間計算步長取為0.05 s。將計算得到的波面位移、波面處水平速度與垂向速度以及沿垂向分布的水平速度與近似解析解（流函數(shù)波浪的數(shù)值解析解）進行了比較，具體結果見圖2 至圖7。圖2、圖4和圖6 給出了計算波面位移、波面處水平速度與垂向速度與近似解析解的比較；圖3、圖5 和圖7 給出了x=4L波峰面下水平速度的垂向分布與解析解的比較。由圖2、圖4 和圖6 可見，波面以及波面處速度與解析解誤差較小，為了明確誤差的具體值，計算了波長與解析解的誤差（從2L～8L范圍內標定），3 種情況下僅第一組誤差是0.52%，其他組幾乎不存在誤差；以2L～8L范圍內變量最大值的平均值為統(tǒng)計，最大波峰面與解析解的誤差、最大水平速度、最大垂向速度與相應的解析解誤差均在1%范圍內。進一步，從x=4L波峰面下水平速度（圖3，圖5，圖7）對比來看，伴隨水深的加大，波峰面水平速度分布的精確度逐漸變差，這與文獻中給出的線性Stokes 波的速度分布與解析解在kh=3.5 的積分誤差是1%的規(guī)律是類似的，這里不再進行速度分布的積分誤差的詳細計算。

圖4 水深為30 m 的計算結果與解析結果的比較Fig.4 Comparisons of analytical solution and numerical simulation for water depth 30 m

4 聚焦波演化的模擬

圖5 水深為30 m 波峰面下的水平速度分布與解析解的比較（x=4L）Fig.5 Comparisons of velocity profile between calculated and analytical solution for water depth 30 m (x=4L)

受水槽長度的限制，實驗室中一般采用線性相位聚焦方式產生較大的聚焦波進行物理模型實驗。在深水情況下，由于受波浪非線性的影響導致實際聚焦位置與設置的線性聚焦位置存在向后偏移，對應的聚焦時間也會向后偏移。這種偏移是因為聚焦波的波幅增大后，線性波的波長不能詮釋聚焦波波幅增加后引起波長增加的現(xiàn)象，其本質就是線性色散關系表達式已經不適用于闡釋強非線性的波長，此時應考慮波幅引起色散關系的變化，即波幅離散引起頻散。關于這一點在上面的算例中已有相關印證，以第3 個工況為例，線性波長為56.192 m，而流函數(shù)波浪的波長為61.679 m，二者誤差為9.7%。聚焦波的非線性增強，則這種偏移量越大[2]。

圖6 水深為40 m 的計算結果與解析結果的比較Fig.6 Comparisons of analytical solution and numerical simulation for water depth 40 m

圖7 水深為40 m 波峰面下的水平速度分布與解析解的比較（x=4L）Fig.7 Comparisons of velocity profile between numerical simulation and analytical solution for water depth 40 m (x=4L)

Baldock 等[2]開展的深水實驗含有不同波浪頻率段范圍，在不同周期內等分周期成29 份，每個頻率的波幅取值相等，其中B 組（T=0.6～1.4 s，寬譜）和D 組（T=0.8～1.2 s，窄譜）作為主要參考用來考察本文數(shù)值模型，B 組和D 組的線性無因次水深kh分別為1.568～7.825 和2.026～4.403。為了模擬聚焦波的生成過程，本文在入射邊界處引入線性解疊加方式，具體做法：給定每一個組成成分的波面位移的時間歷程，利用線性化的Boussinesq 方程確定出其他所有變量的時間歷程值，詳細過程可參見文獻[14,16]。設置計算域為24 m，空間計算步長采用0.04 m，時間計算步長為0.01 s，末端設置4 m 的海綿邊界層消波，計算時間為40 s。對B 組和D 組的3 種不同設定的聚焦波峰面為22 mm、38 mm 和55 mm 工況進行數(shù)值模擬，計算的聚焦位置處的波面位移和最大波峰面下的速度分布與實驗結果進行了比較，具體見圖8 和圖9。整體而論，計算波面位移與實驗結果吻合程度較好，水平速度的分布與實驗結果基本吻合。寬譜情況下（B 組），計算的水平速度出現(xiàn)“S”形，圖9 與圖7 的結果有相類似的規(guī)律，這是因為數(shù)值模型在線性無因次水深kh＞3.5 后誤差將逐漸增大，在kh＜3.5 只有18 個組成波在這個范圍內；窄譜情況下（D 組），線性無因次水深kh＜3.5 范圍內22 個組成波，因此D 組工況比B 組工況模擬效果要好。關于D22 算例數(shù)值結果與實驗結果吻合程度不好的原因，本文也嘗試了采用兩個不同的雙層Boussinesq 水波方程模型[14,16]；Liu 等[14]數(shù)值模擬結果與本文計算結果相差不大，也存在與實驗偏差較大的情況，鑒于這一組的波浪非線性相對較弱，因此可以初步認為有可能是實驗結果存在一定的偏差。

對于設定聚焦波面為55 mm 的情況，固定波浪中心頻率不變，即保持1 s 周期不變，設定線性聚焦位置為8 m，30 s。在周期范圍△T分別為0.8 s、0.7 s、0.6 s、0.5 s 和0.4 s 情況下聚焦波時的波峰面、波峰面速度、聚焦時間和空間差等在表1 中給出，時空偏移是數(shù)值模型（或實驗模型）得到的實際聚焦位置與線性波群設定的聚焦時間和位置之間的偏移量。

5 結論

本文建立了基于單層Boussinesq 水波方程的垂向二維數(shù)值模型，首先與流函數(shù)波浪解析解進行了對比，并將這個數(shù)值模型應用于聚焦波的模擬中，研究得出主要結論如下：

圖8 聚焦的計算波面位移與實驗波面位移的比較Fig.8 Comparisons of surface elevations between modeled and experimental data

圖9 聚焦波峰面下的水平速度與實驗結果的比較Fig.9 Comparisons of velocity profiles between modeled and experimental data

表1 計算結果與實驗結果的比較Table 1 Comparisons of results between modeled and experimental data

（1）數(shù)值模型計算的波面位移、波面位移處的水平速度與垂向速度均與流函數(shù)的解析解有較好的吻合，同時垂向速度在適用范圍內與解析解的比較，也展示出本文建立的數(shù)值模型能夠勝任模擬強非線性波浪。

（2）聚焦波演化中，波浪非線性性起重要作用，本文Boussinesq 型方程具有良好的色散性與非線性性能，保證了方程具有模擬聚焦波的能力；數(shù)值模擬波面與實驗結果的對比，展現(xiàn)了本文Boussinesq 數(shù)值模型能夠模擬聚焦波波面演化。

（3）與精確模擬的波面相比而言，在模擬速度分布的精度方面，在理論方面，方程的色散適用范圍明顯大于水深速度分布的適用范圍。數(shù)值模型受限于線性速度分布的精度，其模擬非線性速度分布的能力也明顯受限于kh＜3.5。盡管數(shù)值模型可以較為精確模擬波面以及波峰面處的水平速度與垂向速度，但對于kh＞3.5 的組成成分較多的聚焦波，基于本文Boussinesq 水波方程的數(shù)值模型無法精準模擬水平速度分布情況。

本文明確了水平速度分布的非線性適用范圍，在此基礎上可期望將模型拓展到研究聚焦波與建筑物相互作用問題。此時只需將速度表達式（公式（9）至公式（12））代入垂向歐拉方程中，進行兩次垂向積分，最終可以得到波浪力，關于這一方面將在下一步進行研究。