魏 凱，姜沫臣，洪 杰

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

隨著我國橋梁事業(yè)的不斷發(fā)展，橋梁建設(shè)不斷向沿海地區(qū)拓展，如杭州灣大橋、東海大橋、港珠澳大橋等。波浪作用是跨海橋梁結(jié)構(gòu)區(qū)別于陸地橋梁的主要荷載之一[1]。波浪由外海向跨海橋梁所在的近岸海域傳播時，隨著水深急劇變淺，波浪波長變短、波陡增大，能量迅速向波峰附近很小的區(qū)域集中。當(dāng)波峰水質(zhì)點(diǎn)的水平運(yùn)動速度超過波浪傳播速度時，波浪會發(fā)生卷曲，并在極短時間內(nèi)破碎，即發(fā)生卷破。一般情況下，各種波浪破碎形式中卷破波對結(jié)構(gòu)產(chǎn)生作用力最大，危害最深，這里僅研究這一形式，因此后文涉及到的破碎波特指卷破波。一個成長完全的破碎波會給結(jié)構(gòu)帶來巨大的沖擊力，其荷載峰值可能達(dá)到波高相同的非破碎波浪的5倍以上[2]，不僅影響海上結(jié)構(gòu)物的動力荷載分布，其瞬態(tài)高頻成分還可能激發(fā)結(jié)構(gòu)振動[3]，放大結(jié)構(gòu)的動力響應(yīng)[4]，威脅橋梁結(jié)構(gòu)安全[5]。研究橋墩受破碎波作用規(guī)律對跨海橋梁建造維護(hù)及橋梁發(fā)展具有十分重要的意義。

國內(nèi)外學(xué)者對海濱圓柱結(jié)構(gòu)的波浪作用開展了大量研究。Wienke等[6]開展了大比例尺圓柱破碎波浪作用荷載試驗(yàn)，討論了圓柱和波浪破碎點(diǎn)相對位置對破碎波浪荷載的影響。Chella等[7]采用淺化破碎波浪沖擊直立圓柱的試驗(yàn)方法，研究了破碎波特性及幾何特征與水深、遠(yuǎn)岸波陡、斜坡坡度以及入射波高的關(guān)系，探究了不同周期下破碎位置波速波壓分布及自由液面形狀的特性。魏凱等[8]開展了極端波浪沖擊高樁承臺尺縮模型的水槽試驗(yàn)，建立了波浪沖擊荷載時程模型，給出了沖擊荷載峰值和沖擊上升段持續(xù)時間的邊緣分布以及基于Copula的聯(lián)合概率分布形式。柳淑學(xué)等[9]通過波浪聚焦的方法在水池中產(chǎn)生多向聚焦波和破碎波，討論了破碎波的波面特性、破碎波的破碎指標(biāo)及傳播過程中的頻譜變化。

隨著數(shù)值模擬技術(shù)的發(fā)展，數(shù)值模擬方法被廣泛應(yīng)用于破碎波。Bihs等[10]利用REEF3D研究了順流向并列排布兩柱體受破碎波浪力大小與其間距及破碎位置的關(guān)系。常爽等[11]利用聚焦波浪理論，建立了生成瞬態(tài)聚焦破碎波浪的數(shù)值水槽。王修亭等[12]基于OpenFOAM模擬近岸單樁結(jié)構(gòu)物所受破碎波浪荷載。高學(xué)平等[13]采用改進(jìn)的標(biāo)記單元法對處于非完全繞流區(qū)的柱上波浪力機(jī)制進(jìn)行了分析和討論。王紅川等[14]建立了拋物型緩坡方程，描述波浪在破波帶內(nèi)能量衰減和波浪破碎系數(shù)之間的關(guān)系，并將數(shù)值與試驗(yàn)進(jìn)行對比驗(yàn)證。

