淺水中質量源造波方法

田正林，孫昭晨，梁書秀

(大連理工大學 海岸及近海工程國家重點實驗室，大連 116024)

淺水中質量源造波方法

田正林，孫昭晨，梁書秀

(大連理工大學 海岸及近海工程國家重點實驗室，大連 116024)

針對質量源造波過程中源區上方水體高度難以確定的問題提出了改進方法，將源區高度值考慮到源強方程中，避免了造波過程中源區上方水體高度值多次試算帶來的不便。該方法彌補了質量源造波方法在水深小于1.5 m的淺水中所造波浪與理論波浪之間差距較大的缺陷，尤其是水平網格步長與源區高度值相比小于10的情況，推廣了質量源造波方法在淺水中的應用。文中驗證了波高對所造波浪的影響以及時間步長對淺水中所造波浪的影響，得出結論：當波高與水深的比值小于0.2時所造波浪波高與理論值的擬合較好；當時間步長小于0.005 s時所造波浪無明顯的上下漂浮現象，波形穩定，符合理論波形要求。

質量源造波；淺水造波；源區高度；改進方法

數值波浪水槽是港口海岸工程數值模擬的基礎，造波系統與消波系統是數值波浪水槽的重要組成部分，造波系統造波能力及消波系統消波能力又是數值波浪水槽可用性的關鍵。大多數值模擬采用加長計算域或是波動沒有達到準穩定狀態時結束模擬避免二次反射波的影響[1]。為了得到長時間穩定的數值模擬，對數值水槽的要求是一能夠降低二次反射波，二所造目標波浪不受反射波的影響，流體域內部造波方法能夠很好地解決該問題。Larsen和Dancy[2]在水域內部沿一定長度的直線段上通過擾動水面生成目標波浪，避免了二次反射波的影響。

目前，流體內域造波主要有動量源造波和質量源造波。動量源造波是將目標波浪作為人工源項加到動量方程中生成所需波浪場，質量方程保持不變，自由水面以下區域均屬造波源區，造波過程中在整個造波源區內施加動量源方程；質量源造波是將目標波浪作為人工源項加到質量方程中生成所需波浪場，造波過程中并不是自由水面以下區域均屬造波源區，源區頂部與自由水面之間存在一定高度無質量源的水體，正因為此部分水體的存在使得源區高度的確定變得繁瑣。Lin和Liu[3]通過試算得出質量源造波源區中心到靜水面之間的距離為13～12水深時所造波浪比較理想。周玲玲，丁全林，蘭慶琳等[4]指出質量源距離水面13水深時的模擬結果較為精確，并將此處定義為質量源的標準位置。余勇飛[5]在深水數值水槽造波時質量源區上方水體高度采用一倍波高時所造波浪精度較高，與理論波浪匹配較好。李宏偉[6]在造孤立波時質量源區上方水體高度采用半倍波高時，所造孤立波波形與理論波形匹配較好。Chen和Hsiao[7]在FLOW-3D中開發了水體內域造波質量源區域的設計方法，詳細描述了質量源位置的確定公式，為質量源造波方法的應用提供了便利。

已有文獻質量源區上方水體高度大多采用試算法確定，只給出經驗高度區間，并且試算法的成功快慢取決于個人經驗，對于經驗不足者尤顯耗時。筆者造波過程中發現質量源區上方水體高度為一倍波高時在水深大于1.5 m的深水情況下所造波浪極好，但在水深小于1.5 m的淺水中所造波浪并不理想，尤其是在水平網格步長與源區高度值相比小于10時所造波浪波高與目標波浪波高相差較大，必須不斷地加大或減小源區上方水體高度以減小差距。而通常情況下需要造多種波高的目標波浪，這樣每次變化波高時都需要不斷調整源區上方水體高度，既費時也是沒有必要的。故本文針對質量源造波源區上方水體高度難以確定的問題，提出改進方法，避免了多次試算帶來的麻煩，進一步完善了質量源造波方法。

1 數值模型

