李東洋,張慶河,焦方騫
(天津大學(xué) 水利工程仿真與安全國(guó)家重點(diǎn)實(shí)驗(yàn)室,天津 300072)
斜坡堤是保護(hù)沿海地區(qū)免受潮、浪襲擊的重要水工建筑物,由于斜坡堤迎浪面較緩,波浪爬高往往較大,很容易發(fā)生越浪,因此準(zhǔn)確模擬斜坡堤越浪過程并確定越浪量,對(duì)于斜坡堤的設(shè)計(jì)具有重要意義。
早期對(duì)斜坡堤越浪的研究主要以物理模型試驗(yàn)為主,王聰?shù)萚1]在總結(jié)前人物理模型試驗(yàn)的基礎(chǔ)上,研究了不規(guī)則波作用下?lián)趵藟Ω叨鹊葘?duì)越浪量的影響。近年來,隨著計(jì)算機(jī)技術(shù)和并行計(jì)算的發(fā)展,利用數(shù)值模型研究斜坡堤越浪成為港口海岸工程學(xué)科關(guān)注的熱點(diǎn)。追蹤自由表面運(yùn)動(dòng)的斜坡堤越浪數(shù)值模型主要分為兩大類,一類基于N-S方程,如Van Gant等[2]較早基于二維不可壓縮的N-S方程研究了波浪與寬肩臺(tái)可滲斜坡堤相互作用的水動(dòng)力特性,Losada等[3]基于體積平均的RANS方程和VOF方法建立了COBRAS(Cornell Breaking Waves and Structure)模型,模擬分析了不規(guī)則波在可滲斜坡堤上的越浪過程。關(guān)大瑋[4]基于CFD軟件Flow-3D建立三維數(shù)值波浪水槽,模擬研究了規(guī)則波和不規(guī)則波作用下不可滲光滑斜坡堤的越浪過程。閆科諦等[5]基于OpenFOAM軟件建立數(shù)值波浪水槽,通過網(wǎng)格細(xì)化,刻畫出柵欄板凹槽的形狀,進(jìn)而模擬了規(guī)則波作用下非可滲堤心柵欄板護(hù)面的斜坡堤越浪過程。Xiang等[6]等通過將Fluidity和Y3D耦合的方法,模擬研究了規(guī)則波作用下扭王塊體護(hù)面斜坡堤周圍的水動(dòng)力特性。另外一類模型為無(wú)網(wǎng)格方法,如葉曉文[7]基于不可壓縮SPH數(shù)學(xué)模型研究了規(guī)則波沿可滲斜坡堤的爬坡和越浪過程, Ren[8]基于弱可壓的光滑粒子流體力學(xué)方法(WCSPH)模擬了波浪與可滲潛堤的相互作用,并分析了堤體內(nèi)部的水動(dòng)力特性。溫鴻杰[9]基于改進(jìn)的SPH方法,模擬研究了規(guī)則波在出水堤上的破碎和越浪過程。
總體上來看,目前關(guān)于斜坡堤越浪數(shù)值模擬的研究還主要集中于規(guī)則波,對(duì)護(hù)面塊體進(jìn)行三維全尺度模擬且考慮護(hù)面塊體墊層以下可滲特性的工作也比較少見,這使得利用數(shù)值模型解決斜坡堤越浪等實(shí)際工程問題仍受到一定的制約。因此,利用數(shù)值模型展現(xiàn)真實(shí)海況下斜坡堤越浪過程,需要進(jìn)一步探索??紤]到SPH類模型在壓力脈動(dòng)等方面還存在有待解決的問題,本文將利用Jesus[10]和Lara等[11]基于OpenFoam開源軟件開發(fā)的ihFoam求解器,通過OpenFoam網(wǎng)格劃分工具snappyHexMesh精細(xì)刻畫網(wǎng)格的方法生成扭王字塊,模擬不規(guī)則波作用下扭王塊體護(hù)面斜坡堤越浪過程,以使斜坡堤越浪的數(shù)值模擬在解決工程問題時(shí)達(dá)到一定的實(shí)用性,為斜坡堤的設(shè)計(jì)提供依據(jù)。
本文通過將滲流運(yùn)動(dòng)的Darcy-Forchheimer方程引入到兩相流求解器interFoam中,生成全新的考慮孔隙介質(zhì)流的ihFoam求解器,采用主動(dòng)吸收式速度入口造波,通過求解體積平均的RANS方程,模擬水體和空氣兩相不可壓縮流體的運(yùn)動(dòng)。紊流模型選用RNGk-ε模型,采用VOF法捕捉自由表面。模型控制方程如下
(1)
(2)
(3)
(4)
式中:D50為組成孔隙介質(zhì)的中值粒徑;KC=umT/(nD50),其中um為最大震蕩速度;T為震蕩周期。α和β分別為Darcy-Forchheimer線性項(xiàng)系數(shù)和非線性項(xiàng)系數(shù),2011年Lara等[13]提出了α和β關(guān)于孔隙率和中值粒徑的經(jīng)驗(yàn)公式如下
(5)
當(dāng)n=1時(shí),即純水體情況時(shí),式(2)轉(zhuǎn)化為標(biāo)準(zhǔn)的雷諾時(shí)均方程。
圖1 數(shù)值水槽和邊界條件示意圖Fig.1 Numerical wave flume and boundary condition diagram
數(shù)值水槽和邊界條件示意圖如圖1所示,左側(cè)入口設(shè)置為主動(dòng)吸收造波邊界。右側(cè)出口設(shè)置為輻射消波邊界。數(shù)值水槽前后側(cè)邊界設(shè)置為可滑移固邊界,法向速度為0。底面及護(hù)面塊體設(shè)置為固壁邊界,數(shù)值邊界條件設(shè)置為無(wú)滑移邊界條件,壁面處的流速為0。
為了驗(yàn)證數(shù)值模型的合理性,我們收集已有物理模型試驗(yàn)數(shù)據(jù)與數(shù)值模型結(jié)果進(jìn)行對(duì)比分析。
圖2 物理模型試驗(yàn)尺寸示意圖[14]Fig.2 Physical model size diagram[14]
陳偉秋等[14]通過物模試驗(yàn)研究了不規(guī)則波作用下扭王塊體護(hù)面斜坡堤的平均越浪量。試驗(yàn)是在長(zhǎng)175 m、寬1.2 m、高1.5 m的波浪水槽中進(jìn)行的,試驗(yàn)?zāi)P偷木唧w尺寸如圖2所示,斜坡堤左右兩側(cè)坡度均為1:1.5,堤頂高度為0.5 m,斜坡堤堤底寬為2.11 m,兩端各有坡度為1:2的拋石棱體作為護(hù)底,海側(cè)斜坡鋪設(shè)一層扭王字塊,單個(gè)扭王字塊的高度是6 cm。試驗(yàn)水深0.4 m,入射波為不規(guī)則波,采用JONSWAP譜進(jìn)行模擬,有效周期為2 s,有效波高分別為6.6 cm、7.1 cm、7.6 cm。
表1 多孔介質(zhì)物理參數(shù)表Tab.1 Physical parameter of porous media
模擬過程中除扭王塊體采用全尺度模擬外,斜坡堤堤體其他組成部分均作為孔隙介質(zhì)處理,各層物理參數(shù)見表1。
試驗(yàn)中采用JONSWAP譜進(jìn)行模擬,JONSWAP譜的表達(dá)形式為
(6)
(7)
(8)
(9)
式中:S(f)為波浪譜密度函數(shù);Hs為有效波高;fp為譜峰周期;f為波浪頻率;γ為譜峰因子,取3.3;σ是波浪密度的標(biāo)準(zhǔn)差,當(dāng)f≤fp時(shí),σ=0.07,當(dāng)f>fp時(shí),σ=0.09。
3-a 波面歷時(shí)圖 3-b 計(jì)算波譜與理論波譜對(duì)比圖圖3 Hs=7.1 cm情況下數(shù)值造波結(jié)果Fig.3 Hs=7.1 cm Simulation results of free surface
不規(guī)則波采用線性波疊加法生成[15],根據(jù)等分頻率法將波浪能量分割成70份。因?yàn)閿?shù)值模擬的側(cè)壁采用的是可滑移邊界條件,避免了壁面對(duì)內(nèi)部流場(chǎng)的影響,所以可以利用寬度較窄的計(jì)算水槽反演試驗(yàn)結(jié)果,具體數(shù)值水槽設(shè)置如下:長(zhǎng)15 m,高0.6 m,寬0.24 m,網(wǎng)格水平空間步長(zhǎng)取L/100,為了準(zhǔn)確追蹤自由表面,在其靜水位上下2倍有效波高范圍內(nèi)取30層網(wǎng)格,上下網(wǎng)格漸變,寬度方向上取15層網(wǎng)格。數(shù)值模擬時(shí)間為220 s,取x=12 m(即堤腳)處測(cè)點(diǎn)進(jìn)行分析,以有效波高為7.1 cm為例,圖3-a顯示了其水位歷時(shí)曲線。取后100個(gè)波,即20~220 s的實(shí)測(cè)數(shù)據(jù)進(jìn)行快速傅里葉分析,得到不規(guī)則波浪的能量譜。圖3-b顯示了計(jì)算模擬波譜與理論波譜的對(duì)比情況,由圖可知兩者比較接近。模擬有效波高與理論有效波高誤差為2.4%,滿足《波浪模型試驗(yàn)規(guī)程》[16]要求。
圖4 AutoCAD三維建模中扭王字塊及扭王字塊排放示意圖 圖5 扭王塊體和胸墻經(jīng)snappyHexMesh處理后局部網(wǎng)格剖分示意圖Fig.4 AutoCAD three-dimensional model of accropode blocks Fig.5 Accropode blocks and parapet local mesh encryption by snappyHexMesh
OpenFOAM提供了強(qiáng)大的網(wǎng)格生成工具,通過自帶的網(wǎng)格處理工具snappyHexMesh將扭王字塊和胸墻轉(zhuǎn)變成貼體網(wǎng)格模擬,可以逼近真實(shí)工況下斜坡堤上的波浪爬坡和越浪。
首先參照《防波堤設(shè)計(jì)與施工規(guī)范》[17]中規(guī)定的扭王字塊的尺寸,使用AutoCAD三維建模繪制出扭王字塊的外部輪廓。周雅等[18-19]通過物理模型試驗(yàn)研究得出扭王字塊護(hù)面塊體的排放型式對(duì)越浪量的影響,因此為了簡(jiǎn)化模型,將扭王字塊規(guī)則排放在斜坡堤上,如圖4所示。圖5顯示了精細(xì)網(wǎng)格刻畫后的扭王塊體和胸墻局部網(wǎng)格形狀。由圖可以看出數(shù)模所采用的網(wǎng)格能較好描述扭王護(hù)面塊體和胸墻的外部形狀。
6-a越浪量歷時(shí)圖6-b累計(jì)越浪量圖圖7 平均越浪量計(jì)算值與試驗(yàn)值比較圖6 Hs=7.1cm情況下越浪量計(jì)算結(jié)果Fig.6Hs=7.1cmSimulationresultsofovertoppingdischargeFig.7Comparisonofnumericalaveragedischargeandexperimentalaveragedischarge
8-a t=174.36 s 9-a t=174.36 s
圖6顯示了有效波高為7.1 cm情況下的波浪越浪量歷時(shí)圖6-a和累計(jì)越浪量圖6-b,在有效周期一定的情況下,越浪量隨有效波高的增大而明顯增大。取20~220 s(100個(gè)波)時(shí)間內(nèi)的越浪量計(jì)算平均值以得到平均越浪量。圖7顯示了平均越浪量數(shù)值模擬結(jié)果與試驗(yàn)值的比較,由圖可知斜坡堤越浪的數(shù)值模擬結(jié)果與試驗(yàn)結(jié)果吻合較好,表明利用snappyHexMesh網(wǎng)格處理工具生成貼體網(wǎng)格的方法對(duì)扭王護(hù)面塊體進(jìn)行全尺度模擬,對(duì)墊層以下塊石利用滲流模型模擬斜坡堤越浪可以獲得較好的效果。
8-b t=174.56 s 9-b t=174.56 s
8-c t=174.76 s 9-c t=174.76 s圖8 最大單波越浪過程三維示意圖 圖9 堤前附近流場(chǎng)圖 Fig.8 3D process of the largest Fig.9 Variation diagram of single-wave overtopping velocity fields
圖8顯示了有效波高為7.6 cm,在174.5 s前后發(fā)生越浪時(shí)波面的三維過程,圖9顯示了Y=0.1 m,XZ剖面扭王字塊附近的流場(chǎng)變化。從圖中可以看出,越浪發(fā)生前,波浪沿著坡面向上爬升,在向上水體的帶動(dòng)下,水面抬高,波浪拍擊胸墻表面,直至越過堤頂形成越浪。在整個(gè)越浪過程中,由于扭王字塊護(hù)面的作用,堤內(nèi)的滲流流速相對(duì)很小。
從越浪過程、越浪波面和流速變化過程看,本文所建立的模型較準(zhǔn)確地反映了波浪與扭王塊護(hù)面斜坡堤的相互作用,具備了為斜坡堤設(shè)計(jì)優(yōu)化提供支持的能力。
表2模型與原型波要素
Tab.2 Wave elements under experimental model and prototype size
類別有效波高(m)有效周期(s)水深(m)實(shí)驗(yàn)室尺度0.07120.40原型尺度0.071×25=1.7752×251/2=100.40×25=10.00
現(xiàn)階段對(duì)于原型尺寸下的波浪與斜坡堤相互作用的研究,大多是通過比尺縮放,進(jìn)行實(shí)驗(yàn)室物模試驗(yàn),再將試驗(yàn)值根據(jù)《波浪模型試驗(yàn)規(guī)程》換算成原型,而直接研究原型尺寸下斜坡堤越浪的工作還比較少。因此,為進(jìn)一步論證模型可否應(yīng)用于描述現(xiàn)場(chǎng)尺度條件下的斜坡堤越浪,首先依據(jù)《波浪模型試驗(yàn)規(guī)程》[16]將文獻(xiàn)[14]中的物理模型試驗(yàn)按照1∶25的比尺放大成原型工況,通過模擬該尺度下的斜坡堤越浪,并與文獻(xiàn)[14]按照重力相似計(jì)算得出的越浪量進(jìn)行比較,從而驗(yàn)證模型是否可用于模擬計(jì)算現(xiàn)場(chǎng)尺度下的斜坡堤越浪。以模型有效周期為2 s,有效波高0.071 m為例,按照重力相似準(zhǔn)則確定原型波要素如表2所示。
圖10 原型工況下累計(jì)越浪量Fig.10 The accumulated overtopping discharge under prototype size
原型模擬數(shù)值水槽長(zhǎng)×寬×高=375 m×6 m×15 m,依然按照上述方法建模。圖10顯示了原型尺度下累計(jì)越浪量的數(shù)值模擬結(jié)果,單個(gè)越浪發(fā)生時(shí)其實(shí)歷時(shí)很短,可以看成是一個(gè)瞬時(shí)過程,圖中折線每次的抬升就代表發(fā)生了一次越浪,“臺(tái)階”的高度即為單個(gè)越浪發(fā)生時(shí)的越浪量,從而可以得出,單波最大越浪量發(fā)生在950 s左右,單波最大越浪量約為1.3 m3/m。取100~1 100 s(100個(gè)波)時(shí)間內(nèi)的越浪量計(jì)算平均值得到平均越浪量,原型尺度下數(shù)值計(jì)算的越浪量為7.89×10-3m3/m/s。根據(jù)模型試驗(yàn)越浪量折算到原型的平均越浪量為6.68×10-5×251.5=8.36×10-3m3/m/s,與原型數(shù)值模擬結(jié)果較為接近。合理的模型試驗(yàn)可以反映原型結(jié)果是比較公認(rèn)的,因此,從原型直接模擬結(jié)果與模型試驗(yàn)折算到原型的越浪量比較吻合來看,利用OpenFOAM直接模擬原型工況下的斜坡堤越浪是可行的,并推薦采用原型尺寸進(jìn)行數(shù)值研究,進(jìn)而優(yōu)化確定設(shè)計(jì)階段中斜坡堤的斷面型式和頂高程等。
本文基于OpenFOAM建立數(shù)值波浪水槽,控制方程采用了體積平均的RANS方程和包含非線性項(xiàng)的滲流方程描述斜坡堤內(nèi)部的水體滲流流動(dòng)。利用AutoCAD三維建模繪制出扭王塊體的外部輪廓,利用snappyHexMesh網(wǎng)格劃分工具對(duì)斜坡堤扭王護(hù)面塊體進(jìn)行全尺度網(wǎng)格精細(xì)劃分。采用建立的模型模擬研究了正向入射不規(guī)則波與扭王塊體護(hù)面斜坡堤的相互作用,模擬結(jié)果表明,平均越浪量的計(jì)算值與實(shí)驗(yàn)結(jié)果吻合較好。通過將實(shí)驗(yàn)室尺度放大成原型尺度,模擬研究了原型條件下的斜坡堤越浪,從計(jì)算結(jié)果可以看出,平均越浪量的計(jì)算值與模型通過重力相似準(zhǔn)則折算成原型的結(jié)果較為吻合。從實(shí)驗(yàn)室尺度和原型尺度模擬結(jié)果看,對(duì)護(hù)面塊體進(jìn)行全尺度模擬的數(shù)值波浪水槽目前已可以較為合理地描述復(fù)雜護(hù)面塊體斜坡堤的越浪過程,可以應(yīng)用于設(shè)計(jì)階段對(duì)斜坡堤斷面型式和頂高程的優(yōu)化確定。
需要指出的是,本文僅給出了波浪正向入射情況下斜坡堤的越浪模擬結(jié)果,有關(guān)護(hù)面塊體本身的穩(wěn)定性以及斜向入射波與斜坡堤的相互作用有待進(jìn)一步研究。
[1]王聰, 陳國(guó)平, 嚴(yán)士常, 等. 不規(guī)則波作用下斜坡堤越浪量試驗(yàn)研究[J]. 水運(yùn)工程, 2017(2): 49-52.
WANG C, CHEN G P, YAN S C,et al. Experimental study on wave overtopping of sloping seawall under irregular wave action[J]. Port & Waterway Engineering, 2017(2): 49-52.
[2]Van Gent M R A. Wave action on and in permeable structures[D]. The Netherlands: Delft University, 1995.
[3]Losada I J, Lara J L, Guanche R, et al. Numerical analysis of wave overtopping of rubble mound breakwaters[J]. Coastal Engineering, 2008, 55(1): 47-62.
[4]關(guān)大瑋. 規(guī)則波與不規(guī)則波的海堤越浪數(shù)值模擬[D]. 廣州:華南理工大學(xué), 2016.
[5]閆科諦, 張慶河. 柵欄板護(hù)面斜坡堤越浪數(shù)值模擬研究[J]. 中國(guó)港灣建設(shè), 2016, 36(2):16-19.
YAN K D, ZHANG Q H. Numerical simulation of overtopping on sloping dike with fence panels[J]. China Harbour Engineering, 2016, 36(2):16-19.
[6]Xiang J, Latham J P, Vire A, et al. Coupled fluidity/Y3D technology and simulation tools for numerical breakwater modelling[J]. Coastal Engineering Proceedings, 2012, 1(33):66.
[7]葉曉文. 斜坡式防波堤越浪過程的SPH模擬[D]. 大連:大連理工大學(xué), 2011.
[8]Ren B, Wen H, Dong P, et al. Improved SPH simulation of wave motions and turbulent flows through porous media[J]. Coastal Engineering, 2016, 107: 14-27.
[9]溫鴻杰. 非線性波浪與可滲結(jié)構(gòu)相互作用的SPH模型[D]. 大連:大連理工大學(xué), 2016.
[10]Jesus M D, Lara J L, Losada I J. Three-dimensional interaction of waves and porous coastal structures. Part I: Numerical model formulation[J]. Coastal Engineering, 2012, 64:57-72.
[11]Lara J L, del Jesus M, Losada I J. Three-dimensional interaction of waves and porous coastal structures. Part II: Experimental validation[J]. Coastal Engineering, 2012, 64: 26-46.
[12]Higuera P, Lara J L, Losada I J. Three-dimensional interaction of waves and porous coastal structures using OpenFOAM?. Part I: Formulation and validation[J]. Coastal Engineering, 2014, 83: 243-258.
[13]Lara J L, Losada I J, Maza M, et al. Breaking solitary wave evolution over a porous underwater step[J]. Coastal Engineering, 2011, 58(9):837-850.
[14]陳偉秋, 陳凌彥, 王登婷,等. 斜坡堤后坡砌石護(hù)面穩(wěn)定厚度的模型試驗(yàn)研究[J]. 水運(yùn)工程, 2016(6): 93-98.
CHEN W Q, CHEN L Y, WANG D T, et al. Experimental study on thickness of inner block revetment of sloped seawall[J]. Port & Waterway Engineering, 2016(6): 93-98.
[15]俞聿修, 柳淑學(xué). 隨機(jī)波浪及其工程應(yīng)用[M]. 大連:大連理工大學(xué)出版社, 2011.
[16]JTJ/T234-2001,波浪模型試驗(yàn)規(guī)程[S].
[17]JTS154-1-2011, 防波堤設(shè)計(jì)與施工規(guī)范[S].
[18]周雅, 林登榮, 李慶銀, 等. 不規(guī)則波作用下斜坡堤越浪量試驗(yàn)研究[J]. 水道港口, 2016(4):331-335.
ZHOU Y, LIN D R, LI Q Y, et al. Experimental research of wave overtopping on sloping dike under irregular waves[J]. Journal of Waterway and Harbor, 2016(4): 331-335.
[19]鐘雄華,陳國(guó)平,嚴(yán)士常,等.不同擺放方式扭王字塊體穩(wěn)定性研究[J].水道港口,2016(5):479-483.
ZHONG X H, CHEN G P, YAN S C,et al. A stability analysis of different stack types of accropode[J]. Journal of waterway and harbor, 2016(5):479-483.