李 杰

(西南交通大學(xué)土木工程學(xué)院，四川成都610036)

隨著水波動力學(xué)和計算機科學(xué)的發(fā)展，數(shù)值模擬由于其易修改、避免了物理實驗的縮比效應(yīng)等特點，成為了國內(nèi)外研究的重點方向。數(shù)值波浪水槽則成為了模擬波浪運動，研究其與結(jié)構(gòu)物之間相互作用的重要工具。其關(guān)鍵在于造波和消波方法的技術(shù)實現(xiàn)。目前較為常用的數(shù)值造波方法包括動量源項造波法、動邊界造波法(推板式或搖板式)、給定入射邊界速度法及造波區(qū)域法；常用消波方法包括動量源項消波法和多孔介質(zhì)消波法。

M Brorsen和 JLarsen［1］根據(jù)源造波理論，通過在計算域內(nèi)設(shè)置的造波源可分別向源區(qū)域兩邊同時生成兩列方向相反的波。Peter Troch和 Julien De Rouck［2］應(yīng)用VOF追蹤數(shù)值波浪水槽的自由表面，首次實現(xiàn)了推板造波且消除邊界反射波浪的技術(shù)，減少水槽長度從而提高了計算效率。劉加海等［3］通過仿照物理水槽造波方法，利用Fluent軟件模擬了二維規(guī)則波浪。嚴(yán)汝建等［4］利用Fluent軟件，基于搖板造波結(jié)合阻尼層消波方法，實現(xiàn)了二維規(guī)則波與不規(guī)則波的模擬。黃華等［5］通過模擬活塞推板運動造波，設(shè)置人工阻尼區(qū)消波，對二階Stokes波進行了仿真。王本龍等［6］推導(dǎo)了在變地形情況下適用的高階Boussinesq波浪模型，并數(shù)值模擬分析了波流相互作用問題。董志等［7］分析了以上幾種造波-消波方法組合的優(yōu)劣性，提出了優(yōu)化的造波-消波配置。文獻［8］中指出，通常在入射邊界處，難以給定理論速度或波高與實際波浪非線性相匹配，而模擬造波板運動能解決此問題。盤俊等［9］測試分析了武漢理工大學(xué)某小型推板式波浪水槽的造波性能。并指出物理波浪水槽造波的推板式適用于淺水造波，而搖板式適用于較深水造波。

因此本文針對梅沃特(Le Mehaute B，1969)提出的波浪理論適用范圍中中等水深區(qū)域波浪，基于推板式動邊界造波法結(jié)合動量源消波法建立數(shù)值波浪水槽進行造波、消波性能測試，并與理論值進行對比。了解中等水深波浪特性差異，掌握造波性能標(biāo)準(zhǔn)尺度，得到推板造波方法模擬此類波浪的適用性量化依據(jù)。

1 數(shù)值水槽及波浪參數(shù)

本文利用FLUENT軟件UDF功能，以不可壓縮黏性流體的N-S方程為控制方程，通過 VOF法捕捉波浪自由面，水槽上部邊界設(shè)為pressure-outlet，兩側(cè)邊界設(shè)為symmetry，其余邊界均設(shè)為wall，采用壓力速度耦合PISO算法求解；設(shè)置動網(wǎng)格以模擬推板運動，并在水槽末端設(shè)置消波段。

1.1 數(shù)值水槽尺寸及網(wǎng)格劃分

水槽規(guī)格1：長(x)×寬(z)×高(y)為(500×3×8)m，水深d為6 m；由于所選用橢圓余弦波波長較長，為得到沿程足夠多的波數(shù)，水槽規(guī)格2：長(x)×寬(z)×高(y)為(1000×3×8)m，水深d為6 m。

在y方向靠近液面處、x方向入口動邊界處進行網(wǎng)格加密(圖1)，消波段區(qū)域進行網(wǎng)格放疏以減少網(wǎng)格數(shù)量(圖2)。

圖1 水槽入口端網(wǎng)格加密

圖2 水槽末端消波段網(wǎng)格放疏

1.2 波浪參數(shù)設(shè)置

本文工況設(shè)置以無因次參數(shù)d/gT2和H/gT2的變化為參考，選用波浪參數(shù)如表1所示。

2 數(shù)值穩(wěn)定判據(jù)及消波效果檢驗

2.1 數(shù)值穩(wěn)定判據(jù)

以波高H=0.4 m，周期T=4 s的StokesⅡ波舉例，x=50 m處的波面時程曲線如圖3所示，當(dāng)波面高程穩(wěn)定時，說明此時數(shù)值結(jié)果已達到穩(wěn)定狀態(tài)。

表1 波浪參數(shù)

圖3 波面時程曲線

2.2 消波效果檢驗

當(dāng)數(shù)值結(jié)果已經(jīng)穩(wěn)定，為排除反射波對入射波的干擾，需對消波效果進行檢驗。取穩(wěn)定后某時刻(t=190 s)波面曲線，如圖4所示，在距離造波入射口x約為475 m處，波高已為0 m，達到消波要求。

圖4 消波效果