基于FLUENT求解器，求解雷諾平均N-S方程，采用雷諾應力粘性模型封閉方程組。采用有限體積法離散方程，PISO壓力速度耦合算法求解離散方程，對流項采用三階MUSCL格式，用改進的QUICK方法跟蹤自由表面。邊界條件為左右側為對稱邊界條件，底部為墻邊界條件，上部為壓力進口邊界條件。以下僅就文中的控制方程、造波方法及消波方法進行描述。

1.1控制方程

張量形式的連續性方程和動量方程如下

(1)

(2)

對于二維問題i,j=1,2；ui,uj為速度；p為壓強；ρ為密度；g為重力加速度。粘性應力張量τij為

(3)

μ為分子運動粘性。

1.2造波方法

數值波浪水槽的質量源設置于計算域內，質量源區距離左側消波區至少一倍波長的距離可以很好地消除波浪二次反射(Lin & Liu[1])，本文造波源區距離左側消波區2L(L為波長)。

造波源區內質量守恒方程為

(4)

s(xi,t)為質量源強方程，其大小由設計的入射波波面方程決定，造波源區內非零，造波源區外為零。理論上，當質量源寬度遠小于波長時，可作點源處理，假定質量源引起的質量增減對應于入射波波面的變化，并且源區內質點的水平速度等于設計波浪的水平速度，因為質量源可以看做均勻分布的點源，所以在造波源域Ω內每個網格中應滿足下列關系[5-6，8-9]

s(xi,t)ΔxΔyΔt=2u(xi,t)ΔyΔt

(5)

化簡式(5)得下式

(6)

式中：xi為“源強”的水平位置；Δx為網格的水平尺度；Δy為網格的垂向尺度；右端系數2表示質量源產生傳播方向相反的兩列波。

二維二階Stocks 波的水平速度方程如下

(7)

將式(7)代入(6)得源強方程如下

(8)

式中：H為波高；k為波數；ω為頻率；D為水深。

以下通過對造波源區域源強方程的改進，可以得到更適合淺水區域的造波方法。造波源區簡圖如圖1所示，造波源區底部為墻邊界并無水體向下噴出，考慮單位時間內造波源區其他3個方向噴出的水體質量守恒可推得

s(t)·Δx·Hs=2u·Hs+Δx·v

(9)

式中：u為造波源區流體質點水平速度；v為造波源區流體質點垂向速度；Hs為造波源區高度。化簡(9)可得

s(t)=2uΔx+vHs

(10)

圖1 源造波區簡圖Fig.1 Sketch of the source wave-maker zone

質量源造波假定源區內各點產生的水平速度等于波生成對應點產生的水平速度，造波源區上方水體依靠水體自重形成的恢復力產生波浪，因此源區上方至自由水面高度的確定比較關鍵。仔細研究方程(10)可以發現：造波過程中u和v的量級相當，由于Δx一般小于波長的5%，非常小，而Hs一般大于23倍的水深，對于水深大于1.5 m的深水Δx相比Hs更為小量，至少小一個量級，式中的第二項對源強的貢獻較小可以忽略，但對于水深小于1.5 m的淺水，第二項會隨著水深的減小而增大，其對造波源的貢獻也會逐漸增大，如果第二項繼續忽略會影響所造波浪的波高。數值波浪水槽造波時如果水深小于1.5 m，Hs小于1.0 m，Δx又小于0.1 m，可以看出此時第一項與第二項相比在10以內，不足一個量級，如果水深繼續減小，Hs繼續減小，式中第二項對源強的貢獻會更大，不應忽略。根據Broren.M.,Larsen.J.[10]描述，源造波的原理類似于“水下爆炸”原理，那么源爆炸瞬間質點各個方向的速度大小應該相等。所以質量源向外噴出的流體速度各方向應該大小相等，因此u=v，最終可得源強方程為

(11)

周玲玲，丁全林，蘭慶琳等[4]指出為避免數值振蕩，質量源初始時刻應設置為0，二階Stocks波令初始相位為-π2，并在計算歷程開始時刻設置軟啟動緩慢施加質量源驅動波面變化。本文采用此種做法，最終施加的源強方程如下

圖2 源區右側數值水槽示意圖Fig.2 Right side schematic diagram of numerical wave tank

(12)

式中：N為緩變控制參數；t0為緩變時間。

