海洋工程中直立圆柱波浪爬升问题的数值研究❋

唐　鹏， 于定勇， Bai Wei， 田　艳

(1.中国海洋大学工程学院，山东 青岛 266100; 2.新加坡国立大学，新加坡 119077)



海洋工程中直立圆柱波浪爬升问题的数值研究❋

唐鹏1， 于定勇1， Bai Wei2， 田艳1

(1.中国海洋大学工程学院，山东 青岛 266100; 2.新加坡国立大学，新加坡 119077)

本文基于Navier-Stokes方程和连续性方程，利用VOF方法捕捉自由液面和松弛造波法实行造波，建立基于黏性流体理论的三维数值波浪水槽，并在该数值波浪水槽内对圆柱周围的波浪爬升效应进行数值模拟。选用两种不同的波浪参数，探讨波陡参数对柱体周围的波浪爬升效应及其所受载荷的影响。与势流模型相比，本文采用的黏性流模型可以较好的捕捉实验中出现的二次波峰现象。计算表明，本文建立的数值水槽可以较好用于直立圆柱波浪爬升问题的计算。

波浪爬升；圆柱；数值水槽；OpenFOAM

引用格式：唐鹏， 于定勇， Bai Wei， 等. 海洋工程中直立圆柱波浪爬升问题的数值研究[J]. 中国海洋大学学报(自然科学版), 2016, 46(10): 116-122.

TANG Peng, YU Ding-Yong, Bai Wei， et al. Numerical simulation of wave run-up on cylindrical offshore structures[J].Periodical of Ocean University of China, 2016, 46(10): 116-122.

波浪在传播的过程中遇到结构物时，会产生绕射效应，波浪爬升可描述为入射波浪在传播过程中遇部分浸入水中的结构物时波浪沿其表面向上爬升的非线性现象。波浪爬升过程的物理解释为: 当入射波浪遇部分浸入水中的柱体时，其周围的自由波面将发生剧烈的变化，其中部分水体绕过立柱继续朝前传播，而受立柱阻碍的另一部分波浪开始实现能量的转化，即由动能转化成势能，导致波浪沿立柱表面迅速向上攀升[1]。对海洋平台，随着波浪爬升往往会造成平台甲板上浪和对平台底板抨击作用，从而引起海洋平台的破坏。为避免波浪爬升对浮式平台的破坏，在设计中，应该充分考虑波浪爬升，对波浪爬升幅度有准确的估计，保证平台甲板有足够的标高。同时，平台甲板也不能过高，甲板过高往往会产生工程成本增加和对平台稳定性方面产生影响，因此对结构物爬升问题研究具有非常重要的现实意义。

Galvin和Hallermeier等[2]首次采用物理模型试验方法对圆柱的波浪爬升进行了研究。其通过采用安装于圆柱壁面附近的多个浪高仪来获取其圆柱周围自由水面的波形分布。Nielsen[3]采用规则波对圆形和带导角的方柱进行了物理模型试验，研究了两种柱体形状对波浪爬升效应和柱体周围波面分布的影响，其结果表明，与圆柱相比，带导角的方柱周围的局部波浪放大率比圆柱更大，非线性效应更强。Contento G[4]等针对波陡较大、深水规则波浪下圆形直立柱表面的波浪爬升效应进行试验研究。当A/D≫1(A为入射波波幅，D为圆柱直径)时候，必须考虑黏性效应的影响，黏性效应对波浪爬升的贡献较大，不宜采用势流理论求解。Donald G. Danmeier[5]基于黏性流模型，利用ComFLOW软件数值模拟了规则波作用下重力式平台的波浪爬升效应。结果表明，入射波陡较小时，线性波浪爬升占主导，随着波陡的不断增大，平台周围的波峰逐渐尖瘦，波谷趋于平坦，波浪爬升的高阶非线性特征越来越明显，且在平台的迎浪面和侧面位置，一个波浪周期内出现二次波峰现象。

