黃瑞啟，黃筱云,2，張少?gòu)?qiáng)，程永舟,2

(1.長(zhǎng)沙理工大學(xué)水利與環(huán)境工程學(xué)院,湖南 長(zhǎng)沙 410114; 2.水沙科學(xué)與水災(zāi)害防治湖南省重點(diǎn)實(shí)驗(yàn)室,湖南 長(zhǎng)沙 410114)

庫(kù)岸邊坡失穩(wěn)、崩落，滑坡體以一定速度沖擊水體，將產(chǎn)生巨大的涌浪。國(guó)內(nèi)外滑坡浪事件屢見(jiàn)不鮮，如1963年意大利的瓦伊昂水庫(kù)滑坡[1]、2007年清江水布埡水庫(kù)大堰塘滑坡[2]、2008年長(zhǎng)江巫峽龔家方滑坡[3]等。按照滑坡發(fā)生位置，可以將滑坡分為水上滑坡(subaerial)、水下滑坡(submarine)和濱水滑坡(partially submerged)3種類(lèi)型[4]。與濱水滑坡、水下滑坡相比，水上滑坡會(huì)沖擊水面，涌浪的初始形態(tài)更加復(fù)雜，且與滑坡體運(yùn)動(dòng)與動(dòng)力特征密切相關(guān)[5]。水上滑坡產(chǎn)生涌浪過(guò)程可分為以下4個(gè)階段[6]：①滑坡體開(kāi)始移動(dòng)、加速、沖擊水面；②滑坡體入水滑行產(chǎn)生涌浪；③涌浪向外傳播；④涌浪抵達(dá)岸線(xiàn)、爬高。按近場(chǎng)首波的非線(xiàn)性特征，滑坡涌浪大致可分為線(xiàn)性波、孤立波、橢圓余弦波和涌潮[7]，具體類(lèi)型與滑坡體體積以及滑坡體與水之間的動(dòng)量交換有關(guān)。

滑坡涌浪的特征研究主要采用理論分析[8-12]、物理模型試驗(yàn)[13-19]和數(shù)值模擬[20-26]方法?，F(xiàn)階段滑坡涌浪的數(shù)值模擬方法有基于N-S方程的數(shù)值模擬[27-28]、基于淺水方程的數(shù)值模擬[29-30]、基于Boussinesq波方程的數(shù)值模擬[31-32]和基于勢(shì)流方程的數(shù)值模擬[33-34]4類(lèi)。若要完全模擬滑坡體沖擊水面造成水花飛濺和涌浪的傳播過(guò)程，基于N-S方程的數(shù)值模型是唯一選擇?；贜-S方程的數(shù)值模型又分為二維和三維兩類(lèi)。在二維模型中，滑坡體寬度與水體寬度一致，涌浪只在一個(gè)方向上傳播；而三維模型可以呈現(xiàn)涌浪的不同徑向傳播?；掠坷说娜S數(shù)值模擬研究遠(yuǎn)少于二維數(shù)值模擬研究，且滑坡涌浪傳播遠(yuǎn)場(chǎng)特征的分析與討論較少見(jiàn)。

本文采用FLOW-3D軟件模擬水上滑坡產(chǎn)生涌浪的三維過(guò)程，分析不同徑向上周期與波高變化過(guò)程以及涌浪傳播近遠(yuǎn)場(chǎng)的界限，并探討了遠(yuǎn)場(chǎng)涌浪能量大小及其能量轉(zhuǎn)換率，以期為滑坡涌浪地質(zhì)災(zāi)害的預(yù)警提供參考。

1 數(shù)值模型及驗(yàn)證

1.1 模型概況

三維滑坡涌浪模型(以下簡(jiǎn)稱(chēng)“本文模型”)采用FLOW-3D軟件建模，其中流體控制方程為不可壓縮流體運(yùn)動(dòng)的Navier-Stokes方程，采用標(biāo)準(zhǔn)RNGk-ε模型描述流體紊動(dòng)特征，GMO模型控制滑坡體運(yùn)動(dòng)，TruVOF技術(shù)捕捉水面起伏和飛濺。通過(guò)調(diào)整GMO模型中的摩擦系數(shù)和碰撞系數(shù)可以改變滑坡體下滑、撞擊水面、觸底等狀態(tài)，以獲得準(zhǔn)確的涌浪結(jié)果，碰撞系數(shù)取1表示碰撞物完全彈性，取0則是完全塑性。

1.2 模型驗(yàn)證

通過(guò)與Liu等[35]的物理試驗(yàn)結(jié)果對(duì)比來(lái)檢驗(yàn)本文模型的有效性。該物理試驗(yàn)在長(zhǎng)104 m、寬3.7 m、高4.6 m的水槽中進(jìn)行，槽內(nèi)水深為3 m，滑坡坡度為27°，滑坡體長(zhǎng)91 cm、寬65 cm、高46 cm，質(zhì)量為475.52 kg，初始時(shí)刻滑坡體頂面距水面高度為0.46 m。本文模型采用漸變網(wǎng)格，最小網(wǎng)格尺寸是0.02 m，時(shí)間步長(zhǎng)根據(jù)柯朗條件自動(dòng)調(diào)整，摩擦系數(shù)和碰撞系數(shù)均取0.1。

