基于RANS的跨海桥梁高桩承台波浪作用数值模拟

2019-11-29徐博,魏凯

铁道标准设计 2019年12期

徐博,魏凯

(西南交通大学土木工程学院桥梁工程系,成都 610031)

跨海桥梁作为我国沿海交通网络中的关键性结点与重要组成部分[1]，促进了我国海上交通基础建设的进一步发展，具有重要的战略意义。由群桩和承台组成的高桩承台结构，因具备地质适应性好、造价低、施工便捷、风险小等特点，广泛应用于我国跨海桥梁工程[2]。如我国的东海大桥[3]、杭州湾跨海大桥[4]、港珠澳大桥[5]、平潭海峡公铁两用大桥[6]的下部基础形式都采用高桩承台结构。不同于陆地桥梁，海洋桥梁处于复杂多变的海洋环境，通常水深、浪大，波浪作用成为跨海桥梁受力和设计的关键荷载之一[7]。极端波浪荷载给跨海桥梁带来了严峻的挑战，许多跨海桥梁曾在海洋自然灾害下受到严重损坏，给修复重建工作带来了巨大代价。2004年飓风Ivan、2005年飓风 Katrina 和2008年飓风Ike形成的极端波浪和风暴潮先后造成美国多座跨海桥梁被毁[8-10]。2013年，寒潮引发的极端波浪和高潮位造成青岛栈桥发生坍塌。历次事故表明，极端波浪是造成跨海桥梁破坏的主要因素之一，跨海桥梁高桩承台在水面处设有巨大的承台，当极端波浪的波峰作用于承台上时，会产生极大的瞬时波浪冲击荷载。因此，确保跨海桥梁高桩承台在极端波浪下的结构安全，成为我国跨海桥梁发展面临的重要挑战。

对于承台这类大尺度结构物的波浪荷载计算问题，国内外学者开展了大量的研究工作。最早由MacCamy和Fuchs基于势流理论提出了线性绕射理论计算大尺度垂直圆柱的波浪力[11]。Bhatta和Rahman[12]基于绕射理论推导了针对圆形承台波浪力计算的解析表达式。李忠献等[13]提出了考虑绕射效应的圆截面承台波浪作用解析计算方法。刘桦等[14]基于势流理论和莫里森方程提出了计算桩基承台结构波浪力的实用计算方法。杜修力等[15]提出了考虑水体可压缩性的大直径深水圆柱结构动水压力时域和频域的计算方法。势流理论虽然计算效率高，但基本假设中忽略了流体黏性，无法模拟波浪在结构物附件的剧烈变形、湍流、旋涡、破碎等现象，对于此类问题，一些学者提出了基于不可压缩黏性流体理论的CFD计算方法。康啊真等[16]利用浸没边界法处理不规则结构物表面，对大型圆柱波浪力进行数值模拟。谭长建[17]采用大涡模拟方法考虑湍流的影响，通过建立三维波浪与桥墩相互作用数值模型，对桥墩在波浪作用下的受力进行了详细研究。向宝山等[18]通过建立三维数值波浪水槽，研究了承台淹没深度及波浪入射角度对海洋环境中复合桩基水平向波浪荷载的影响。而现有的研究大多针对大型圆形结构，关于跨海桥梁高桩承台极端波浪冲击作用的研究还未充分完善，现有研究成果还无法针对跨海桥梁高桩承台的荷载设计提供有效的指导。因此，研究跨海桥梁高桩承台极端波浪冲击作用具有重要的意义。

通过求解流体运动的雷诺时均Navier-Stokes方程(RANS)和k-ω湍流模型，基于FLOW-3D软件建立了波浪与典型高桩承台桥梁深水基础相互作用的三维数值模型，分析了作用在承台上的波浪力特性、波浪作用下承台底面的压强分布以及周围流场的变化，并研究了承台底面净空对承台波浪力的影响。

1 数值水槽模拟方法

1.1 控制方程

在波浪与结构相互作用研究中模拟波浪运动时，假设水体是不可压缩黏性流体，采用连续性方程和雷诺时均Navier-Stokes方程(RANS方程)作为波浪运动的控制方程

(1)

(2)

式中，i=1，2，3，xi表示x，y，z坐标；ui表示流场各方向的时均速度；Ai表示各方向可流动流体的面积分数；VF表示流体的体积分数；t表示时间；p表示流体压强；ρ表示流体密度；Gi表示流体在各方向上的重力加速度；fi表示各方向上黏滞力引起的加速度，表示如下

(3)

(4)

其中，τb，i表示各方向面上流体剪切应力；Sij表示应变率张量，v表示动力黏度，vT表示湍流黏度。

1.2 湍流模型

选用k-ω湍流模型来计算湍流黏度

vT=k/ω

(5)

式中，ω=ε/k，k和ε由下式计算得到

(6)

(7)

(8)

其中

(9)

式(7)中，σω=2.0，α=13/25，β=β0fβ，其中β0=9/125，fβ为

(10)

其中

(11)

式(11)中

(12)

1.3 边界条件