OpenFOAM(Open Source Field Operation and Manipulation)是一个在Linux下运行的计算流体力学(CFD)类库，基于有限体积法FVM(Finite Volume Method)求解偏微分方程，是一个针对不同的的流动编写不同的C++程序集合。本文利用有限体积法求解N-S方程、VOF方法捕捉自由液面和松弛造波法实行造波，建立三维数值波浪水槽。通过修改OpenFoam中已有两相流求解器interFoam中的interFoam.C文件编译适合本文的求解器waveFoam，进行波浪与直立圆柱爬升问题的数值研究。

1　数值方法

1.1 控制方程

本文中，固液两相均为不可压缩流体，不可压缩黏性流的连续方程和N-S方程表示为：

连续方程：

(1)

N-S方程：

(2)

式中：ui为速度，i,j=1,2,3代表三个不同维度，即u1=u，u2=v，u3=w；p为压力；ρ为流体密度；μ为分子黏性系数；t为时间；gi为重力。

VOF方法[6]是目前常用捕捉自由液面的方法之一，通过求解体积分数函数来捕捉自由液面。

(3)

γ为流体单元的体积分数，0≤γ≤1。在空气中γ=0，在水中γ=1，自由表面位置0<γ<1。

1.2 数值水槽

本文采用Jacobsen N G[7]提出的松弛造波法实现数值造波和消波，建立三维数值波浪水槽。松弛造波法首先是由Mayer等[8]在利用海绵层(Sponge layer)消波方法中提出的。松弛造波法的基本原理是在造波区和消波区内每一时刻对速度以如下形式进行更新。

u=(1-αR)utarget+αRumodel，

(4)

式中；umodel是通过求解N-S方程得到的速度值；utarget是期望得到的目标波浪速度值；αR是与空间位置有关的加权函数, 在边界位置αR=0。可以通过下式来确定：

(5)

2　数值计算和结果分析

2.1 三维数值波浪水槽的建立和验证

为了验证本文数值计算方法是否正确，将利用新建立的三维数值波浪水槽模型来实现线性波的模拟。该数值波浪水池长20m，宽度为3m，高度为1m，水深h=0.65m。x轴正向与波浪传播方向一致。反射波往往会影响到实验的结果，因此消波效果的对于数值水槽模拟非常重要，本算例在出口处3m设置为消波区，模型相关参数如图2所示。边界条件设置为：入口边界(inlet)patch；出口边界(outlet)patch；底部边界(bottom)wall；顶部边界(atmosphere)patch；前后壁面(front和back)边界wall。采用的波浪参数为：波高H=0.08m，周期T=1.5s，波长λ=3.058m。

图1　入口和出口区域松弛层中αR和χR变化图Fig.1　Sketch of the variation of αR(χR) for both inlet and outlet relaxation zones

图2　数值波浪水池示意图Fig.2　Schematic diagram of numerical wave tank

合理的网格分布对生成理想的波浪场是非常重要的，本文的数值计算中，在水平方向上采用均匀网格,在竖直方向上采用非均匀网格，自由液面附近对网格进行加密处理从而更加准确的获得自由表面位置。表1给出了数值计算中采用的3种网格参数。

表1　数值计算网格参数Table 1　Mesh information for simulation

建立数值波浪水槽是为了实现对波浪条件的模拟，评价和验证一个数值波浪水槽是否可行的一个重要参数为波高，图3给出了x=7.25m处位移历程曲线与理论解的对比图，可以看出生成波形稳定，数值计算波高接近理论波高值，说明本模型可以很好的实现对波浪条件的模拟。图4给出了t=15s时刻下的3种网格生成的波形与理论解的对比图，从图中可以看出，当波浪传递到水槽的最右边界时，波高已经趋于零，说明本模型设置的消波区有很好的消波作用，消波区能够很好地吸收波能。随着网格的加密，三维数值水槽生成的波形图更加符合理论解。同时，由于波浪的非线性产生的紊动耗散产生的原因，产生波形与理论解也存在一定差异。

