0.   引　言

随着海洋开发逐渐向深海领域发展，海洋环境更加恶劣，在大波幅波浪作用下，海洋结构物载荷与运动的非线性特性显著增强。目前，常用的商用软件（例如AQWA，Seasam，Hydrostar，Wamit等）都是基于线性或弱非线性的假定进行数值模拟，这必然会导致计算精度不足。因此，针对波浪和海洋结构物相互作用的问题，开展完全非线性数值模拟研究，对于提高海洋工程装备的安全性具有重要意义。

1976年，Longuet-Higgins等^[1]基于势流理论提出了一种完全非线性方法，采用混合欧拉—拉格朗日方法追踪自由水面质点，自由水面和物面边界条件得到实时满足，每一时刻都对网格进行重构。近年来，随着计算机运算能力的不断提高，完全非线性方法得到了快速发展。周斌珍^[2]采用完全非线性方法建立了开敞水域模型，并对漂浮圆柱的受迫振动、绕射问题以及Ringing现象进行了深入研究，揭示了非线性作用的机理。Bai等^[3-6]基于高阶边界元法（HOBEM），运用区域分解技术，建立了波浪水池的完全非线性三维数值模型，并对圆柱的受迫振动以及非线性波浪作用下圆柱的运动特性及波浪爬高特性等进行了分析研究。宁德志^[7]采用快速多极子去奇异边界元法（Fast Multi-pole Boundary Element Method）分别建立了非线性三维开敞水域模型及非线性数值波浪水槽模型，并对圆球的强迫纵荡、垂荡和转动问题以及波浪对圆柱的绕射等问题进行了探讨。Feng等^[8-9]利用完全非线性数值模型对并列双箱在波浪作用下的水动力共振问题进行了研究。然而，上述研究都仅局限于简单几何形状物体，鲜有考虑船舶这种带有复杂曲面形状的海洋结构物。

本文将针对某浮式生产储卸油轮（Floating Production Storage and Offloading，FPSO）模型进行数值模拟研究。首先，综合考虑自由水面、物面以及入射波浪的非线性影响，基于势流理论，结合高阶边界元法，建立波浪与海洋结构物相互作用的三维完全非线性开敞水域模型。然后，采用混合欧拉—拉格朗日（MEL）方法追踪瞬时自由水面，利用四阶Runge-Kutta法预报更新下一时间步的水面、物体位移和速度，每一时刻都对自由水面和物体湿表面网格进行更新重构。最后，通过对自由漂浮的圆柱、无航速的Wigley船模型进行模拟，验证所提模型的准确性，并模拟某FPSO模型的水动力特性。

1.   数学模型

1.1   初边值问题

图1所示为建立的完全非线性开敞水域模型。图中：S_F为瞬时自由水面边界；S_B为瞬时物体湿表面；S_D为水底边界；d为水深；n为物面单元单位法向矢量，以指出流体方向为正方向。定义2组坐标系，分别描述流体域和物体的运动。一个是固定坐标系O-xyz，其中z = 0位于静水面上，垂直向上为正；另一个是船体坐标系O'-x'y'z'，其中原点O'在船体质心上。当物体静止处于平衡位置时，2组坐标系相互重合。

图  1  坐标系和计算域定义

Figure  1.  Definition of coordinate system and computational domain

下载: 全尺寸图片 幻灯片

假定流体为无黏性、不可压缩、无旋运动的理想流体，则整个流域可用速度势函数描述，总的流场速度势ϕ满足Laplace方程：

由物质导数和伯努利方程推导可得瞬时自由水面的完全非线性运动学和动力学边界条件。本文采用混合欧拉—拉格朗日方法追踪瞬时自由水面，其边界条件在拉格朗日系统下可以写为

式中：${\rm{d}}/{\rm{d}}t = {\rm{\partial}} /{\rm{\partial}} t + \nabla \phi \cdot \nabla$，为物质导数；X(x, y, z)为瞬时自由水面流体任意质点的位置矢量；η为自由水面高程；g为重力加速度。

在瞬时物体湿表面S_B上，满足如下边界条件：