1.3消浪層

為了吸收二次反射波和不需要的源造波，數值波浪水槽左右端應布置至少一倍波長的海綿吸波層耗散波能。本文采用的吸收方程是由Niels and David[11]提出的橢圓方程(13)，x0表示消波層起始點x坐標；xs表示消波層長度，方程從海綿層起點1逐漸減小到海綿層終點0。

圖3 實際波形與理論波形的對比Fig.3 Comparisons of present wave elevation and theory wave elevation

海綿層消波方程μ(X)

(13)

數值波浪水槽源區右半側示意圖，如圖2所示。

2 數值模擬結果

本文數值波浪水槽長36 m，高0.7 m，水深0.5 m，質量源高度Hs為水深的23，△x小于波長的5%。以下對所造波浪結果、波形沿程變化、波高對目標波浪的影響以及時間步長對波形的影響進行描述。

圖4 波形沿程變化圖Fig.4 Wave elevation displacement of present model with distance from the source zone

2.1實際波形與理論波形對比

為了驗證所造波浪的有效性，給出水深D=0.5 m，波高H=0.06 m，周期T分別為0.8 s、1.0 s、1.2 s、1.4 s時的實造波浪波形與理論波浪波形的對比(圖3)。從圖3中可以看出：不同周期時的實際波形與理論波浪波形匹配較好，雖波谷處振幅值比理論值略小但仍都可以滿足精度要求。

2.2波浪沿程變化情況

為了說明波浪沿程變化情況，給出了水深D=0.5 m，波高H=0.06 m，周期T分別為0.8 s、1.0 s、1.2 s、1.4 s時造波源區(約5倍波長處)右側波浪10倍波長范圍內的波浪沿程變化情況(圖4)。從圖4可以看出：所造波浪沿程衰減較小，滿足目標波浪精度要求。

圖5 所造波浪波高與理論波高的對比Fig.5 Comparisons of the targeted wave height and theoretical wave height

2.3波高對目標波浪的影響

為了檢驗改進的質量源造波方法所造不同波高的波浪與理論波浪的差距，本部分計算了D=0.5 m，T=1.2 s，波高H分別為0.06 m、0.08 m、0.10 m、0.12 m造波源區(約5倍波長處)右側波浪10倍波長范圍內的波浪隨時間的變化(圖5)。從圖5可以看出：隨著波高的增大，所造波浪波高與理論波高的差距在增大；波高H=0.10 m時所造波浪波高只是在波峰處略小于理論值，波谷處仍與理論值擬合較好，但是當波高H=0.12 m時所造波浪波峰和波谷處均與理論值有偏差，所造波浪效果變差。由此可知當波高與水深的比值超過0.2時所造波浪波高與理論值的擬合較差，所以此種造波方法適合于波高與水深的比值在0.2以下的情況。

2.4時間步長的影響

圖6 波形隨時間變化圖Fig.6 Time step effects on wave surface

為了檢驗時間步長對改進的質量源造波方法所造波浪的影響，給出了水深D=0.5 m，周期T=1.2 s，波高H=0.06 m時所造波浪波形隨時間的變化。時間步長驗證之前需要對網格進行反復調整驗證，當網格連續3次調整后的波形無較大變化后停止，進行時間步長驗證。時間步長分別取0.01 s、0.005 s、0.002 s、0.000 2 s進行驗證。對比所得波面由圖6可以看出：當dt=0.01 s、dt=0.005 s時，隨著時間的增長波面出現了幅度較大的上下飄動現象，而當dt=0.002 s、dt=0.000 2 s時，波面已無明顯的上下飄動現象，波面較穩定。

3 結論

將質量源區高度值考慮到造波源強方程之中，量化源區上方水體高度，克服了水深小于1.5 m時造波過程中源區上方水體高度難以確定的困難，得到了適合于淺水的質量源造波方法。該方法可以在水深小于1.5 m的淺水中造出比較理想的波浪，更適合于水平網格步長相比于源區高度值小10倍的情況。文中驗證了波高對所造波浪的影響、時間步長對淺水中所造波浪的影響及波浪沿程變化，得出結論：當波高與水深的比值小于0.2時所造波浪波高與理論波高值的擬合較好；當時間步長小于0.005 s時所造波浪無明顯的上下漂浮現象，波形穩定，符合理論波形要求，波浪沿程衰減較小。