Oumeraci等[15]研究發(fā)現(xiàn)，結(jié)構(gòu)物形狀特征會對其所受破波力產(chǎn)生一定的影響。而跨海橋梁橋墩多為圓端形，且因?yàn)闃蚨蘸奢d不同，截面長寬比差異較大。但現(xiàn)有研究多限于圓柱、直防波堤等，對不同長寬比圓端形橋墩所受波浪破碎作用的特性尚無清晰認(rèn)識。因此，通過建立數(shù)值水槽，采用波浪淺化破碎的方式，考慮k-ω湍流模型，計算了破碎波對圓端形橋墩的破波力，研究了破波力與橋墩長寬比的關(guān)系及破波力組成成分準(zhǔn)靜態(tài)力和沖擊力的變化規(guī)律。

1 數(shù)值模擬

1.1 控制方程

數(shù)值模擬流體運(yùn)動時，假設(shè)流體不可壓縮，通過雷諾平均N-S方程描述流體運(yùn)動，得到流體連續(xù)性和動量守恒控制方程為：

(1)

(2)

式中：i=1，2，3時，xi分別表示x，y，z坐標(biāo)；ui表示流場各方向的時均速度；Ai表示流體所占的面積分?jǐn)?shù)；VF表示流體的體積分?jǐn)?shù)；t表示時間；p表示流體壓強(qiáng)；ρ表示流體密度；gi表示各方向流體加速度。其中，fi表示各方向黏滯力引起的加速度，可表示為：

(3)

(4)

式中：τb,i表示流體在x，y，z方向面上的剪切應(yīng)力；Sij表示應(yīng)變率張量，ν表示動力黏度，νT表示湍流黏度，由湍流模式理論計算得到。

1.2 湍流模型

雷諾平均運(yùn)動方程雖很好地描述了湍流運(yùn)動，但方程中增加了未知的雷諾應(yīng)力項(xiàng)，從而形成基本方程的不封閉問題，為此需引入湍流模型尋找附加條件，選用k-ω湍流模型[7]，湍流黏度根據(jù)式(5)進(jìn)行計算：

νT=k/ω

(5)

其中，k為湍動能；ω為單位湍動能耗散率，由下式計算得到：

(6)

(7)

其中，i=1，2，3時，xi分別表示x，y，z坐標(biāo)；VF表示流體的體積分?jǐn)?shù)；t表示時間，Ai表示流體所占的面積分?jǐn)?shù)；σk、β*、σω、α、β為閉合系數(shù)，σk=2.0，β*=9/100，σω=2.0，α=5/9，β=3/40。

1.3 邊界條件

數(shù)值水槽采用前端造波邊界生成斯托克斯五階波，水槽邊壁采用對稱邊界，邊界上梯度為0，不會對波浪進(jìn)行反射。底部采用墻邊界，上部采用壓力邊界模擬大氣，壓強(qiáng)為101 325 Pa。水槽末端消波段取為6 m，尾部采用出流邊界，使用Sommerfeld輻射條件動態(tài)估計邊界條件，在保持水槽水量不變的條件下允許波浪出流。即認(rèn)為在邊界上：

(8)

其中，φ為要輻射的變量，C為波浪的傳播速度，n為輻射邊界的法向向量。同時，在水槽末端加設(shè)一定長度的消波段，進(jìn)一步消除波浪在出流邊界的反射。采用Flow-3D軟件建立數(shù)值水槽模型，采用立方體網(wǎng)格，輸入網(wǎng)格尺寸dx后，自動在三維水槽范圍內(nèi)生成均一化網(wǎng)格。

1.4 模型設(shè)置

在開展橋墩破碎波浪作用模擬之前，分別以前人試驗(yàn)和數(shù)值模擬為參考，驗(yàn)證了文中數(shù)值水槽在破碎波浪作用下結(jié)構(gòu)受力模擬方面的準(zhǔn)確性。