圖1為6個(gè)測(cè)點(diǎn)(圖2)波面高度η實(shí)測(cè)結(jié)果與模擬結(jié)果的比較?？傮w上，本文模型模擬結(jié)果較Liu等[35]模擬結(jié)果更接近實(shí)測(cè)結(jié)果，說(shuō)明本文模型各參數(shù)設(shè)置合理。

圖1 本文模型模擬結(jié)果與Liu等[35]模擬和試驗(yàn)結(jié)果比較

圖2 Liu等[35]物理試驗(yàn)測(cè)點(diǎn)布置(單位：cm)

2 數(shù)值模擬試驗(yàn)

2.1 試驗(yàn)設(shè)置

整個(gè)試驗(yàn)區(qū)域縱向長(zhǎng)度設(shè)為200 m，水面橫向?qū)挾葹?00 m，水深為10 m，滑坡坡度為45°。滑坡體為前緣角度為45°的四棱柱，長(zhǎng)9 m、寬8 m、厚3 m，置于區(qū)域中線(xiàn)上，如圖3所示。為節(jié)省計(jì)算資源，將整個(gè)區(qū)域按中線(xiàn)分開(kāi)，并將與中線(xiàn)重合的邊界設(shè)為對(duì)稱(chēng)邊界，其余邊界設(shè)置為出流邊界，并布置多孔介質(zhì)用于消波。

圖3 三維滑坡涌浪數(shù)值模擬試驗(yàn)測(cè)點(diǎn)布置(單位：m)

2.2 數(shù)值模擬結(jié)果

滑坡涌浪的產(chǎn)生與傳播區(qū)域可以分為飛濺區(qū)、近場(chǎng)和遠(yuǎn)場(chǎng)。一般而言，入水點(diǎn)至滑坡體停止點(diǎn)的水域?qū)儆陲w濺區(qū)；從停止點(diǎn)開(kāi)始，波浪出現(xiàn)色散的區(qū)域稱(chēng)為近場(chǎng)；呈現(xiàn)振蕩波特點(diǎn)的區(qū)域則為遠(yuǎn)場(chǎng)[5]。就本文數(shù)值模擬試驗(yàn)而言，坡腳線(xiàn)可以看作飛濺區(qū)與近場(chǎng)的界限。為獲取滑坡涌浪近遠(yuǎn)場(chǎng)涌浪三維特征，在0°徑向上等間距布置6個(gè)測(cè)點(diǎn)，30°徑向上布置3個(gè)測(cè)點(diǎn)，60°徑向上布置2個(gè)測(cè)點(diǎn)，如圖3所示，其中，2、8號(hào)測(cè)點(diǎn)和10號(hào)測(cè)點(diǎn)徑距一致，3、9號(hào)測(cè)點(diǎn)和11號(hào)測(cè)點(diǎn)徑距相等。表1為滑坡體以4 m/s初始速度下滑時(shí)各測(cè)點(diǎn)的首、次波波峰高度ηc與周期T模擬結(jié)果。

表1 各測(cè)點(diǎn)首、次波波峰高度與周期

圖4 0°徑向各測(cè)點(diǎn)波面時(shí)程變化

圖5 坡腳各測(cè)點(diǎn)波面時(shí)程變化

圖6 不同徑向等徑距測(cè)點(diǎn)波面時(shí)程變化

3 滑坡涌浪特征分析

為分析滑坡涌浪特征，本文通過(guò)改變滑坡體厚度s、寬度b與入水速度v以獲得不同的涌浪形態(tài)，數(shù)值模擬試驗(yàn)具體工況見(jiàn)表2。

表2 數(shù)值模擬試驗(yàn)工況與結(jié)果

3.1 滑坡體運(yùn)動(dòng)

(1)

式中l(wèi)為滑坡體長(zhǎng)度。由于滑坡體下滑主要受到水體的正面阻力，若不考慮側(cè)向阻力，有

(2)

根據(jù)27組工況的滑坡體滑坡特征時(shí)間數(shù)值模擬結(jié)果(表2)，相對(duì)滑坡特征時(shí)間的擬合公式為

(3)

式(3)的擬合度為0.991，均方誤差為3.47×10-4。

3.2 近遠(yuǎn)場(chǎng)

Walder等[5]將近場(chǎng)定義為從飛濺區(qū)邊界至涌浪動(dòng)能與勢(shì)能趨近于某一漸近值的位置，Watts[13]則采用如下能量積分表示通過(guò)位置x斷面的涌浪總勢(shì)能：

(4)