3 數(shù)值波浪水槽性能分析

參考文獻［10］中對造波設(shè)備產(chǎn)生波浪的穩(wěn)定性及重復(fù)性要求，由于本文采用數(shù)值模擬造波，不受自然環(huán)境影響，省略重復(fù)性檢驗。故對周期、波長的誤差及穩(wěn)定性，波高的誤差、穩(wěn)定性及沿程衰減情況進行分析。

3.1 周期誤差及穩(wěn)定性分析

為將數(shù)值模擬產(chǎn)生波浪(模擬波)的周期與相應(yīng)理論周期值定量對比，引入周期相對誤差值［11］ErT，其計算公式見式(1)，ErT值越小表示誤差越小。其中T0表示周期理論值，T′表示模擬波周期平均值(圖5～圖7)。

圖5 StokesⅡx=50m H=0.4m

圖6 Airy x=50m H=0.04m

圖7 Cnoidal x=50m H=0.6m

周期穩(wěn)定性誤差值［10］RT，其計算公式見式(2)，RT值越小表示模擬波周期越穩(wěn)定。其中T′表示模擬波平均周期，Tmax表示最大周期，Tmin表示最小周期。

各模擬波穩(wěn)定后的時程曲線均接近重合于理論值，周期相對誤差值ErT及周期穩(wěn)定性誤差值RT均極接近于0。

3.2 波長誤差及穩(wěn)定性分析

圖8～圖10所示為三類模擬波浪穩(wěn)定后，某時刻波面沿程分布與理論值的對比，為分析波長的誤差及穩(wěn)定性，引入波長相對誤差值ErL、波長穩(wěn)定性誤差值RL，其計算方法同ErT、RT。波長誤差指標(biāo)結(jié)果見表2。

圖8 StokesⅡH=0.4m

圖9 Airy H=0.04m

圖10 Cnoidal H=0.6m

表2 穩(wěn)定后波長誤差指標(biāo)

3.3 波高誤差

對比波浪相關(guān)參數(shù)模擬值與理論值，差異最顯著的是波高，參考文獻［12］引入以下誤差指標(biāo)評價波高誤差。

(1)波高穩(wěn)定性誤差值［10］RH，其計算方法同上。

(2)峰值誤差［12］ErP表示模擬波穩(wěn)定后的某時刻波面升高值Y’(i)的峰值與谷值相對于理論值Y0(i)的峰值與谷值的相對誤差求和平均值，其計算公式見式(3)。

(3)波峰的平均誤差［12］EF表示模擬波穩(wěn)定后的波峰值W′F(i)與理論波峰值WF之差的求和平均值與理論波峰值WF的比值，其計算公式見式(4)。

(4)波谷的平均誤差［12］EG表示模擬波穩(wěn)定后的波谷值W′G(i)與理論波谷值WG之差的求和平均值與理論波谷值WG的比值，其計算公式見式(5)。

(5)波高沿程衰減度［9］DH表示單位距離相對波高沿程衰減程度，其計算公式見式(6)。H1、H2分別表示不同測點所監(jiān)測平均波高。

穩(wěn)定后波高誤差指標(biāo)結(jié)果見表3。

橢圓余弦波模擬波，其波高沿程未出現(xiàn)衰減，取數(shù)值水槽前段(即為除去消波段)，紅線表示理論波高峰值連線位置(圖11)。

圖11 Cnoidal H=0.6m

4 結(jié)論

(1)本文基于推板造波數(shù)值方法結(jié)合動量源消波方法建立的數(shù)值水槽消波效果良好。

(2)在周期方面，模擬波與理論波對比，幾乎無誤差、周期穩(wěn)定性好。

(3)在波長方面，Airy模擬波誤差最小且波長最穩(wěn)定。由于非線性影響，波長穩(wěn)定性StokesⅡ模擬波＜Cnoidal模擬波＜Airy模擬波，即豎向整體呈現(xiàn)穩(wěn)定性隨H/gT2增大而降低趨勢。且在水深較大區(qū)域，非線性影響顯著，即對于StokesⅡ波，其波高越大，波長越不穩(wěn)定。

表3 穩(wěn)定后波高誤差指標(biāo)

(4)在波高方面，Cnoidal模擬波誤差最小。由于非線性影響，波高穩(wěn)定性StokesⅡ模擬波＜Cnoidal模擬波＜Airy模擬波，即豎向整體呈現(xiàn)穩(wěn)定性隨H/gT2增大而降低趨勢。且在水深較大區(qū)域，非線性影響顯著，即對于StokesⅡ波，其波高越大，波高越不穩(wěn)定。

(5)沿程波高衰減程度Cnoidal模擬波＜Airy模擬波＜StokesⅡ模擬波，即橫向整體呈現(xiàn)波高衰減度隨d/gT2增大而增大趨勢。且在水較深區(qū)域，非線性影響顯著，即對于StokesⅡ波，其波高越大，衰減越嚴(yán)重。

(6)在峰值方面，Cnoidal模擬波衰減弱，但“毛刺”現(xiàn)象突出。