以Irschik等[16]開展的破碎波沖擊圓柱結(jié)構(gòu)試驗(yàn)和Bihs等[10]基于該試驗(yàn)開展的數(shù)值模擬為參照，開展對比研究。試驗(yàn)采用的水槽寬5 m、高7 m，試驗(yàn)段斜坡坡度為1∶10，試驗(yàn)水深3.8 m，試驗(yàn)中結(jié)構(gòu)模型為直徑0.7 m的圓柱，入射波為規(guī)則波，波高為1.3 m，周期為4 s，波浪由深向淺傳播。

根據(jù)上述試驗(yàn)建立如圖1所示的三維數(shù)值水槽，通過CFD數(shù)值模擬，分析了破碎位置和結(jié)構(gòu)受力。圖2表示波浪破碎的過程，分別表示將要破碎、正在破碎和已破碎的時刻，圖2(b)描述了波浪在墩前破碎時浪舌頂部與靜水面之間的距離——波浪表面高程ηb。在橋墩前端豎向位置z=3～5 m內(nèi)每隔0.1 m布置一測點(diǎn)，設(shè)置如圖3所示測點(diǎn)群α，負(fù)責(zé)監(jiān)控橋墩前側(cè)的壓強(qiáng)分布情況。

圖1 破碎波沖擊橋墩CFD數(shù)值模型Fig. 1 Numerical model of pier impacted by breaking wave

圖2 波浪破碎過程Fig. 2 Wave breaking process

圖3 圓端形橋墩測點(diǎn)群αFig. 3 Monitor point group α of round-ended pier

1.5 參數(shù)率定

為研究網(wǎng)格大小對計算結(jié)果的影響，分別采用不同網(wǎng)格尺寸dx=0.05 m、0.10 m、0.20 m進(jìn)行數(shù)值模擬，波浪、水深參數(shù)與2.1節(jié)中參考試驗(yàn)相同，計算得到的圓柱結(jié)構(gòu)破波力如圖4所示。同時將計算結(jié)果與試驗(yàn)結(jié)果[16]進(jìn)行對比可知，隨著網(wǎng)格尺寸由0.20 m減小到0.05 m，波浪力時程與試驗(yàn)波浪力曲線擬合得越來越好，破波力計算峰值與試驗(yàn)結(jié)果誤差為10.2%、5.6%、0.6%(圖5)。計算時間上，dx=0.20時僅用時2min，dx減小到0.10 m及0.05 m時用時分別為2.5 h、49.0 h。dx=0.10 m時數(shù)值模擬破波力波峰、波谷與試驗(yàn)偏差較小，整體波形吻合良好，并且和前人的數(shù)值模擬結(jié)果[10]相比誤差小于1%；由于后續(xù)數(shù)值模擬設(shè)計大量工況，綜合考慮模擬結(jié)果捕捉波浪破碎形狀需求以及計算資源問題，后續(xù)驗(yàn)證及分析均采用dx=0.10 m下的計算結(jié)果。

圖4 波浪力時程對比Fig. 4 Time history comparison of wave force

圖5 網(wǎng)格單元大小對波浪力峰值影響Fig. 5 The effect of the cell size on breaking wave force peak

2 模擬結(jié)果

基于驗(yàn)證后的數(shù)值水槽，研究3個工況即波高分別為H=1.1 m、1.3 m、1.6 m，水深d=3.8 m，周期T=4.0 s下不同長寬比(L/D)圓端形橋墩破波力受力特性，各工況中L分別為0、0.5D、1.0D、1.5D、2.0D、2.5D、3.0D、4.0D、5.0D、10.0D，其中L、D表示如圖3，橋墩迎水面圓端的圓心設(shè)在坡頂。

2.1 破碎波浪特性

首先研究不同波高和周期對破碎位置的影響，并對所選工況合理性進(jìn)行研究，在數(shù)值水槽中將波高H、周期T作為輸入變量，計算得到破碎位置與輸入變量間的關(guān)系如圖6所示。

