深吃水圆筒型浮式核能平台涡激运动数值模拟

2024-02-26曹辰泽何炎平刘亚东

船舶力学 2024年2期

曹辰泽，何炎平，王梓，刘亚东

（上海交通大学a.海洋工程国家重点实验室；b.船舶海洋与建筑工程学院；c.海洋装备研究院，上海 200240）

0 引言

海上浮式平台已经使用了一个多世纪，它们最常用于石油和天然气开采，但近几年来此类平台被设计用于石油和天然气生产、海上风力发电和海上火箭发射等，部分难以在陆地上进行的工程作业开始越来越依赖海上浮式平台[1]。在图1所示的三种海上浮式平台中，深海Spar平台是一种重心低于浮心的浮式平台，由圆柱形主船体和上部模块两大部分组成；圆筒形FPSO与Spar平台相似，但具有更小的长径比，通常应用于深海和危险的海洋环境；海上浮式核能平台（offshore floating nuclear platform,OFNP）概念由麻省理工学院（Massachusetts Institute of Technology,MIT）首次提出[2]，因与圆筒形FPSO外形类似但吃水更深，也可称为深吃水圆筒型浮式核能平台。该平台创造性地将先进的轻水反应堆技术（light water reactors,LWRs）与海上浮式平台技术相结合，可布置在中国南海岛礁附近并通过相应管线向岛礁输送电能等能源，克服了传统船舶运输补给方式在恶劣海况下难以开展的局限性，可有效解决南海岛礁开发和建设过程中的能源供给问题。

图1 三种海上浮式平台Fig.1 Three types of offshore floating platforms

与深海Spar 平台类似，深吃水圆筒型浮式核能平台的圆筒形主船体结构使得平台在均匀来流作用下，尾流区两侧会形成交替性的漩涡脱落，该现象会导致平台在顺流方向和垂直顺流方向上分别受到周期性的拖曳力和升力作用，进而诱导平台发生周期性的涡激运动，持续性的涡激运动将严重加速系泊和立管系统疲劳损害[3],并对平台内部核反应堆设备的稳定性和安全性带来不利影响。深吃水圆筒形浮式核能平台的概念在国际上提出较晚，国内外对该平台涡激运动方面的研究较少涉及，但已有学者对相似平台进行了一定的涡激运动研究，研究方法主要分为模型试验和数值模拟两类。Fujarra等[4]通过模型试验的方法研究了不同来流角、吃水深度等参数对单柱式生产、储卸油平台（monocolum production,storage and offloading system,MPSO）涡激运动的影响；王颖等[5]结合模型试验和数值模拟两种方法，对TrussSpar平台硬舱部分的尾流流场、涡脱激励及相应的涡激运动特征进行了较为深入的对比研究；何佳伟等[6]采用数值模拟的方法对比分析了有无螺旋侧板Spar平台涡激运动响应的不同。

IDDES 方法是一种基于分离涡模拟（detached eddy simulation,DES）并经过两次改进的湍流模拟方法，该方法解决了分离涡模拟中的模化应力损耗（modeled stress depletion,MSD）现象和延迟分离涡模拟（delayed detached eddy simulation,DDES）中的对数边界层不匹配现象（log layer mismatch），提高了近壁面湍流的求解质量，可以对海上浮式平台涡激运动这类大分离流动问题进行有效的数值模拟。近几年来逐渐有学者采用IDDES 方法对海上浮式平台涡激运动问题展开了数值模拟研究，如胡晓峰等[7]采用IDDES 方法研究了不同流向角、吃水和立柱倒角对深吃水半潜平台纵荡、横荡和艏摇运动响应的影响；杨敏等[8]采用IDDES方法模拟了深海环境下深吃水半潜平台从“预锁定”到“锁定后”范围内折合速度下的涡激运动响应以及平台主尺度参数对涡激运动响应的影响。本文将采用流体力学仿真软件STAR CCM+，首次将IDDES 方法应用于深吃水圆筒型浮式核能平台涡激运动数值模拟研究；在研究平台横荡、纵荡、艏摇运动响应的同时，还将对平台表面及尾流区三维流场细节特征展开具体分析，从而揭示平台涡激运动特性产生的内在机理，研究成果可为平台的结构设计提供参考依据。