式中，V为物面单元速度矢量。针对固定边界，满足V = 0；针对运动的物体，V可以表示为

式中：U为结构物的平动速度；Ω为结构物的转动速度；r_b为位置矢量。

1.2   速度势分离

本文采用Ferrant等^[10-11]提出的入射势和散射势分离技术，将总波浪$\left( {\phi ,\eta } \right)$分解为入射波浪$\left( {{\phi _{\rm{I}}},{\eta _{\rm{I}}}} \right)$和散射波浪$\left( {{\phi _{\rm{S}}},{\eta _{\rm{S}}}} \right)$，即

式中的入射波浪通过高阶Stokes理论解给定。将式(6)代入式(1)~式(4)，整理可得关于散射波浪$\left( {{\phi _{\rm{S}}},{\eta _{\rm{S}}}} \right)$的初边值问题。

式中：r为流体质点与空间固定坐标系原点O的水平距离；x₀和y₀为初始时刻流体质点的水平坐标；μ(r)为阻尼系数函数，其作用是对流体出口处的散射波浪进行消波，以防止波浪反射对结构物造成二次作用，具体表示为

式中：ω为波浪圆频率；λ为入射波波长；r₀为阻尼层起点半径；α₀和β₀为关系数，本文中，α₀=β₀=1。

为了避免初始效应导致数值不稳定，干扰模拟进程，本文在与(ϕ_I, η_I)有关的项上作用一个缓冲函数：

式中，T_m为缓冲时间，一般取入射波周期的整数倍，这里取为2倍。

时域分析理论是一种完全模拟实际作用过程的理论，其控制方程中不仅含有空间项，还含有时间项，因此还需要给出初始条件：

1.3   边界积分方程的建立

在整个计算域内应用格林第二定理，将上述关于散射波浪(ϕ_S, η_S)的边界值问题转化成如下的边界积分方程进行求解：

式中：p= (x₀, y₀, z₀)，为源点坐标；q= (x, y, z)，为场点坐标；α(p)为固角系数；S₁为整个计算域表面；G(p, q)为格林函数。

由于整个计算域是关于O-xz平面对称，且海底是水平的，所以选取简单Rankine源及其关于对称面（y = 0）和海底（z = −d）的镜像为格林函数。鉴此，在计算时只需考虑一半计算域，而且海底被排除在外，式（17）则为G ( p , q )格林函数。

其中，

这里，采用HOBEM法，将计算域离散为8节点的四边形单元，并将其表示为参数坐标下的等参单元。在每个单元内，则采用二次形函数来描述任意点的几何坐标、速度势和速度势的法向导数，即

式中：$\left( {\xi ,\varsigma } \right)$为参数坐标；X_k，ϕ_k，${\left( {\partial \phi /\partial n} \right)_k}$，h_k分别为单元上第k个节点的坐标、速度势、速度势法向导数和形函数；K为单元内节点的个数。将式(19)代入式(16)中，积分方程可以离散成如下形式：

式中：$J\left( {\xi ,\varsigma } \right)$为雅各比行列式；N_e1和N_e2分别为自由水面和物面边界上离散的曲面单元数。上述积分采用标准的Gauss-Legendre积分，将8节点四边形单元离散为4×4个样本点进行求解。最终积分方程转换为如下线性方程组：

式中：A₁₁，A₁₂，A₂₁，A₂₂均为系数矩阵A的元素；B₁和B₂均为向量B的元素。上述积分边界每一时刻都在变化，因此每一时间步的系数矩阵A和向量B都需重新建立再求解。计算过程中，假定当前时刻物面上的速度势导数和自由水面上的速度势均已知，由积分方程(20)计算当前时刻的物面速度势和自由水面上的速度法向导数后，采用四阶Runge-Kutta法，基于自由水面边界条件式(8)~式(11)，计算下一时刻水质点位置的和速度势。自由水面流体质点的速度分量可由式(22)计算：

式中：n_x，n_y和n_z为法向向量n沿x，y和z方向的分量；$\partial {\phi _{\rm{S}}}/\partial \xi$和$\partial {\phi _{\rm{S}}}/\partial \varsigma$可直接通过式(19)求解。需要注意的是，一个节点可能由几个单元共用，因此需对节点周围单元速度取平均值，得到该节点的速度分量。

