2. 中海油能源发展股份有限公司天津分公司;
3. 韩国海洋大学船舶与海洋系统工程系
2. Tianjin Company, CNOOC Energy Technology & Services Limited;
3. Department of Naval Architecture and Ocean Systems Engineering, National Korea Maritime & Ocean University
0 引言
随着海洋石油工业技术的发展和海上油田水深的增加,水下生产系统已成为深水油气田开发的必然选择。水下生产装备等大型结构物,如水下采油树、水下管汇及PLET等在下放安装过程中所处作业环境恶劣,同时由于其自身重力较大,在穿越飞溅区阶段结构物底部将受到巨大的波浪砰击载荷,有可能导致结构物发生变形和弯曲等结构受损,引起安装失败甚至重大工程事故[1-6]。鉴于上述情况,亟需提出一种有效的数值模拟工具,可以准确预测结构物所能承受的临界波浪砰击载荷,以判断水下装备穿越飞溅区时是否存在潜在危险。
结构物冲击入水时,不但固体受到水面的强烈砰击载荷,而且因固体的运动与变形也会导致流场形态发生改变,即流体与固体之间存在着强非线性耦合作用。针对具有复杂边界条件的入水冲击问题,很难通过解析法求解,目前研究该问题的主要方法为水池试验与数值模拟。黄震球等[7]采用改变翼缘高度与翼缘间距的试验方案,研究了翼缘高度和翼缘间距对平底结构砰击压力的影响。邹丽等[8]通过自由落体的入水方式,分别在静水和规则波中开展了不同横剖面的曲面楔形体入水砰击问题试验研究。彭晟等[9]对新一代江海直达船艏部利用落水试验的方法,进行了一系列不同落水高度及不同入水角度的入水砰击试验。
随着计算机技术及计算流体力学(CFD)的飞速发展,数值模拟方法以其经济性、可靠性与可重复性的优点成为更多学者研究入水冲击问题的重要手段,而且物理试验的一些局限性,如试验设备的可重用性较低、复杂模型的制作、试验成本较高及试验周期性较长等因素均可通过建立数值水槽轻易克服。卢炽华等[10]通过数值方法研究了弹性结构入水冲击中的弹性效应,通过引入射流单元及最大射流厚度处理了冲击引起的射流问题。孙琦等[11]基于有限元软件LS-DYNA建立了三维弹性弹丸体在多物质流体介质中运动的有限元动力分析模型,采用ALE方法对弹丸体撞水过程进行了流固耦合模拟。杨力等[12]采用ALE算法对二维弹性平底结构等速冲击入水过程进行了数值模拟和分析。王易君等[13]应用动网格技术和VOF方法模拟了平底结构的自由入水过程和气-液两相流动,分析了结构物自由落水过程的运动特性、砰击载荷特性及入水空泡演化过程。
结构物冲击入水过程中由于其运动幅度较大,采用传统的网格更新方法(如光顺法、动态层法及网格重构法)计算易出现计算收敛困难、耗时较长及网格负体积等数值计算问题。本文通过重叠网格技术模拟二维平底结构等速冲击入水运动,基于有限体积法与VOF模型求解气液两相流动的雷诺时均方程(RANS),研究平底结构入水过程中砰击载荷特性,并对结构物入水过程中气液两相流动变化情况进行分析。所得结果可为判断水下结构物能否安全穿越飞溅区提供指导。
1 控制方程及重叠网格技术 1.1 控制方程本文采用的控制方程为质量守恒方程及动量守恒方程(N-S方程),并通过VOF模型计算自由面的位置。对于不可压缩的两相流动,将气体与液体密度均视为常数,只考虑水的体积分数的质量守恒方程为:
(1) |
式中:α为水的体积分数,α∈[0, 1],水与空气的体积分数之和为1个单位;u为气-液混合物的速度矢量。
考虑气-液两相为均匀混合,则动量方程可表示为:
(2) |
式中:ρ为气-液混合物密度,ρ=αρw+(1-α)×ρair;ρw为水的密度;ρair为空气的密度;pd为流体动压力,pd=p-ρgx;g与x分别表示重力与位置向量;μ为气-液混合物动力黏度,μ=αμw+(1-α)μair;μw、μair分别为水和空气的动力黏度。
对非稳态的N-S方程做时均化处理,产生了时均化的Reynolds方程,即RANS方程,如式(3)所示。
(3) |
式中:ui与uj为速度分量,i, j分别取1、2、3;
笔者通过引入湍流模型使方程封闭,采用求解RANS方程的方法可以同时满足节约计算成本及得到较高精度流场信息的要求。在此湍流模型选取k-ω SST模型,该模型结合了近壁面边界层区域k-ω模型与高雷诺数区域k-ε模型两者的优势,在提升数值解精度的同时也保证了求解稳定性。
1.2 重叠网格技术与传统的动网格技术不同,重叠网格技术通过简化网格拓扑结构提高了网格适用性。采用重叠网格技术解决具有复杂边界条件的动网格问题(如结构冲击入水等)时,需把整个计算区域分为前景区域及背景区域(至少两个区域),其中前景区域即固体结构位于的子区域。分别对各子区域网格进行离散,利用插值技术实现各区域间流场信息的交互[14]。重叠网格示意图如图 1所示,首先生成背景网格B覆盖全计算区域,并通过前景网格A处理固体部分。在流场计算之前,需对背景网格B进行扫描,标识位于结构内部的网格节点(洞内点)并通过挖洞技术屏蔽此处节点,在之后的CFD计算中这部分洞内点不参与流场计算。节点f、g、h位于固体内部,将此处节点标识为洞内点,同时洞内点周边的e、i节点也需从计算域中除去。在CFD流场计算中,需对各节点进行标识,将全体节点分为有效节点、无效节点及边缘节点,从而进行空间插值。对于网格A,n、s节点为边缘节点,对于网格B,与洞内点相邻的有效节点d、j为边缘节点。
对于每一类网格的边缘单元,应选择另一类网格的相关单元作为贡献单元进行辅助空间插值。对于网格A中的节点n,应选择网格B中的单元b、c作为其贡献单元;对于网格B中的节点j,则选择网格A中的节点r、s作为其贡献网格,以此类推。在流场计算之前,需计算及储存这些重要的节点信息。综上所示,重叠网格建立背景网格与前景网格之间的耦合关系主要为挖洞、寻点及插值3步骤。
2 数值计算模型 2.1 几何模型本文研究水下结构物入水冲击过程中其底部所受的砰击载荷特性,为减小数值模拟的计算量,在不影响计算分析的基础上忽略平底结构宽度方向的影响(默认其为1个单位),将实际的三维结构简化为二维模型。
二维平板、空气域、水域及前景网格区域等几何模型的详细尺寸如图 2所示,流体空间包含了固体结构的全部运动区域。为消除结构物入水冲击时流场壁面边界产生的反射波与固体结构发生耦合作用导致的计算误差,结合无界流场概念[15]及本案例计算的工作需要,设置流场边界尺寸与结构物尺寸之比不小于10:1,结构物模型尺寸与实际工程中水下结构物底部尺寸之比约为1:1。为保证数值计算精确性,前景网格与背景网格应保证足够大的重叠面积,两套网格的网格重叠数至少为4层,且网格尺寸应尽量一致。此处设置前景网格区域尺寸为30 m×20 m。
实际工程中水下结构物穿越水面通常以匀速下放的方式进行,故为降低数值模拟时的计算成本及缩短仿真时间,平板下端面与水面之间无需设置过大距离,此处距离设为1 m。
2.2 数值方法及边界条件通过有限体积法对基于VOF多相流模型的RANS方程进行空间及时间上的离散,采用基于压力的瞬态求解器求解。压力-速度耦合方式采用Couple算法,数值计算过程中压力项离散格式采用PRESTO!,动量方程的离散格式采用二阶迎风格式,湍流方程的离散格式采用一阶迎风格式,通过Geo-Reconstruct对自由液面进行重构。
由于固体结构运动速度较低,所以将水及空气均可视为不可压缩流体,二维平板与流域的下部、左部及右部均设为无滑移壁面,流域上部边界条件为压力出口,其出口压力为标准大气压。
2.3 网格收敛性分析及数值计算验证为保证数值计算网格的收敛性,笔者以平板垂直方向等速(速度1.0 m/s)入水作为分析案件,分别建立3种不同的网格尺寸(0.6、0.5和0.4 m)分析模型,监测平板底面受到的砰击压力随时间的变化情况。图 3为网格尺寸0.4 m时的网格划分图。
图 4为3种不同网格尺寸下结构物等速入水冲击的压力时历图。由图 4可知,关于峰值压力3种网格预测结果基本一致,并且砰击压力与时间的变化关系基本相同,但是随着网格精度的提高,预测精度也逐渐提高,对关键峰值的获取则更加准确。本文考虑数值计算的准确性与可靠性,选择0.4 m的网格尺寸进行以下的数值计算。
为验证本文提出的数值方法的正确性及可靠性,在此模拟1.0 m×0.6 m的二维平板从水平面上方0.3 m处做自由落体运动,并设置流体域的长、宽分别为平板的10倍。监测平板底面中点位置受到的砰击压力随时间的变化情况,结果如图 5所示。在T=0.243 s时刻平板与水平面接触,此时产生峰值压力,将此压力值与相应的试验值[7]对比可知,两者之间的相对误差为8.34%,计算值略小于试验值。另外,与水池试验方法相同,在二维平板底面布置11个监测点,计算平板在撞水时刻其底面受到的砰击压力分布情况,结果如图 6所示。监测点砰击压力以比值的形式给出,即监测点砰击压力与平板中心处砰击压力的比值,pc为平板底面中心处的砰击压力值。由图 6可知,数值计算与水池试验得到的压力分布规律相同。
3 计算结果及讨论 3.1 入水速度对砰击压力的影响
水下生产装备下放时,通常以恒定的速度快速通过飞溅区,根据实际海上作业工况,不同水下结构物穿越飞溅区速度基本保持在0.4~0.8 m/s范围内,不超过1.0 m/s。为了研究工作的可靠性与广泛性,本文分别取下放速度v为0.4、0.6、0.8和1.0 m/s进行研究。不同下放速度入水冲击时平底结构所承受砰击压力随时间的变化曲线如图 7所示。
由图 7可知,随着下放速度的增大,结构受到的最大砰击压力(峰值压力)随之增大,且峰值压力均发生在结构撞水时刻,随后砰击压力逐渐减小并趋于稳定。
以平板底面中点为原点建立一维坐标系,取18个监测点记录平板在撞水时刻底面压力分布曲线,监测点布置如图 8所示。
图 9为不同下放速度下结构物撞水瞬间下端面压力分布情况。由图 9可知,当平板做等速垂向入水冲击时,以结构底部中心为中轴线,左右两侧所受压力对称分布,峰值压力出现的位置为平板底面中点,并从中点向两侧砰击压力逐渐减小,在平板左、右边缘处达到最小值。
3.2 入水横纵速度比对砰击压力的影响
实际海上装备下放时,并不能保证结构物做横向速度为0 m/s的垂直入水运动,本节主要研究结构物下落时横向速度与垂直速度(下放速度v)的比值R对流体砰击压力的影响。选择R分别为0.0、0.2、0.4和0.6进行分析,分析时保持垂向速度(下放速度v)为1.0 m/s不变。经数值计算,得到不同横纵速度比时平板底面砰击压力时历图,如图 10所示。由图 10可知,横纵速度比的增大对结构所受的砰击压力影响较小,4种不同速度比下的砰击压力曲线接近重合,只有峰值压力略有不同。这是因为横向速度的增大使得下降合速度(横向速度与下放速度v的合成速度)增大,砰击峰值压力与入水速度成正比,于是随着横向速度的增加,峰值压力逐渐增大。
与3.1节相同,这里研究不同横纵速度比情况下,结构撞水时刻下端面砰击压力的分布情况,结果如图 11所示。4种不同速度比条件下,平板底面砰击压力分布与图 9相似,随着合速度的增大,砰击压力逐渐增大。但需注意,随着横纵速度比的增大,下端面砰击压力的分布发生了变化,当横向速度不为0时,峰值压力的出现位置偏离平板中轴线,且偏离程度随着速度比的增大而增大。另外,当结构物具有横向平移速度入水时,平板左、右边缘的压力分布不再关于中轴线对称,出现一侧边缘压力增大,另一侧边缘压力减小的现象,且随着速度比值的增大,两侧边缘处的压力差值也逐渐增大。
3.3 入水角度对砰击压力的影响
实际海上作业过程中,结构物受到风浪流的联合作用通常发生绕中轴线的旋转运动,即入水时平板底面与水平面之间具有一定的倾斜角度,本节主要研究入水角度β对结构所受砰击压力的影响。因水下结构物质量较大,入水角度通常保持在较小范围内,这里选取β分别为0°、1°、2°、3°、4°、5°和6°进行比较,其中入水速度均保持垂直方向1.0 m/s(即下放速度v)不变。经数值计算,得到不同入水角度时平板底面砰击压力时历图,如图 12所示。由图 12可知,当固体结构做无倾角的垂直入水时,峰值压力最大,随着入水角度的增大,峰值压力逐渐减小。但应注意,当β大于4°时,此规律将发生变化,此时砰击压力随着β的增大逐渐增大,入水临界倾角为4°。随着入水时间的延长,即固体在水中深度的增加,不同入水角度的平板受到的砰击压力基本收敛于同一数值。
3.4 入水全过程气液两相流动变化情况
固体结构入水全过程是一个具有复杂边界条件的流固耦合问题,在此过程中,固体不仅受到流体的砰击载荷作用,同时固体的运动也会对流体区域产生显著影响,气-液交界面在此过程中会发生复杂变化。本节观察并分析从平板与水接触时刻至被水全包围整个过程气液两相流动的变化情况,流场形态的改变会对入水后结构的运动姿态与运动稳定性产生明显影响。这里选取无水平速度、无入水倾角的等速垂直入水(下放速度v=1.0 m/s)案例进行相关分析。
图 13展示了平板入水全过程中气液两相流体速度的变化规律。平板下落至水面附近过程中压缩空气在其下端面形成“空气泡”,“空气泡”在平板及水域的相互挤压下逐渐破碎,并沿着平板左右两侧往斜上方流动,在平板左、右上角位置形成空气旋涡,平板上方的空气垂直冲击其上端面。从平板上端面与水平面平齐至完全位于水平面之下的过程中(T=3.96~5.00 s时刻),平板下方及其左右两侧被气液两相混合物包围,垂直冲击上端面的空气被溅起的水花向平板中轴线处挤压。在T=6.50 s时刻附近,从平板左右两侧溅起的水流在其中轴线处汇合撞击,使空气沿平板中轴线方向产生垂直向上的射流,并在中部射流与左右两侧水域之间产生空气旋涡。在T=8.00 s时刻附近,左右两侧溅起的水流撞击程度进一步加大,使得水与空气均产生大量垂直向上的射流,此时平板四周的空气较少。在T=10.00 s时刻附近,液体射流在重力作用下落至水平面,并在平板上端面位置形成水流漩涡,此时平板周围逐渐产生气泡。在T=12.60 s时刻附近,水平面基本恢复稳定,产生的气泡逐渐与平板脱离,并在浮力的作用下升向水面。
4 结论
基于VOF模型及重叠网格技术,对二维平底结构等速入水冲击过程中的砰击载荷进行系统分析,数值计算结果与试验数据的对比验证了所提出方法的正确性及可靠性。重点研究了下放速度、入水横纵速度比及入水角度对平板所受砰击压力的影响,分析了不同条件下结构撞水时刻底面砰击压力的分布情况,以及平板入水全过程气液两相流动变化情况,得到如下结论。
(1) 平板垂直等速入水时,结构受到的砰击压力随下放速度的增大而增大,所产生的峰值压力也逐渐增大。在撞水时刻,平板底面砰击压力分布规律为中轴线位置最大,向左右两侧砰击压力逐渐减小,且左右两侧关于中轴线对称分布。
(2) 当平板具有横向移动速度时,不同横纵速度比对砰击压力的影响较小,砰击压力随时间的变化曲线基本相同,仅峰值压力略有差异,平板撞水时刻底面砰击压力不再关于中轴线对称分布,峰值压力的出现位置偏离中轴线,且结构左右边缘处砰击压力的差值随横纵速度比的增大而增大。
(3) 入水角度对平板落水时所受的砰击压力影响较大,当入水角度小于4°时,砰击压力随入水角度的增大而减小,当入水角度大于4°时,砰击压力随着入水角度的增大而增大,入水临界角度为4°。
(4) 平底结构入水全过程中,气液两相流动变化情况较为复杂,主要分为形成空气泡、产生气液射流及水面恢复稳定3个阶段。
[1] |
BUCHNER B, BUNNIK T, HONIG D, et al. A new simulation method for the installation of subsea structures from the splash zone to the ultra deep[M]. Netherlands: [s.n.], 2003.
|
[2] |
WANG Y Y, TUO H H, LI L W, et al. Dynamic simulation of installation of the subsea cluster manifold by drilling pipe in deep water based on OrcaFlex[J]. Journal of Petroleum Science and Engineering, 2018, 163: 67-78. |
[3] |
ZHANG K, DUAN M L, LUO X L, et al. A fuzzy risk matrix method and its application to the installation operation of subsea collet connector[J]. Journal of Loss Prevention in the Process Industries, 2017, 45: 147-159. |
[4] |
ZHANG X H, DUAN M L, MAO D F, et al. A mathematical model of virtual simulation for deepwater installation of subsea production facilities[J]. Ships and Offshore Structures, 2016, 12(2): 182-195. |
[5] |
脱浩虎, 李华, 张赢今, 等. 基于OrcaFlex的深水采油树下放安装过程仿真研究[J]. 石油机械, 2017, 45(3): 58-62. TUO H H, LI H, ZHANG Y J, et al. Simulation of running and installation of deep water X-tree based on OrcaFlex[J]. China Petroleum Machinery, 2017, 45(3): 58-62. |
[6] |
王涛.深水采油树下放安装关键技术研究[D].青岛: 中国石油大学(华东), 2014. WANG T. Research on the key technology of the installation of deepwater subsea tree[D]. Qingdao: China University of Petroleum(Huadong), 2014. http://cdmd.cnki.com.cn/Article/CDMD-10425-1016750593.htm |
[7] |
黄震球, 张文海. 减小平底体砰击的试验研究[J]. 华中工学院学报, 1986, 14(5): 725-730. HUANG Z Q, ZHANG W H. Experimental investigation on the reduction of flat-bottom body slamming[J]. Journal of Huazhong University of Science and Technology, 1986, 14(5): 725-730. |
[8] |
邹丽, 李振浩, 孙铁志, 等. 楔形物体在静水与波浪中自由入水砰击试验研究[J]. 海洋工程, 2019, 37(1): 1-11. ZOU L, LI Z H, SUN T Z, et al. Experimental study on free fall slamming of wedge-shaped objects under still water and waves conditions[J]. The Ocean Engineering, 2019, 37(1): 1-11. |
[9] |
彭晟, 吴卫国, 夏子钰. 江海直达船艏部结构入水砰击试验[J]. 中国舰船研究, 2016, 11(4): 14-21. PENG S, WU W G, XIA Z Y. Experimental investigation on the slamming load of the flat bow of a sea-river link container ship[J]. Chinese Journal of Ship Research, 2016, 11(4): 14-21. |
[10] |
卢炽华, 何友声. 二维弹性结构入水冲击过程中的流固耦合效应[J]. 力学学报, 2000, 32(2): 129-140. LU C H, HE Y S. Coupled analysis of nonlinear interaction between fluid and structure during impact[J]. Acta Mechanica Sinica, 2000, 32(2): 129-140. |
[11] |
孙琦, 周军, 林鹏. 基于LS-DYNA的弹体撞水过程流固耦合动力分析[J]. 系统仿真学报, 2010, 22(6): 1498-1501. SUN Q, ZHOU J, LIN P. Dynamic analysis of fluid-structure interaction for water impact of projectile using LS-DYNA[J]. Journal of System Simulation, 2010, 22(6): 1498-1501. |
[12] |
杨力, 周勇, 李华峰, 等. 基于ALE方法的平底结构入水冲击动力特性研究[J]. 舰船科学技术, 2017, 39(23): 125-133. YANG L, ZHOU Y, LI H F, et al. Research on impact characteristics of water entry of a flat-bottom structure based on ALE method[J]. Ship Science and Technology, 2017, 39(23): 125-133. |
[13] |
王易君, 李明海, 张中礼, 等. 基于VOF法的平底结构自由落体入水砰击载荷模拟[J]. 振动与冲击, 2017, 36(2): 185-189. WANG Y J, LI M H, ZHANG Z L, et al. Numerical simulation on the slamming load in the water-entry process of flatted-bottom body based on the method VOF[J]. Journal of Vibration and Shock, 2017, 36(2): 185-189. |
[14] |
MA Z H, QIAN L, MARTÍNEZ-FERRER P J, et al. An overset mesh based multiphase flow solver for water entry problems[J]. Computers and Fluids, 2018, 172: 689-705. |
[15] |
姜苗苗.组合柱体绕流场模拟与水动力干扰分析[D].哈尔滨: 哈尔滨工程大学, 2017. JIANG M M. Simulation of the flow around combined columns and analysis of hydrodynamic interference[D]. Harbin: Harbin Engineering University, 2017. |