1 数值模型概况

数值模拟对象采用MIT设计的OFNP-300[9]，由于在中国南海及大部分真实海域中，兴波引起的波浪力对该平台水平面内涡激运动的影响较小，后续数值模拟忽略自由液面的影响，仅考虑平台吃水线以下的结构，如图2 所示。其中过渡三角体的高度与宽度值均为HA，三个柱靴结构呈120°等间隔布置，采用1:100缩尺比得到数值模拟模型，平台主尺度参数及具体数值见表1。

表1 OFNP-300主尺度Tab.1 Main particulars of OFNP-300

图2 OFNP-300水线面下数值模型Fig.2 Numerical model below waterplane of OFNP-300

2 数值模拟方法

2.1 IDDES方法

在STAR CCM+中进行公式调用可实现IDDES 方法的使用，该方法混合DDES 和壁面大涡模拟（wall modeled LES,WMLES）两个分支[10]，湍流长度尺度定义为

当来流中没有湍流脉动时，IDDES中DDES分支激活，fe= 0，此时湍流长度尺度为

该方法可对来流信息进行准确判断并通过相应湍流长度尺度在近壁区采用RANS 方法求解，在远离壁面区域采用LES方法求解，平台表面及尾流区三维流场细节特征可得到有效捕捉。

2.2 平台三自由度运动方程

研究海上浮式平台涡激运动问题时主要考虑横荡、纵荡、艏摇三个自由度的运动响应，采用STAR CCM+中的动态流固耦合模块（dynamic fluid-body interaction,DFBI）对平台运动进行求解，三个自由度的运动方程分别为

式中，M为平台总质量，Cx、Cy、Cz为平台结构阻尼系数，Kx和Ky为平台有效平移刚度，Kz为有效旋转刚度，FD(t)、FL(t)和Mz(t)分别为作用在平台上的流体拖曳力、流体升力和流体力矩，Iz为平台绕Z轴惯性矩，α(t)为平台绕Z轴旋转角位移。

此外，采用无量纲参数对涡激运动响应结果进行分析，如横荡最大振幅比、横荡标称振幅比、艏摇最大振幅比、艏摇标称振幅比、最大升力系数(CL)max和标称升力系数(CL)nominal，具体表达式分别为

式中，σ为标准差。

2.3 流体域和边界条件设置

为减小侧壁效应对数值模拟结果的影响，流体域设置为长20DH、宽10DH、高3T的长方体区域，如图3 所示。流体域前后两侧分别设置为均匀速度入口边界条件和压力出口边界条件，流体域底部和左右两侧均设置为对称平面边界条件。本文数值模拟忽略自由液面对平台涡激运动的影响，故顶部自由液面设置为对称平面边界条件，保证该面上法向速度和压强梯度为0。

图3 流体域、边界条件和系泊方式示意图Fig.3 Schematic of fluid domain,boundary conditions and mooring

采用三根水平线性弹簧来等效实际系泊方案，三根处于相同预拉伸状态的等长弹簧间隔120°布置并将刚度均设置为10.3 N/m。弹簧一端连接平台表面导缆孔，另一端约束在锚泊点上，三个连接点的垂向坐标值均与平台质心高度值保持一致。

2.4 重叠网格法

基于重叠网格方法，将复杂的流体域分成背景区域和重叠区域，如图4所示。各区域中的网格均由切割体网格生成器独立生成，重叠区域中平台表面处的边界层网格则采用棱柱层网格生成器生成，流场信息在重叠边界处通过插值方法进行匹配和耦合，并采用“挖洞”方式对重叠区域内不参加流场计算的网格单元进行删除。该方法容易设置，并且允许平台发生更加灵活的运动，例如六自由度运动[11]。