2.   数值实现

因为完全非线性方法的自由水面和物面边界条件都均实时满足，所以每个时间步都要重新捕捉及更新自由水面和物体湿表面网格。自由水面和物面的网格重构是完全非线性模型工作的重点及难点。

2.1   自由水面网格更新

本文采用弹簧近似法^[12-13]对自由水面动网格进行处理。弹簧近似法的本质是把网格节点看作由弹簧连接起来的，整个网格像弹簧系统一样变形。因此，只需要在初始时刻形成一次网格即可。网格节点位移可由式(23)确定：

式中：δ_i和δ_j分别为节点i和j的水平位移矢量；N_i为网格节点总数；k_ij为节点i和j之间的弹簧刚度。k_ij的计算式为

式中，${l_{ij}}$为节点i与j之间的距离。更新后，网格节点i的水平位置矢量${{X}}_i^{{\rm{new}}}$可表示为

式中，${{X}}_i^{{\rm{origin}}}$为节点i在初始网格中的水平坐标。

至此，新的均匀网格节点位置已经完成更新。下面将通过二次形函数插值求解该网格节点上的速度势以及波面高程。为方便描述，将通过上述方法求得的网格记为网格1，将通过混合欧拉—拉格朗日方法（MEL）方法求得的网格记为网格2。为进行插值求解，需首先确定网格1中的节点在网格2中所属的单元，可以通过式(26)来确定：

式中：S_m为网格2中第m个单元的面积；S_j为三角形子元素的面积；该三角形由网格1中所求节点和第m个单元的节点构成；M₁为三角形子元素的个数。当式(26)成立时，网格1所求节点一定在第m个单元的内部，此时可以通过式(27)来求解节点的参数坐标$\left( {\xi ,\varsigma } \right)$：

式中：(x_i, y_i)为网格1上所求节点的坐标；(x_k^m, y_k^m)为网格2上第m个单元节点的坐标。

通过上述处理得到的自由水面的计算域网格相对均匀，可避免因网格变形、扭曲而导致的数值不稳定性。

2.2   物面网格更新

对于船舶这种具有复杂曲面的海洋结构物来说，物面网格的更新较为困难。为了降低网格更新难度及工作量，本文将物面网格分为靠近水线区域和远离水线区域，如图2所示。

图  2  物面区域划分示意图

Figure  2.  Schematic diagram of body surface domain division

下载: 全尺寸图片 幻灯片

因为远离水线区域始终处于瞬时水线以下，所以这部分网格可以不用更新。针对靠近水线区域的网格，通过背景网格插值技术予以更新。初始时刻通过读入船体各站的型线数据进行网格划分，包括水线以上的干舷部分和水线以下的部分，其中水线以上干舷部分的网格只用来更新物面网格，不参与计算。

物面网格均在船体坐标系下更新。假设船体x坐标不变，根据更新后自由水面的水线位置对船体z坐标进行重新划分，y坐标通过背景网格插值来求解，即：将船体网格投影到O-xz平面，通过式(26)确定当前所求节点在背景网格中所在的单元m。当式(26)成立时，所求节点一定在第m个单元的内部。此时，首先通过式(27)求得所求节点的参数坐标$\left( {\xi ,\varsigma } \right)$，然后利用二次形函数插值求得节点的y坐标。当所有节点的坐标均完成更新后，需要将船体坐标系下的物面节点坐标转换到固定坐标系下，2组坐标系下的坐标满足如下转换关系：

其中，