[1]Lin P, Liu L F. Discussion of "Vertical variation of the flow across the surf zone" [Coast. Eng. 45 (2002) 169-198][J]. Coastal Engineering, 2004, 50(3):161-164.

[2]Larsen J, Dancy H. Open boundaries in short wave simulations-A new approach[J]. Coastal Engineering, 1983, 7(3):285-297.

[3]Lin P, Liu L F. Internal Wave-maker for Navier-stokes Equations Models[J].Journal of Waterway, Port, Coastal and Ocean Engineering, 1999,125(4): 207-215.

[4]周玲玲,丁全林,蘭慶琳,等.基于LB方法的質量源造波的數值波浪水槽[J].水電能源科學，2015,33(3):99-103. ZHOU L L,DING Q L,LAN Q L, et al.Numerical Wave Tank with Mass Source Generating Waves Based on LB Method[J]. China Civil Engineering Journal, 2015, 33(3): 99-103.

[5]余勇飛.波浪的數值模擬及其流場結構分析[D].哈爾濱：哈爾濱工業大學，2013.

[6]李宏偉.造波理論與方法研究[D].哈爾濱：哈爾濱工程大學，2013.

[7]Chen Y L, Hsiao S C. Generation of 3D water waves using mass source wavemaker applied to Navier-Stokes model[J]. Coastal Engineering, 2016, 109:76-95.

[8]房卓,張寧川,臧志鵬.透空式梳式防波堤的數值模擬和波浪透射系數的研究[J].水道港口,2011,32(2)：86-93. FANG Z, ZHANG N C, ZANG Z P. Numerical simulations of open comb-type breakwater and research on its wave transmission coefficient[J].Journal of Waterway and Harbor,2011,32(2):86-93.

[9]FANG Z, ZHANG N C, ZANG Z P. Experimental and Numerical Study on Hydrodynamic Performance of Impermeable Comb-Type Breakwater[J]. Journal of Ship Mechanics, 2012, 16(6):632-645.

[10]Brorsen M, Larsen J. Source generation of nonlinear gravity waves with the boundary integral equation method[J]. Coastal Engineering, 1987, 11(2):93-113.

[11]Jacobsen N G, Fuhrman D R, Fredsφe J. A wave generation toolbox for the open-source CFD library: OpenFoam[J]. International Journal for Numerical Methods in Fluids, 2012, 70(9):1 073-1 088.

The method of mass source wavemaker in shallow water

TIANZheng-lin,SUNZhao-chen,LIANGShu-xiu

(StateKeyLaboratoryofCoastalandOffshoreEngineering,DalianUniversityofTechnology,Dalian116024,China)

In view of the problem that the distance between the top of the source region and the still water level is difficult to determine in the process of mass source wavemaker, the improved method was proposed in this paper. In the mass source functions took the height of the source region into account, the inconvenience of many trials to get the height of the mass source region can be eliminated in the process of wavemaker. This method remedies the defect that the difference between targeted wave and theoretical wave in the shallow water whose depth less than 1.5 m, especially horizontal grid step is an order of magnitude 10 smaller than the source region height, the mass source wave-maker method is promoted in the shallow water.The effects of wave height and the time step on simulated waves were studied. The conclusions show that the fitting between made wave and theoretical wave is good when the ratio of wave height and water depth is less than 0.2;when the time step is less than 0.005 s, the floating phenomenon of the simulated waves is well controlled, waves surface is stable and conforms to the requirements of the wave theory.

mass source wavemaker; shallow-water wavemaker; height of source region; improved method

TV 139.2

：A

：1005-8443(2017)04-0325-05

2017-01-13；

：2017-04-06

國家自然科學基金資助項目(51279028)；國家海洋局公益性行業科研專項(201405025-1)

田正林(1985-),男，黑龍江省哈爾濱人，博士研究生，主要從事波浪與結構物相互作用的研究。

Biography： TIAN Zheng-lin(1985-),male, doctor student.