圖6中工況分別為周期4.0 s時波高從1.1 m增加到1.6 m、波高1.3 m時周期從3.8 s增加到5.6 s，破碎位置零點(diǎn)為斜坡頂點(diǎn)。該箱型圖中，橫線自上到下表示最大值、上四分位數(shù)、中位數(shù)、下四分位數(shù)及最小值，加號代表異常值?？梢钥闯?，隨著入射波高的增加，波浪在更前的地方破碎。其中H=1.3 m時，波浪剛好在斜坡頂端破碎。而破碎位置隨周期的變化規(guī)律與波高相反，周期增大時，波浪破碎位置向斜坡頂點(diǎn)后方移動。T=4.0 s時，波浪破碎位置處在斜坡頂端，考慮到分析模擬時需要區(qū)分捕捉波浪破碎形式，采用H=1.1 m，1.3 m，1.6 m的3個工況進(jìn)行分析研究。

圖6 入射波高及周期對破碎位置的影響Fig. 6 The effect of incident wave height and wave period on the breaking location

2.2 墩前流場壓強(qiáng)隨橋墩長寬比變化情況

圖7表示H=1.3 m時計算得破波力峰值時刻測點(diǎn)群α壓強(qiáng)分布，圖8為H=1.3 m時該時刻的破碎點(diǎn)高程?？梢郧逦乜闯鲈诰嚯x平均液面(z=3.8 m)0.70～0.98 m處，即圖8中顯示的高程4.50～4.78 m處，橋墩前端壓強(qiáng)達(dá)到極值，此處即為波浪破碎點(diǎn)，該高度在長寬比變化過程中基本保持不變。測點(diǎn)群壓強(qiáng)大小隨長寬比增加先增大，在L/D=2.0時基本保持平穩(wěn)，不再上升。根據(jù)圖7，在波浪破碎點(diǎn)高度處橋墩圓端前出現(xiàn)低壓的原因?yàn)椴ㄉ嘞聻榭諝鈪^(qū)(如圖2(b)所示)，破碎時波浪首先沖擊橋墩上下兩側(cè)，使沖擊范圍內(nèi)區(qū)域壓強(qiáng)較高。此時中間區(qū)域水浪較少，因此無法產(chǎn)生較大壓力。波浪沖擊橋墩破碎后波高迅速降低，并繞射向墩后匯合，使橋墩兩側(cè)壓強(qiáng)等高線向后下方傾斜。

圖7 H=1.3 m時破波力峰值時刻測點(diǎn)群α豎向壓強(qiáng)分布Fig. 7 Vertical pressure distribution at monitoring point α when breaking wave force peaks under H=1.3 m

圖8 H=1.3 m時破波力峰值時刻最大壓強(qiáng)點(diǎn)高程Fig. 8 Pressure maximal height when force peaks under H=1.3 m

2.3 橋墩長寬比對破波力的影響

數(shù)值計算結(jié)果表明，橋墩所受波浪力受長寬比、波高影響。為找出其相關(guān)性，將3個工況算得的破波力分別進(jìn)行對比，研究每個工況破波力峰值隨長寬比變化情況，并將力隨長寬比的增速表示如圖9。由圖9知，3個工況內(nèi)，隨橋墩長寬比的增加，破波力整體呈先增長后趨于平緩的態(tài)勢，增速逐漸變低。在長寬比到達(dá)一定程度后，破波力會有所降低，但幅度基本不超過5%。圖7中墩前壓強(qiáng)同樣呈現(xiàn)這一規(guī)律，壓強(qiáng)與破波力同趨變化，因此認(rèn)為該壓強(qiáng)是破波力變化的表征因素。并且隨著入射波高的增加，初始破波力增速越來越小，增速變緩會在L/D更小的情況下發(fā)生。入射波高H從1.1 m增加到1.6 m時，初始破波力增速從18.7%降到17.1%及9.8%。以5%為邊界，H=1.1 m時，L/D=2.5后增速變緩；H=1.3 m及H=1.6 m時，長寬比分別在達(dá)到2.0和1.0后破波力便基本穩(wěn)定。不同波高的破碎位置不同，為探究該位置對上述變化趨勢的干擾性，在3個波高對應(yīng)的各工況內(nèi)更改結(jié)構(gòu)位置，使其前端剛好位于波浪破碎點(diǎn)。由2.1節(jié)可知H=1.6 m時破碎點(diǎn)向入射端前移2.4 m，H=1.1 m時破碎點(diǎn)向消波端后移2.0 m，分別記為Δx=-2.4 m，Δx=2.0 m。移動位置后更改結(jié)構(gòu)長寬比，所得破波力結(jié)果如圖10所示，發(fā)現(xiàn)破波力仍先隨長寬比的增大而增大，后逐漸趨于穩(wěn)定。說明所設(shè)工況中結(jié)構(gòu)位置對于破波力的影響遠(yuǎn)小于長寬比變化產(chǎn)生的影響。