本文模型中，水槽左侧入流边界定义为波浪边界，采用速度入口法造波；右侧为出流边界，采用Sommerfeld辐射边界条件进行消波；前后两侧设定为对称边界；底部设定为壁面边界条件；液面以上设为压力边界，压强设置为1个标准大气压。采用流体体积法(VOF法)[19]追踪运动液体的自由表面，其基本思想是在计算域中的每个单元内定义一个流体体积函数F，F定义为单元内流体所占有的体积与该单元内可容纳流体体积之比。F是空间和时间的函数，满足以下输运方程

(13)

当F=0时，表示单元内不含流体；当F=1时，表示单元被流体占满；当F介于0和1之间时，表示单元包含有自由表面。

为了更好地消除波浪的反射，在水槽的末端段设置阻尼区域，阻尼层消波的原理是在流体运动控制方程中引入阻尼力，阻尼力和阻尼层中水体流速和阻尼系数有关，阻尼系数沿波浪传播方向线性增大。很多研究表明采用辐射出流边界与阻尼层结合的方式能够大大提高消除波浪反射的效果[20]。

2 数值水槽造波验证

图1 波面时程对比

为了验证数值模型对非线性波的模拟效果，将数值模型计算的波面时间历程与理论值进行对比，波浪参数为波高H=0.15 m，周期T=1.15 s，水深d=0.9 m，符合斯托克斯五阶波浪理论。对比结果如图1所示，图中黑色实线为根据斯托克斯五阶波浪理论计算得出的理论解，红色实心点为数值计算结果。由图1可知，本文数值模型与理论解析结果无论在波峰和波谷都较为吻合，波浪传播的稳定性较好，基本不存在波高的衰减和抬升，说明本模型消波效果良好，可以较为准确地模拟波浪的运动。

3 波浪-承台数值建模

高桩承台桥梁深水基础由承台和桩基组成，与波长相比，承台属于大尺度结构物，而桩基属于小尺度结构物，承台会对结构周围的流场产生影响。承台通常位于波浪主要作用范围内，承担着较大的波浪冲击荷载，承台上承受的波浪荷载是结构安全设计考虑的关键荷载之一。为了研究承台波浪荷载作用，基于上述验证的数值模型，建立如图2所示长14 m，宽2 m，高1.3 m的数值波浪水槽。

图2 数值波浪水槽示意(尺寸非等比例)(单位：m)

考虑到承台高程对高桩承台基础波浪力有较大影响，为了研究承台位置对波浪力的影响，定义承台底面距静水面的距离为净空s(图2)。s为正，说明承台底面位于静水面以上，反之，承台底面位于静水面以下。所建立的高桩承台模型上部是40 cm×26.67 cm×15 cm的矩形承台，下部是由12根直径为5 cm的圆柱组成的群桩，几何模型和详细尺寸如图3所示。

设置了5种不同的净空条件，净空s分别取为4，2，0，-7.5，-15 cm。波浪参数选取如下：水深0.9 m，波高0.18 m，周期1.25 s，选用斯托克斯五阶波浪理论。为了较精确地模拟波浪的运动同时兼顾计算效率，在划分网格时，在水面附近波高范围内网格尺寸设置为Δx=0.025 m，Δy=0.025 m，Δz=0.008 m，对于远离波峰和波谷的区域内在Z向采用稀疏的渐变网格，同时满足一定的渐变率和纵横比以保证计算的收敛性和稳定性。为了能够较为精确地模拟波浪与结构物的相互作用，更好地描述承台与桩基的几何形状、捕捉承台周围的流场形态，在承台周围采用较精细的网格进行局部加密，加密网格的尺寸为Δx=0.012 5 m，Δy=0.012 5 m。同时为了避免由于网格大小突变而导致的数值计算不稳定问题，在网格加密区周围设置渐变网格，使网格从粗网格逐渐过渡到细网格。网格划分如图4所示。

4 结果与讨论

4.1 波浪力结果

图5 承台水平波浪力时程曲线

图6 承台竖向波浪力时程曲线

4.2 压力场分析

在高桩承台结构分析设计时，不仅需要考虑承台承受的波浪合力，还需要了解局部波浪冲击压强。为了分析波浪动水压强在承台底面的分布，从图6中取出一个波周期内几个典型时刻承台底面处的动水压力场云图，如图7所示。由图7可知，波浪动水压强沿水槽宽度方向呈对称分布。时刻t1是波浪与承台底面接触的瞬间，随着波浪的前进，最大波压强作用点也向前推移，到时刻t2时最大冲击压强达到峰值，承台底面受力面积也迅速增大，此时的竖向波浪冲击力达到最大值。到达时刻t3时，承台底面已完全被水体覆盖，而竖向波浪力数值回归到正常的准静力水平。到达时刻t4时，波浪开始脱离承台，承台底面的前部是正压区，后部出现了负压区，此时竖向波浪力为零。随着波浪继续前进，时刻t5承台底面全部承受负向波压强，此时向下的竖向波浪力绝对值达到最大。到达时刻t6时，波浪完全脱离，承台底面上的压强也变为零。