图4 流体域、平台表面网格与重叠网格法Fig.4 Grids of platform surface,fluid domain and overlapping grid method

为了更好地捕捉平台尾流区三维流场特性，对整个流体域采用四级网格加密，其中加密区域1、2、3 均为长方体，加密区域4 设置为圆柱体，四个加密区域的高度均设置为1.4T。平台表面边界层网格在划分时，棱柱层设置为15层，棱柱层延伸率为1.3，保证y+值小于1。

2.5 数值模拟方法验证

深吃水圆筒型浮式核能平台属于比较新的概念，国内外在该平台涡激运动研究方向上暂无公开文献发表，不过由于该平台主体部分类似Spar型结构[12]，为验证本文数值模拟方法可行性，对文献[13]中一小型Spar 模型的涡激运动响应开展数值模拟并以横荡运动振幅比作为验证参数，验证结果如图5所示。分析可知，在本文研究的折合速度范围内，数值模拟结果同模型试验结果吻合较好，证明本文数值模拟方法具备一定可行性。但在折合速度为8、9时，数值模拟结果相比模型试验有不同程度的偏大，这是由于数值模拟没有考虑模型试验中弹簧系统的阻尼。

图5 试验结果[13]与本文数值模拟结果对比Fig.5 Comparison of CFD simulation results with experimental results

3 数值模拟结果与分析

3.1 自由衰减试验与网格无关性验证

为验证网格无关性，采用粗糙、基础、精细三套网格进行平台静水自由衰减试验，设置均匀来流速度为0，在水平弹性系泊条件下给定平台一个横向的初速度0.14 m/s 并释放，开放平台的横荡和艏摇运动自由度，无量纲时间步长（UΔt/DH）设置为0.016（Δt为时间步长），结果见表2。以基础网格算例作为参照，结果表明用粗糙网格和精细网格算例计算的横荡运动固有周期Tny、艏摇运动固有周期Tnz与用基础网格算例计算结果的相对误差绝对值整体在0.5%以内，因此可以认定网格粗细不会对本文数值模拟结果产生较大影响。为保证计算精度兼顾计算效率，后续数值模拟将采用基础网格进行。

表2 网格无关性验证Tab.2 Grid independence study

图6给出了平台横荡和艏摇静水自由衰减曲线，计算得到平台模型尺度下的横荡、艏摇运动固有周期分别为20.422 s 和6.268 s，根据傅汝德数相似准则换算得到平台实尺度下的横荡、艏摇运动固有周期分别为204.22 s和62.68 s。

图6 OFNP-300自由衰减试验Fig.6 Hydrostatic damping test of OFNP-300

3.2 横荡、纵荡运动时历曲线和水平面内运动轨迹

本文数值模拟采用中国南海海域的海流环境，其最大流速不超过2.5 m/s[14]，表3 给出实尺度下六组工况的流速以及模型尺度下对应的来流速度、折合速度和雷诺数。

表3 流速工况Tab.3 Inlet flow conditions

图7给出OFNP-300横荡、纵荡位移时程图和水平面内质心运动轨迹，当Ur= 4.54时，平台几乎不发生纵荡且横荡振幅较小，运动轨迹未呈现明显特征；当5.45 ＜Ur＜7.26时，纵荡、横荡振幅均逐渐增加，纵荡运动频率近似为横荡运动的两倍，运动轨迹类似“香蕉”形；当7.26 ＜Ur＜9.08 时，纵荡、横荡振幅进一步增加，纵荡运动呈现出强烈的不规则性，运动轨迹在顺流方向上变宽。综合来看，平台纵荡振幅在整个折合速度范围内保持在较小范围，后文流场分析将给出进一步解释。