图3　x=7.25m处位移历程曲线与理论解的对比图Fig.3　Comparison of the elevation between simulation and analytical result at x=7.25m

图4　t=15s时刻下的3种网格的波形与 理论结果的对比图Fig.4　Comparison of the wave profile between simulation and analytical results at t=15s

因此，本文建立的三维数值波浪水槽可以较好地模拟波浪水槽的基本特性，可以用于模拟波浪水槽中进行的断面及三维实验。

2.2 波浪爬升数值计算

本算例模拟的是规则波绕过直立圆柱的爬升现象。几何模型采用类似上节数值水槽的几何模型。长度x方向为20m，其中x=17～20m设置为消波区，同时在x=0～4m设置入口消波区，用来消除波浪与结构物相互作用产生的反射波，消除反射波对入口造波的影响。宽度z方向为7m，高度y方向为1m，水深h=0.65m，圆柱直径d=2r=0.5m，距离入口边界为7.5m。采用两种波浪参数，算例一：波高H=0.08m，周期T=1.5s，波长λ=3.058m，波陡H/λ=0.026；算例二：波高H=0.19m，周期T=1.5s，波长λ=3.058m，波陡H/λ=0.062。

数值计算中，在水平方向上采用非均匀网格,靠近圆柱处采用加密网格，Δxmin=Δzmin=0.005m，在竖直方向上采用非均匀网格，自由液面附近对网格进行加密处理从而更加准确的获得自由表面位置，Δymin=0.01m。为考虑波浪的紊流对结构物的影响，本文采用基于k-ε两方程湍流模型，其中k为湍动能，ε为湍动能耗散率。图5给出了圆柱周围网格细部图。数值计算中采用的最大柯朗数为0.25，采用可变时间步长保证柯朗数满足要求。

图6给出了算例一和算例二圆柱在波浪作用下的无量纲化的水平波浪力历程曲线。其中：ρ为流体密度；g为重力加速度；r为圆柱半径；a为入射波的波幅。可以看到由于波浪的周期性运动，圆柱也受到水平的周期振荡力。

图6　水平波浪力历程曲线Fig.6　Time history of the horizontal force on the cylinder

图7给出了算例一中圆柱表面在x=7.25m处位移与图3中水槽处位移历程对比图，即有结构物和无结构物相同位置处的自由表面处位移随时间变化对比曲线。从图中明显可以看出，有结构物的自由表面处位移增大，在t=9.5、11、12.5、14和15.5s时波峰作用该位置波幅达到最大，最大波幅为a=0.0477m，意味着波浪爬升位移达到0.0077m，爬升幅度与波幅的比值为19.25%。图8给出了算例二中圆柱表面在x=7.25m圆柱自由表面处位移图，其中最大波幅为a=0.231m，意味着波浪爬升幅度达到0.136m，爬升幅度与波幅的比值为143.15%，接近1.5倍的入射波波幅。Isaacson[9]和Kriebel[10]模型试验表明在极端情况下，波浪作用在圆柱的爬升幅度最大可以超过1.5倍的波幅。随着波陡H/λ的增大，圆柱周围的波峰逐渐尖瘦，有更多的动能转化为势能，导致波浪沿立柱的爬升幅度增大，尤其在圆柱的迎浪面，因此在波浪作用下形成的波浪爬升对于结构物的影响是不能忽视。同时发现当波面下降到最低位置时，波面出现了小的波动，也就是二次波峰(Secondary crest)现象。

图7　自由表面处位移历程对比图Fig.7　Comparison of the surface elevation time history

图8　算例二 圆柱表面在x=7.25m处位移图Fig.8　Time history of the surface elevation at x=7.25m in case 2

图9　算例二 圆柱表面在x=7.5m处位移图Fig.9　Time history of the surface elevation at x=7.5m in case 2