圖9 橋墩受破波力隨長寬比變化情況Fig. 9 Breaking wave force change with aspect ratio

圖10 改變橋墩位置對破波力的影響Fig. 10 The effect of the pier’s location on breaking wave force

3 討論與分析

非破碎波通常通過莫里森方程計算拖曳力和慣性力之和得到其波浪力，其中拖曳力系數(shù)及慣性力系數(shù)與KC數(shù)、雷諾數(shù)、粗糙度有關(guān)，而波浪破碎時產(chǎn)生的力可以看作是準(zhǔn)靜態(tài)力(拖曳力和慣性力)與沖擊力[17]之和。從前設(shè)計中通常將準(zhǔn)靜態(tài)力乘以系數(shù)得到破波力，目前多通過Von Karman方法或Wagner方法計算沖擊力。

為研究橋墩所受破波力隨長寬比變化的原因，將數(shù)值計算得到的結(jié)果(合力)通過濾波方式分解，探究準(zhǔn)靜態(tài)力及沖擊力與橋墩長寬比的關(guān)系。數(shù)值模擬中假設(shè)結(jié)構(gòu)為剛體，因此無需考慮結(jié)構(gòu)自振影響。分解步驟如下：1)采用三階Butterworth低通濾波器，截斷頻率設(shè)置為波頻fc=0.25 Hz[18]，低通濾得結(jié)果作為準(zhǔn)靜態(tài)力時程；2)總力時程減去準(zhǔn)靜態(tài)力時程即得到?jīng)_擊力時程。以H=1.3 m下L/D=1.0工況為例，分解過程如圖11所示。

圖11 H=1.3 m且L/D=1.0時濾波后準(zhǔn)靜態(tài)力及沖擊力Fig. 11 Quasi-static force and slamming force filtered when H=1.3 m & L/D=1.0

采用同樣方法對各工況破波力進(jìn)行濾波，得到結(jié)果如圖12所示。

圖12 破波力及分解后準(zhǔn)靜態(tài)力和沖擊力隨長寬比變化趨勢Fig. 12 Breaking wave force and decomposed quasi-static & slamming force tendency with aspect ratio

根據(jù)圖12可知，橋墩破波力峰值隨長寬比增大最大可增加40%，準(zhǔn)靜態(tài)力隨長寬比增大而增大，但增長較為緩慢。當(dāng)長寬比達(dá)到2.0后，準(zhǔn)靜態(tài)力隨L/D變大的增量幾乎為0。其中，除H=1.1 m、H=1.6 m時在L/D從0增長到2.0時準(zhǔn)靜態(tài)力增速大于5%及L/D=10情況外，其余情況下準(zhǔn)靜態(tài)力增速均小于6%。沖擊力隨長寬比變化的趨勢較準(zhǔn)靜態(tài)力更接近總破波力，L/D開始增大時，沖擊力迅速增大，在長寬比達(dá)到2.5后沖擊力開始穩(wěn)定并有輕微下降趨勢。且通過計算可知，模擬中準(zhǔn)靜態(tài)力僅相當(dāng)于分解出的沖擊力的20.3%～41.2%，因此準(zhǔn)靜態(tài)力對總力的影響較小。橋墩所受破波力的變化主要與沖擊力有關(guān)。由于橋墩為圓端形，因此認(rèn)為在長寬比變化過程中，波浪沖擊時pile-up效應(yīng)會受到影響，從而改變橋墩整體受力。