當(dāng)上述積分沿程趨近一定值時(shí)，便可確定近場(chǎng)的外邊界。圖7～9為不同徑向上沿程各測(cè)點(diǎn)的涌浪總勢(shì)能(圖中s*=s/h為滑坡體相對(duì)厚度，b*=b/h為滑坡體相對(duì)寬度)，橫坐標(biāo)x*=x/h?？梢钥闯觯倓?shì)能沿程減小，其中，弗勞德數(shù)對(duì)近場(chǎng)總勢(shì)能的沿程變化影響較小，而相對(duì)厚度與相對(duì)寬度的影響明顯。在1～3倍水深之間，總勢(shì)能的下降幅度較大，而后下降幅度開(kāi)始減少，當(dāng)涌浪達(dá)到6～7倍水深區(qū)間時(shí)，各徑向涌浪總勢(shì)能變化較小，且在7倍水深位置不同徑向的涌浪總勢(shì)能幾乎相等。因此，7倍水深距離可看作三維滑坡涌浪近場(chǎng)的外邊界。

圖7 不同滑坡體弗勞德數(shù)下各徑向涌浪能量積分沿程變化

圖8 不同滑坡體相對(duì)厚度下各徑向涌浪能量積分沿程變化

圖9 不同滑坡體相對(duì)寬度下各徑向涌浪能量積分沿程變化

3.3 波峰高度

(5)

(6)

(7)

(8)

在遠(yuǎn)場(chǎng)，涌浪幾乎以固定頻率傳播，遠(yuǎn)場(chǎng)涌浪最大波幅將出現(xiàn)在近遠(yuǎn)場(chǎng)交界處。取徑距為7倍水深位置作為遠(yuǎn)場(chǎng)起始位置，最大相對(duì)波峰高度出現(xiàn)在次波，其試驗(yàn)值與預(yù)測(cè)值關(guān)系如圖10所示，其擬合公式為

圖10 近遠(yuǎn)場(chǎng)相對(duì)波峰高度試驗(yàn)值與預(yù)測(cè)值關(guān)系

(9)

式(9)的擬合度為0.984 5，均方誤差為7×10-7。

3.4 周期

圖11 遠(yuǎn)場(chǎng)次波相對(duì)周期試驗(yàn)值與預(yù)測(cè)值的關(guān)系

(10)

式(10)擬合度為0.85，均方誤差為7×10-6。

3.5 波形

遠(yuǎn)場(chǎng)涌浪的相對(duì)波長(zhǎng)為

(11)

根據(jù)式(9)(10)(11)，推導(dǎo)出遠(yuǎn)場(chǎng)涌浪非線(xiàn)性Ursell數(shù)公式為

(12)

式中波高H近似為波幅的2倍。27組工況數(shù)值模擬試驗(yàn)得到的Ursell數(shù)如圖12所示,不難發(fā)現(xiàn)，遠(yuǎn)場(chǎng)涌浪均屬于線(xiàn)性波。

圖12 各徑向遠(yuǎn)場(chǎng)涌浪的Ursell數(shù)

3.6 遠(yuǎn)場(chǎng)涌浪能量

計(jì)算遠(yuǎn)場(chǎng)涌浪能量對(duì)滑坡涌浪災(zāi)害的預(yù)測(cè)評(píng)估具有重要的影響。通過(guò)徑距為R處的半圓斷面的涌浪總勢(shì)能Ep(R)為

(13)

根據(jù)線(xiàn)性波理論，遠(yuǎn)場(chǎng)涌浪的勢(shì)能與動(dòng)能相同，故通過(guò)遠(yuǎn)場(chǎng)斷面(R=7h)的涌浪總能量Ew=2Ep，27組工況的近遠(yuǎn)場(chǎng)分界斷面涌浪總能量見(jiàn)表2。涌浪的總能量來(lái)自于滑塊的能量，包括動(dòng)能和勢(shì)能，故定義能量傳遞率ε為遠(yuǎn)場(chǎng)涌浪能量與滑坡體能量的比值(ε=Ew/Es，其中Es為滑坡體的能量)。計(jì)算結(jié)果表明，27組工況的能量傳遞率在8.0%～19.7%之間。能量傳遞率主要與波峰高度時(shí)程變化有關(guān)，根據(jù)數(shù)值模擬結(jié)果，擬合的能量傳遞率公式為

ε=0.295 1(t*/A*)-0.503 4

(14)

式(14)的擬合度為0.860 1，均方誤差為2.08×10-4。能量傳遞率的試驗(yàn)值與預(yù)測(cè)值關(guān)系見(jiàn)圖13。

圖13 能量傳遞率的預(yù)測(cè)值與試驗(yàn)值的關(guān)系

4 結(jié) 論

a.三維滑坡涌浪在不同徑向上均出現(xiàn)衰減，但在7倍水深距離后保持穩(wěn)定，該距離外的波峰高度基本與滑坡體弗勞德數(shù)、厚度以及寬度無(wú)關(guān)，故可認(rèn)為沖擊滑坡涌浪的近遠(yuǎn)場(chǎng)界線(xiàn)位于7倍水深處。

b.相對(duì)滑坡特征時(shí)間是滑坡體弗勞德數(shù)、斜坡角度以及相對(duì)滑坡體尺寸的函數(shù)，決定了近遠(yuǎn)場(chǎng)最大波峰高度和遠(yuǎn)場(chǎng)次波周期的大小，也對(duì)遠(yuǎn)場(chǎng)次波的非線(xiàn)性特征有一定影響。

c.定義了滑坡體的能量傳遞率，本文工況中能量傳遞率在8.0%～19.7%之間。