图9给出了算例二中圆柱表面在x=7.5m处自由表面位移图，相比图8可以更加明显发现，当波面下降到最低位置时，波面出现了较大波动，波浪爬升的高阶非线性特征越来越明显，且在圆柱的迎浪面和侧面位置，一个波浪周期内出现二次波峰现象。Morris Thomas[11]指出二次波峰现象是由于入射波的波前在绕过柱体两侧之后在柱子背后相遇并相互干扰，由此产生的自由波会向上游方向回流并与上游过来的波相互叠加，从而导致二次波峰的出现。Trulsen[12]基于完全非线性势流模型数值计算了圆柱爬升问题，将数值计算结果与实验结果对比中，实验结果中也出现了二次波峰现象，而基于完全非线性势流模型和二阶势流结果没有捕捉到二次波峰现象，说明二次波峰现象确实存在，势流模型无法捕捉到该现象，而本文基于黏性流模型可以很好的捕捉这一现象。

图10和11分别给出了算例一和算例二圆柱不同点在一个周期上的位移图。从图中可以看到，图形以α=180°呈对称变化，在t=(7+5/15)T时，圆柱迎浪面处(α=180°)表面位移达到最大。同时，可以明显看出，算例二(图11)圆周上位移变化幅度比算例一(图10)大，说明在波陡较大的情况下，圆柱的自由表面变化更加剧烈。

图10　算例一圆周不同点在一个周期上的位移图Fig.10　Wave profiles on the waterline of the cylinder at a whole wave period in case 1

图11　算例二 圆周不同点在一个周期上的位移图Fig.11　Wave profiles on the waterline of the cylinder at a whole wave period in case2

图12给出了算例二在一个完整周期圆柱附近的自由表面位移等值量图。从图中可以清楚的看到在一个完整波浪周期作用下圆柱上的位移变化。与图11对应，在t=(7+5/15)T时，圆柱迎浪面(α=180°)处自由表面位移达到最大，即在此处波浪爬升幅度最大。

图12　算例二 在一个完整周期圆柱附近的 自由表面位移等值量图Fig.12　Wave profiles contour of the cylinder at a whole wave period in case2

3　结论

(1)基于Navier-Stokes方程和连续性方程，利用VOF方法捕捉自由液面和松弛造波法实行造波，建立基于黏性流体理论的三维数值波浪水槽。通过数值计算结果与理论结果对比，对比结果符合性较好，可以证明本文建立的三维数值波浪水槽可以很好的实现数值造波和消波，较好地模拟了波浪水槽的基本特性。

(2)利用建立的数值水槽对圆柱周围的波浪爬升效应进行数值模拟。选用两种典型工况条件，其中当入射波陡较小时，波浪爬升以线性特征为主，波浪爬升波幅较小。随着波陡的不断增大，圆柱周围生成的波峰逐渐尖瘦，波谷开始出现波动，出现了与实验相似的二次波峰现象，波浪爬升波幅增大，接近1.5倍的入射波波幅，此时波浪爬升的高阶非线性特征占主导。与势流模型相比，本文采用的黏性流模型可以较好的捕捉实验中出现的二次波峰现象。

(3)在入射波波陡较大的情况下，通过比较圆周不同点在一个周期上的位移可以发现，圆柱附近的自由表面变化幅度也较大。在圆柱迎浪面(α=180°)处自由表面位移达到最大。

[1]单铁兵, 杨建民, 李欣, 等. 深海平台立柱周围波浪非线性爬升研究进展[J]. 海洋工程, 2012,30(1):151-160.

Shan T, Yang J, Li X, et al. Review of the research on non-linear wave run-up around columns of deepwater platform[J]. The Ocean Engineering,2012, 30(1):151-160.

[2]Galvin C J,Hallermeier R J. Wave runup on vertical cylinders[J]. Coastal Engineering Proceedings, 1972,1(13):65-71.

[3]Nielsen F G. Comparative study on airgap under floating platforms and run-up along platform columns[J]. Marine structures, 2003,16(2):97-134.