4 結(jié) 語

基于CFD軟件Flow-3D，詳細(xì)闡述了強(qiáng)非線性破碎波發(fā)展下的流體運(yùn)動方程、湍流模型控制方程及數(shù)值模擬方法，建立合理數(shù)值模型和網(wǎng)格劃分方式，針對斯托克斯五階波的3種入射波高，保持橋墩阻水寬度相同，開展了不同長寬比下破碎波對圓端形橋墩的沖擊力大小研究。探討了破碎波浪特性，從橋墩前壓強(qiáng)變化及分解后的準(zhǔn)靜態(tài)力、沖擊力等方面揭示了長寬比對破波力的影響規(guī)律。得到如下主要結(jié)論：

1) 建立的數(shù)值水槽網(wǎng)格單元尺寸大小dx=0.10 m，模擬所得破波力峰值誤差為5.6%，與前人數(shù)值結(jié)果誤差在1%以內(nèi)。破波力波峰、波谷與試驗(yàn)偏差較小，整體波形吻合良好。以兩倍大小縮小網(wǎng)格單元，計算時間相差約一個數(shù)量級，計算精度僅提高約5%，考慮文中數(shù)值模擬涉及大量工況，因此所建尺寸大小0.10 m網(wǎng)格滿足精度和效率需求。

2) 圓端形橋墩破波力隨橋墩長寬比增大，先增大后趨于穩(wěn)定，破波力峰值隨著入射波高的增大而增大。以5%為邊界，H=1.10 m時，L/D=2.5后破波力峰值增速變緩；H=1.3 m及H=1.6 m時，長寬比分別在達(dá)到2.0和1.0后破波力峰值便基本穩(wěn)定。

3) 破波力隨長寬比變化主要與沖擊力有關(guān)。破波力增大時準(zhǔn)靜態(tài)力及沖擊力均增大，但準(zhǔn)靜態(tài)力增速較?。划?dāng)長寬比達(dá)到2.0后，準(zhǔn)靜態(tài)力隨L/D變大的增量幾乎為0。沖擊力隨長寬比變化的趨勢較準(zhǔn)靜態(tài)力更接近總破波力，L/D開始增大時，沖擊力迅速增大，在長寬比達(dá)到2.5后沖擊力趨于穩(wěn)定并有輕微下降。

為提高計算效率，文中研究數(shù)值模擬采用k-ω湍流模型，該模型雖然無法精確地定量描述橋墩周圍的湍流耗散，但可以滿足參數(shù)化定性分析的需要。因此，文中得出的橋墩長寬比對破波力的影響規(guī)律具有工程指導(dǎo)價值。這里僅考慮卷破這一波浪破碎形式，在文中僅初步考慮結(jié)構(gòu)位置與破波力的關(guān)系，后續(xù)會對更多形式的波浪破碎以及結(jié)構(gòu)位置進(jìn)行研究。在研究過程中采用的結(jié)構(gòu)均為剛性，未考慮結(jié)構(gòu)剛度對破波力的影響，未來將會對考慮流固耦合的破波力進(jìn)行探究。此外，現(xiàn)有文獻(xiàn)中尚無全尺寸橋墩的破碎波浪作用試驗(yàn)，這里所驗(yàn)證的數(shù)值水槽也僅為實(shí)驗(yàn)室尺寸，縮尺效應(yīng)對結(jié)果的影響不能忽略。未來有待進(jìn)一步開展足尺橋墩的破碎波浪作用試驗(yàn)和數(shù)值模擬研究。