图7 不同折合速度下横荡、纵荡位移时程图和水平面内质心运动轨迹Fig.7 Time history diagrams of sway and surge and the motion trajectory

3.3 横荡、艏摇运动和升力系数幅值

图8给出了横荡、艏摇以及升力系数响应幅值随折合速度的变化曲线，其中横荡幅值在整个折合速度范围内逐渐增加，当7.26 ＜Ur＜9.08时，横荡幅值增加速度开始放缓并趋于0，运动表现出一定的“锁定”特性，但横荡运动是否具有“锁定区间”需结合后文频率分析再给出定论；艏摇运动中，当7.26 ＜Ur＜9.08 时艏摇运动响应幅值与折合速度表现出近似线性递增的关系，呈现一定的驰振特征，这一现象不同于典型的涡激运动特性；升力系数在整个折合速度范围内呈现出先增加再减小的变化趋势，当Ur= 7.26 时升力系数幅值达到最大，当Ur＞7.26 时升力系数突然呈现下降趋势，这可能与该折合速度下三维流场特性的改变有关。

图8 不同折合速度下横荡、艏摇和升力系数幅值Fig.8 Amplitudes of sway,yaw and lift coefficients at different reduced speeds

3.4 横荡、艏摇运动响应频率

对OFNP-300 在不同折合速度下的横荡运动和艏摇运动响应进行快速傅里叶变换获得频率随折合速度变化的三维曲线，如图9 所示，fv、fn分别为平台横荡（艏摇）运动频率和固有频率。分析可知，平台横荡运动频率一直保持在一个较为稳定且接近于横荡固有频率的区间，不同于传统Spar 平台横荡运动先上升后锁定再脱离的典型“锁定”现象。艏摇运动频率变化曲线表明，在Ur＞7.26时，主频附近出现了多个峰值并且峰值对应的频率与艏摇运动固有频率接近。由此可见，在高折合速度下，尾流区漩涡的泄放频率不再单一。对比横荡和艏摇运动频率变化可以看出，横向运动频率峰值和艏摇运动频率峰值一致，表明艏摇运动也受交替泄涡现象激励。

3.5 流场水平截面速度云图

同时水平截取z= 2HS和z= 0.5HS处的流场速度云图，如图10所示。分析发现平台主体结构表面能明显观察到流动分离现象，而柱靴结构因为“槽”的存在，其表面流动分离形态受到较大干扰。由此可见，平台主体结构与柱靴结构所在平面内的流场特性存在较大差别，这也验证了后文研究流场三维特性的必要性。

图10 水平截面速度云图Fig.10 Horizontal section velocity contours

3.6 平台受力与流场三维特性分析

通过观察平台表面压力分布情况可大致判断出边界层分离点的位置，如图11所示。柱靴表面边界层分离点位置大致分布于“槽”内靠近背流面一侧，迎流面相比于背流面具有更大范围的高压力区域，而且槽内壁面也存在一定范围的高压力区，这样的压力分布规律对平台的结构设计而言具有一定的工程意义。

图11 平台表面压力（Ur = 7.26）Fig.11 Pressure on platform surface(Ur = 7.26)

对三维流场展开研究更加能体现出流场的复杂性和不稳定性，研究平台附近三维流线结构和三维垂向旋涡结构具有重要意义。从图12 中可以发现柱靴结构打乱了平台表面处的流线，“槽”内流线结构非常复杂，柱靴表面规则泄涡过程被破坏，结合前文z= 0.5HS处水平截面速度云图结果可推测柱靴结构具有一定减涡效果。图13 中，三维漩涡结构用Q准则表示，用来流速度在顺流方向上的分量着色，可观察到尾流区漩涡泄放形态十分复杂，高强度漩涡结构主要分布在平台背流面附近；尾流区不同方向上的漩涡结构相互作用，相互牵制，导致平台表面受力具有较大的不规律性，解释了前文艏摇运动频率在主频附近出现了多个峰值这一现象。