[4]Contento G, Francescutto A, Lalli F. Nonlinear wave loads on single vertical cylinders:Pressure and wave field measurements and theoretical predictions[C]. San Francisco： International Society of Offshore and Polar Engineers, 1998.24-32.

[5]Danmeier D G, Seah R K, Finnigan T, et al. Validation of wave run-up calculation methods for a gravity based structure[C]. San Francisco： International Conference on Offshore Mechanics and Arctic Engineering, 2008.265-274.

[6]Hirt C W, Nichols B D. Volume of fluid (VOF) method for the dynamics of free boundaries[J]. Journal of Computational Physics, 1981,39(1):201-225.

[7]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):1073-1088.

[8]Mayer S, Garapon A, Soerensen L S. A fractional step method for unsteady free-surface flow with applications to non-linear wave dynamics[J]. International Journal for Numerical Methods in Fluids, 1998,28(2):293-315.

[9]Isaacson MdSQ. Wave runup around large circular cylinder[J]. Journal of the Waterway Port Coastal and Ocean Division, 1978, 104(1): 69-79.

[10]Kriebel D L. Nonlinear wave interaction with a vertical circular cylinder. Part II: Wave run-up[J]. Ocean Engineering, 1992, 19(1): 75-99.

[11]Morris-Thomas M T. Wave run up on vertical columns of an offshore structure[D]. Perth :The University of Western Australia, 2003.

[12]Trulsen K, Teigen P. Wave scattering around a vertical cylinder: Fully nonlinear potential flow calculations compared with low order perturbation results and experiment[C]. San Francisco： International Conference on Offshore Mechanics and Arctic Engineering,2002.359-367.

责任编辑陈呈超

Numerical Simulation of Wave Run-up on Cylindrical Offshore Structures

TANG Peng1, YU Ding-Yong1, Bai Wei2, TIAN Yan1

(1.College of Engineering, Ocean University of China, Qingdao 266100, China； 2.National University of Singapore, 119077, Singapore)

Wave run-up and wave impact can cause unexpected damage to offshore structures. Therefore, the design of offshore structures requires accurate predictions of the maximum wave elevation to maintain sufficient airgap below the platform deck. Accurate prediction of wave run-up can both help reducing building costs and avoid the risk of wave impact and damage to the platform. For the increasing number of offshore platforms built for ocean oil and gas exploration, the investigation of wave run-up becomes more and more significant for the design of fixed offshore structures. In this study, the interaction between wave and a surface-piercing vertical cylinder is investigated numerically based on waveFoam solver, which is developed under the framework of the OpenFOAM, an open source code library writing with fully objectoriented C++ language. The governing equations are discretized by using the finite volume method. Two fluids are considered and the volume of fluids (VOF) technique is employed to capture water-air interface. A 3D numerical wave tank of the viscous fluid theory was created based on Navier-Stokes equations and continuity equations, used by VOF method to trace free surface and relaxation wave method to generate wave. Meanwhile make use of the 3D wave tank to simulate wave run-up on a vertical cylinder. The model is validated by compare the extreme wave profile with analytical results. Then the simulation of extreme wave-cylinder interactions is carried out. Two different wave parameters was adopted, wave steepness parameters were discussed for the wave effect and run-up around the vertical cylinder. Compared with the potential flow model, viscous flow model used in this paper can better capture the secondary crest phenomenon appeared in the experiments. Calculations show that the numerical tank established in this paper can be better used to solve the problem of wave run-up on cylindrical offshore structures.

wave run-up; vertical cylinder; numerical wave tank; OpenFOAM

国家留学基金委项目(CSC 201306330022)资助

2014-11-12；

2015-02-15

唐鹏(1985-)，男，博士生。tangpeng851008@163.com

TV 139

A

1672-5174(2016)10-116-07

10.16441/j.cnki.hdxb.20140342

Suppoted by China Scholarship Council(CSC 201306330022)