式中：(x, y, z)为固定坐标系下的坐标；(x', y', z')为船体坐标系下的坐标；(x_c0, y_c0, z_c0)为物体的质心位置；($\zeta_1$, $\zeta_2$, $\zeta_3$ )为物体沿坐标轴的平动位移；(θ, β, γ)为物体绕坐标轴的转动角度。

3.   数值计算及讨论

3.1   模型验证

3.1.1   漂浮圆柱模型

为验证本文模型的准确性，对文献[10]中给出的漂浮圆柱模型进行模拟。参数设置如下：水深d=1 m，圆柱半径R/d=1/3，吃水深度D/d=1/6，波幅A/d=0.05，阻尼阵为0，波浪周期$T\sqrt {g/d}=7$。入射波浪采用Stokes 5阶波。然后，将本文所建完全非线性模型的模拟结果与文献[10]的模拟结果进行对比，结果如图3所示。包括图3及以下各图结果经过无量纲化处理，以便对比。

图  3  垂向位移、速度和加速度的时历曲线

Figure  3.  Time history curves of vertical displacement, velocity and acceleration

下载: 全尺寸图片 幻灯片

由图3可见，两者模拟结果吻合良好，验证了本文模型的准确性。另外，由垂向速度时间历程曲线可以看出，本文结果中的二阶波峰高度并不总是相等，存在能量传递的过程。这可能是因为本研究采用了混合欧拉—拉格朗日方法来追踪瞬时流体质点，相较于文献[10]使用的半混合欧拉—拉格朗日方法，流体质点不仅有垂向运动，还有水平运动，其水平运动对垂向运动产生了一定影响，进而对物体的运动特性也有所影响。

分别对垂向位移、速度和加速度的时间历程曲线进行傅里叶分频处理，可以得到垂向位移、速度以及加速度的幅值谱，如图4所示。从图中可以看出：物体的运动除了主频之外还有很大的二阶成分，加速度甚至还有其他更高阶的成分，因此物体在非线性波浪中的运动具有很强的非线性效应。而实际海浪大多为非线性波浪，因此，采用完全非线性模型对海洋结构物的运动特性进行模拟，有利于提高海洋结构物的安全性。

图  4  垂向位移、速度和加速度幅值谱

Figure  4.  Amplitude spectrum of vertical displacement, velocity and acceleration

下载: 全尺寸图片 幻灯片

3.1.2   零航速Wigley船

为验证本文模型对不同形状物体的模拟效果，对文献[14]中采用的无航速Wigley-Ⅲ船模型进行模拟，并将模拟结果与该文献的实验结果以及文献[15]的数值计算结果进行对比。

该模型的船长L=3 m，船宽B=0.3 m，吃水深度D=0.187 5 m，重心距船底0.17 m，纵向惯性半径R_yy=0.75 m，排水量0.078 m³。船型函数由式(30)确定：

本文只考虑Wigley船的升沉和纵摇2个自由度运动，其他方向都受到限制。

图5所示为Wigley船的网格划分示意图，基于船体对称性，将Wigley船的一半划分为14×6个单元。

图  5  Wigley船网格划分示意图

Figure  5.  Schematic diagram of mesh generation on Wigley hull

下载: 全尺寸图片 幻灯片

图6是Wigley船的升沉运动和纵摇运动的振幅响应（RAO）结果，该结果经过无量纲化处理。从图中可以看出：1）本文模型模拟结果与文献[14]实验结果以及文献[15]数值计算结果基本吻合；2）在入射波波长较小时，本文的数值计算结果没有出现实验结果中的波动现象，这可能由模型实验的精度误差或者是流体黏性导致；3）当$\lambda /L > 2.0$时无实验结果，这可能是因为当入射波波长较长时，物理模型实验中的波浪反射现象较为严重，对实验结果有较大影响。

图  6  Wigley船运动响应对比

Figure  6.  Comparison of motion responses for Wigley hull

下载: 全尺寸图片 幻灯片

模型验证结果表明，本文的模型对不同形状的物体均能得到较好的模拟结果。

3.2   零航速FPSO水动力特性研究

利用本文所建完全非线性模型对某FPSO模型的水动力特性进行模拟，将本文结果进行无量纲化处理，并与AQWA计算结果以及文献[16]的实验结果进行对比。其中，本文实验是在水池中完成的，水池的有效使用尺寸为40 m，最大工作水深为0.7 m，配有实验室自制的可移动式多向不规则波造波机以及微机控制的数据采集系统。

表1所示为FPSO的主尺度参数，为与实验结果^[16]进行对比，选取缩尺比为1∶45，本文对7.5和9 m这2种吃水深度的模型尺度进行模拟计算，网格划分示意图如图7所示。一般情况下，船舶在顶浪或者顺浪时的纵摇和升沉运动幅度最为剧烈，因此，本文主要研究FPSO在顶浪和顺浪时的升沉及纵摇运动响应特性。根据对称性，在一半的FPSO上划分了175个单元，来浪方向为0°（顺浪）和180°（顶浪），所模拟的波浪频率为0.35~1.09 Hz。

表  1  FPSO主尺度参数

Table  1.  Particulars of the FPSO

参数	9 m吃水深度			7.5 m吃水深度
	实尺度模型	缩比模型		实尺度模型	缩比模型
总长/m	112	2.48		112	2.48
垂线间长/m	109	2.42		109	2.42
船宽/m	28	0.62		28	0.62
型深/m	16	0.36		16	0.36
排水量/t	24 610	0.263		20 347	0.218
重心纵向位置/m （相对于船艉 ）	51	1.133		51.2	1.138
重心垂向位置/m （相对于船基线）	7.2	0.16		8.2	0.182
横摇惯性半径/m	8.037	0.179		7.9	0.176
纵摇惯性半径/m	28.149	0.626		27.954	0.621
艏摇惯性半径/m	29.027	0.646		28.905	0.642

下载: 导出CSV 

| 显示表格

图  7  FPSO网格划分示意图

Figure  7.  Schematic diagram of the mesh generation on FPSO

下载: 全尺寸图片 幻灯片

图8和图9分别为7.5 m吃水深度时FPSO在顺浪和顶浪时的运动时历曲线。从图中可以看出，在短波情况下（0.90 Hz），其升沉和纵摇运动响应幅值都较小；在中等波长情况下（0.70 Hz），升沉运动幅值并未达到最大值，但是纵摇固有周期约等于波浪对船体的扰动周期，并发生谐摇，而纵摇幅值达到最大；在最大波长情况下（0.35 Hz），其升沉运动幅值比中等波长时大，纵摇运动幅值比中等波长时小，此时已经处于随波逐流状态。

图  8  波浪0°（顺浪）入射工况下7.5 m吃水深度时FPSO的运动响应时历曲线

Figure  8.  Time histories of motion responses of FPSO with 7.5 m draft and 0° incident angle wave

下载: 全尺寸图片 幻灯片

图  9  波浪180°（顶浪）入射工况下7.5 m吃水深度时FPSO的运动响应时历曲线

Figure  9.  Time histories of motion responses of FPSO with 7.5 m draft and 180° incident angle wave

下载: 全尺寸图片 幻灯片

由此可见，针对升沉和纵摇运动，波浪对船体扰动的大小与波长和船长的比值有较大关系：波浪频率越大，升沉和纵摇运动越缓和，故大船在小波浪中的升沉和纵摇运动幅度都不会很大；最大纵摇角度出现在中等波长下，此时波长略大于船长，最大波长的纵摇幅值介于小波长和中等波长之间。

图10和图11分别为9 m吃水深度时FPSO在顺浪和顶浪时的运动时历曲线。从图中可以看出，其运动规律与7.5 m吃水深度时的运动规律保持一致，表明吃水深度的改变对FPSO的运动特性影响较小。

图  10  波浪0°（顺浪）入射工况下9 m吃水深度时FPSO的运动响应时历曲线

Figure  10.  Time histories of motion responses of FPSO with 9 m draft and 0° incident angle wave

下载: 全尺寸图片 幻灯片

图  11  波浪180°（顶浪）入射工况下9 m吃水深度时FPSO的运动响应时历曲线

Figure  11.  Time histories of motion responses of FPSO with 9 m draft and 180° incident angle wave

下载: 全尺寸图片 幻灯片

图12和图13分别为顺浪和顶浪情况下FPSO的升沉运动响应结果。由图可以看出，本文模拟结果与AQWA计算结果以及文献[16]的实验结果趋势一致，升沉运动响应随着入射波浪周期的增大而增大。在顺浪时，不论是长周期波还是短周期波，本文模拟结果均与AQWA计算结果吻合良好；但是在顶浪时，本文模拟结果与AQWA计算结果存在一定差异，尤其是对高频波浪。从图中可以看到：在高频波浪区段，本文模拟结果仍出现了局部小峰值，与文献[16]的实验结果吻合得更好，而AQWA计算结果未能体现这一现象，表明本文完全非线性模型的模拟结果更加符合实际。

图  12  波浪0°（顺浪）入射工况下FPSO升沉RAO

Figure  12.  The heave RAO of FPSO with 0° incident angle wave

下载: 全尺寸图片 幻灯片

图  13  波浪180°（顶浪）入射工况下FPSO升沉RAO

Figure  13.  The heave RAO of FPSO with 180° incident angle wave

下载: 全尺寸图片 幻灯片

图  14  波浪0°（顺浪）入射工况下FPSO纵摇RAO

Figure  14.  The pitch RAO of FPSO with 0° incident angle wave

下载: 全尺寸图片 幻灯片

图14和图15分别为顺浪和顶浪情况下FPSO的纵摇运动响应结果，从图中可以看到，本文模拟结果与AQWA计算结果以及文献[16]实验结果的趋势一致，FPSO纵摇运动响应随着入射波浪周期的增大先增大，达到谐摇状态，而后开始减小。从图中还可以看到：本文完全非线性模型在共振区段的模拟结果要大于AQWA计算结果，与文献[16]的实验结果吻合得更好，且共振区段的频带宽度相较于AQWA计算结果变宽了，与文献[16]的实验结果所表现出来的现象一致。这是因为当波浪的频率与FPSO的固有振荡频率接近时，FPSO会产生共振，导致运动的非线性效应加强；另外由于考虑了非线性效应，当入射波浪与FPSO相互作用时，可能会诱发产生其他频率的波浪，从而导致了共振区的频带宽度变大。传统的线性方法无法考虑这些非线性因素的影响，对海洋结构物运动的预报偏于危险，不利于海洋结构物的安全性。此外，在低频区段，本文模拟结果以及AQWA计算结果和文献[16]的实验结果相差较大，这可能是因为在水池的有限范围内进行模型实验时，难以完全忽略低频波长波浪反射等因素对实验的影响（模型实验水池长40 m，最大低频波长约为7 m。虽然水池设有斜坡和消浪材料，仍难以完全消除反射波的影响），而数值模拟方法可以很好地排除这些因素对结果的影响。

图  15  波浪180°（顶浪）入射工况下FPSO纵摇RAO

Figure  15.  The pitch RAO of FPSO with 180° incident angle wave

下载: 全尺寸图片 幻灯片

表2给出了不同吃水深度和不同波浪入射角时FPSO的纵摇最大响应周期，可以看到本文数值结果与AQWA计算结果以及文献[16]的实验结果均有所差异，但误差都在可接受的范围内，表明本文模型能够得到较为精确的结果。

表  2  不同入射波时纵摇最大响应周期

Table  2.  Maximum response period of pitch with different incident angle wave

吃水 深度/m	方法	最大响应周期
		0°	180°
7.5	文献[16]	3.195	3.195
	AQWA	2.876	2.876
	本文模型	3.050	3.050
9	文献[16]	3.195	3.195
	AQWA	2.796	2.796
	本文模型	2.876	2.876

下载: 导出CSV 

| 显示表格

图16和图17分别为顺浪和顶浪工况下FPSO的横摇运动响应结果，文献[16]的实验结果与其他2种计算结果出现了较大偏差：模型实验中出现了一定的横摇角，而其他结果中基本没有横摇响应。针对此问题，大多数学者^[17]认为是模型实验产生的误差，主要原因有：1）模型难以保证绝对的左右对称；2）波浪由造波机产生，在较宽的水池中并排分布多块造波板，造波板之间可能产生微小的相位偏差；3）其他可能的系统误差或人为误差。本文采用开敞水域模型，入射波浪直接由解析解给定，并采用对称模型计算，避免了产生上述误差。

图  16  波浪0°（顺浪）入射工况下FPSO横摇RAO

Figure  16.  The roll RAO of FPSO with 0° incident angle wave

下载: 全尺寸图片 幻灯片

图  17  波浪180°（顶浪）入射工况下FPSO横摇RAO

Figure  17.  The roll RAO of FPSO with 180° incident angle wave

下载: 全尺寸图片 幻灯片

图18和图19为7.5 m吃水深度、周期t/T=5时，不同波浪频率下FPSO附近的流场图。从图中可以看到，远场区域自由水面相对规则，表明本文使用的人工阻尼层消波方法效果良好。FPSO对短波产生了较大的扰动，这是因为船长大于波长，FPSO易对短波的传播产生影响。由前面的分析可知：最大纵摇角度出现在中等波长下（0.70 Hz），此时波长略大于船长，从图中可以更直观地得到该结论。由图还可以看出，在顶浪时FPSO周围的扰动相较于顺浪时要更剧烈一些，会对船体的运动产生一定影响，非线性效应会更强一些，这与图12~图15所呈现的现象是一致的。

图  18  波浪0°（顺浪）入射工况下7.5 m吃水深度时FPSO附近流场图(t/T=5)

Figure  18.  Flow field near the FPSO with 7.5 m draft and 0° incident angle wave (t/T=5)

下载: 全尺寸图片 幻灯片

图  19  波浪180°（顶浪）入射工况下7.5 m吃水深度时FPSO附近流场图(t/T=5)

Figure  19.  Flow field near the FPSO with 7.5 m draft and 180° incident angle wave (t/T=5)

下载: 全尺寸图片 幻灯片

图20和图21为9 m吃水深度、周期t/T=5时，不同波浪频率下FPSO附近的流场图。从图中可以看到，其流场特性与7.5 m吃水深度时的流场特性相似，再次表明吃水深度对FPSO的运动特性及流场特性影响较小。

图  20  波浪0°（顺浪）入射工况下9 m吃水深度时FPSO附近流场图(t/T=5)

Figure  20.  Flow field near the FPSO with 9 m draft and 0° incident angle wave (t/T=5)

下载: 全尺寸图片 幻灯片

图  21  波浪180°（顶浪）入射工况下9 m吃水深度时FPSO附近流场图( t/T =5)

Figure  21.  Flow field near the FPSO with 9 m draft and 180° incident angle wave (t/T=5)

下载: 全尺寸图片 幻灯片

4.   结　论

本文利用高阶边界元（HOBEM）方法建立了波浪与海洋结构物相互作用的时域三维完全非线性开敞水域模型，采用混合欧拉—拉格朗日方法追踪瞬时自由水面，利用四阶Runge-Kutta法更新下一时间步的自由水面、物体位移和速度，且每一时刻都对自由水面和物体湿表面网格进行更新。基于完全非线性模型的数值模拟结果，得到如下结论。

1) 通过对非线性波浪作用下自由漂浮的圆柱的模拟，验证了本文模型的准确性，证明了物体在非线性波浪中的运动具有很强的非线性。因此，在对物体运动进行预报时，采用线性波代替非线性波不利于物体的安全性。

2) 通过对Wigley船和某FPSO模型的模拟，表明了本文模型对不同形状的复杂几何物体均有较好的模拟结果，表明了本文模型在实际工程中应用的可行性。

3) 通过对某FPSO模型在顺浪和顶浪中升沉及纵摇运动响应的模拟，结果表明，本文模拟结果与AQWA软件计算结果以及文献实验所得到的运动特性基本吻合，这进一步验证了本文模型的准确性。此外，在考虑非线性效应的情况下，本文模型的升沉运动受到的影响较小，基本上仍与线性理论所得结果一致。纵摇运动表现出明显的非线性特性：共振区段的幅值明显增大且频带宽度变大，这与实验结果的现象一致，说明传统的线性方法尤其是在共振区段，低估了FPSO的运动响应。由于本文采用对称模型计算，避免了横摇运动模拟产生的误差。因此，在实际工程应用中有必要考虑非线性效应，不仅要关注共振点，还要对其附近区域予以足够的重视。