船舶在海上航行遭遇大幅海浪时易产生大幅运动，其可能导致结构因承受的载荷剧增而严重破坏，并对船舶的稳定性和安全性造成影响^[1]，尤其是大幅的垂向运动将导致甲板上浪现象的发生，而甲板上浪是一种对结构物极具破坏力的强非线性现象，朱仁传等^[2]研究表明，该现象主要是海上巨浪和船体的垂向运动(垂荡和纵摇运动)组合作用的结果。甲板上浪也受船型参数(如船首外飘等)、航速的有无及高低的影响，该现象导致的问题一直是船舶结构强度和航行性能的重要影响因素之一。开展甲板上浪问题的研究对船体性能预估、船舶的设计、制造和营运有着深远的意义。对于在微幅波浪中作微幅运动的船舶，可以将船舶在波浪中的受力和运动作为线性问题处理，因而可以应用叠加原理。但船舶作大幅运动如剧烈的垂荡、纵摇等运动时，现象复杂且强非线性特征显著，线性理论较难适用，因而具有较大的难度和挑战性，尤其对于数学解析和势流模型。因此，早期的研究大多以试验研究为主，后期有学者基于线性理论通过概率统计的方法对甲板上浪现象进行分析^[3]。这种传统的理论分析方法仅能描述船舶整体的运动状态，其概率分析无法计算波浪涌上甲板后对甲板及上层建筑等结构物的冲击载荷^[4]。后来的研究大多采用势流理论来对甲板上浪过程进行数值模拟，但是势流理论没有考虑到流体的粘性作用，其假设流体为不可压缩且流动无旋，难以处理非线性的影响，且其不能处理波浪破碎、翻卷、水花飞溅等强非线性现象。关于对甲板上浪问题更详细的研究进展可参见文献[5-6]等。舰船在实海域航行时，其周围的流场是非常复杂的，且航行中的非线性现象对船体及上层建筑的影响也不能忽视。段文洋等^[7]通过试验研究，分析了船舶艏部外飘对甲板上浪的影响以及翻卷状波浪和涌浪对甲板的砰击压力。梁修锋等^[8]通过对浮式生产储油卸油装置(floating production storage and offloading, FPSO)的模型试验和数值模拟，表明了采用并行运算解决甲板上浪等非线性问题的可能性；并采用CFD (computational fluid dynamics)软件FLUENT中的动网格以及VOF (volume of fluid)法对固定甲板情形下的单个极大波上浪问题进行速度剖面上变化规律的研究^[9]。基于求解非结构自由液面粘性流场的新型CFD软件FINE/Marine建立了有航速船舶的甲板上浪模型，并对该数值模型的可靠性进行了验证和研究，然后利用该模型重点对波浪中航行Wigley船甲板上浪的动态过程进行了数值模拟；最后，分析了海上波浪(波长和波幅)对甲板上浪的影响。在建立的固-液-气三相模型中，虽然气相参与模拟计算，但因气体影响与液体相比甚小，故在文中不单独进行分析。

1 控制方程及自由液面捕捉法

基于粘性流理论，采用SST κ - ω(Menter)湍流模型^[10-13]。求解器基于非结构化动网格技术，通过计算不可压缩的非定常雷诺平均方程模拟流场变化，基于有限体积法对输运方程进行空间离散。速度场和压力场分别由以下动量守恒和质量守恒方程^[14]求解得到

$ \frac{\partial }{{\partial t}}\int_V {\rho {\rm{d}}V} + \int_S {\rho \left( {\mathit{\boldsymbol{U}} - {\mathit{\boldsymbol{U}}_d}} \right) \cdot \mathit{\boldsymbol{n}}{\rm{d}}S} = 0 $

(1)

$ \begin{array}{*{20}{c}} {\frac{\partial }{{\partial t}}\int_V {\rho {U_i}{\rm{d}}V} + \int_S {\rho {U_i}\left( {\mathit{\boldsymbol{U}} - {\mathit{\boldsymbol{U}}_d}} \right) \cdot \mathit{\boldsymbol{n}}{\rm{d}}S} = }\\ {\int_S {\left( {{\mathit{\boldsymbol{\tau }}_{ij}}{\mathit{\boldsymbol{I}}_j} - p{\mathit{\boldsymbol{I}}_i}} \right) \cdot \mathit{\boldsymbol{n}}{\rm{d}}S} + \int_V {\rho {\mathit{\boldsymbol{g}}_i}{\rm{d}}V} } \end{array} $

(2)

式中：ρ为密度，V为控制体，S为控制面，U_d为S上n方向的速度，U和p分别为速度和压力，U_i为在x_i方向上的平均速度分量，τ_ij和g_i分别为粘性应力张量和重力矢量，I_i和I_j分别为方向向量。

采用BRICS (blended reconstruction interface capturing scheme)可压缩型离散格式的自由液面捕捉法来处理多相流的交界面，通过引入体积分数c_i表征网格单元中第i相流体的占比，0 < c_i < 1代表流体发生混合，因此将c_i=0.5定义为多相流的交界面，以此来重构固-液-气多相流模型中的自由液面。通过以下方程联立求解c_i：

$ \frac{\partial }{{\partial t}}\int_V {{c_i}{\rm{d}}V} + \int_S {{c_i}\left( {\mathit{\boldsymbol{U}} - {\mathit{\boldsymbol{U}}_d}} \right) \cdot \mathit{\boldsymbol{n}}{\rm{d}}S} = 0 $

(3)

$ \rho = \sum\limits_i {{c_i}{\rho _i}} $

(4)

$ \mu = \sum\limits_i {{c_i}{\mu _i}} $

(5)

$ 1 = \sum\limits_i {{c_i}} $

(6)

2 数值模型及其验证 2.1 船体模型

模型采用船长2.5 m的修正Wigley船型，船型公式为

$ \begin{array}{*{20}{c}} {\left[ {1 - {{\left( {\frac{z}{d}} \right)}^2}} \right]\left[ {1 - {{\left( {\frac{{2x}}{L}} \right)}^2}} \right]\left[ {1 + {a_2}{{\left( {\frac{{2x}}{L}} \right)}^2} + {a_4}{{\left( {\frac{{2x}}{L}} \right)}^4}} \right] = }\\ {\frac{{2y}}{B} - \alpha {{\left( {\frac{z}{d}} \right)}^2}\left[ {1 - {{\left( {\frac{z}{d}} \right)}^8}} \right]{{\left[ {1 - {{\left( {\frac{{2x}}{L}} \right)}^2}} \right]}^4}} \end{array} $

(7)

式中：a₂=0.6, a₄=1.0, α=1.0，L、B、d分别表示船长、船宽和吃水。修正Wigley船型主尺度在表 1中给出。

2.2 网格划分

由于船型和计算问题均具有左右对称性，因此数值模拟采用半船模型，计算域入射面采用造波边界，对称面采用对称边界，流出面和侧面采用远场边界，计算域底面和顶面采用规定压力边界，即在初始化计算过程中施加压力。经过对计算域大小、网格数量和时间步长等进行系统的收敛性验证^[15]，计算域尺寸最终设置为12.5 m×3.75 m×5 m，船艏距入射面2.5 m，先后采用无上层建筑和上层建筑高0.05 m(船型相同)的两个船体进行模拟。有上层建筑的船体计算域网格如图 1所示。

2.3 可靠性验证

为了验证该数值模拟的可靠性，首先设置了与物理模型试验相同的参数，对Wigley船的辐射问题进行了模拟^[16]。在傅汝德数为Fr=U/(gL)^1/2=0.2的迎浪航行条件下进行模拟，强制垂荡幅值为0.01 m，无因次化波数KL为4.0 ~ 30.0。先对该船强制垂荡运动的附加质量和辐射阻尼系数进行模拟求解并与试验数据进行对比验证，然后研究和分析了波长和波幅对甲板上浪的影响。

对比验证结果如图 2所示，图中实线是通过EUT (enhanced unified theory)方法计算得到的结果，A_ij和B_ij分别为附加质量和辐射阻尼系数。从图 2中可以看出整体模拟结果与试验数据吻合较好，因此可利用该模型进行接下来的大幅运动与甲板上浪的时域模拟。

3 数值模拟与结果分析 3.1 波长选取

对于迎浪航行的船舶，垂荡和纵摇等垂向运动对甲板上浪具有主要影响，为了确定Wigley船模型更易发生甲板上浪现象的波长，首先进行了在固定波幅的不同波长条件下船体运动响应的模拟。从若干组模拟中选取了具有代表性的四组波长数据，其垂荡和纵摇的时历曲线如图 3所示，其中t为计算时间；T_e为遭遇周期，艉倾角为正。

图 3(a)所示相对波长(λ/L)为0.8的海况下，船体的垂荡和纵摇幅值与其他三组波长相比甚小，故此处不进行赘述。随着波长的增加，运动幅值也逐渐增大。图 3(b)相对波长1.2的海况下，船体产生了大幅度的垂荡和纵摇运动，当船体中心位于波谷时，船体处于艏倾状态，由于该海况下波长3 m，船艏处于约90%波幅处。当波长达到图 3(c)所示的4 m海况时，其运动幅值略大于3 m海况，当船体中心位于波谷时，船体也处于艏倾状态，船艏处于约40%波幅处，综合其幅值变化考虑，3 m波长更易发生甲板上浪。之后随着波长的增加，运动幅值出现减小的趋势，如图 3(d)所示的5 m海况，其幅值与3 m海况大致相同，当船体中心位于波谷时，船体已基本处于水平状态，且船艏处于零波幅附近，因此较难产生甲板上浪现象。根据上述分析，选取最具代表性的相对波长为1.2和2.0的两组波浪进行系统的船体运动与上浪模拟。

3.2 波浪设置

根据实际海域中的波浪等级划分，选取不同海况波幅的波浪进行船体大幅运动与甲板上浪现象模拟；海况选取按船长为200 m的实船进行比尺缩放，波浪等级表与模拟参数设置如表 2和表 3所示。因波幅0.03 m的中浪海况已发生非常明显的甲板上浪现象，所以重点展示和分析了轻浪、中浪以及大浪海况的模拟结果。

3.3 波浪对船体运动的影响

甲板上浪是一种对结构物极具破坏力的强非线性现象，主要是海上大幅波浪和船体大幅运动组合作用的结果；甲板上浪也受船型参数(如船艏外飘等)、航速有无及高低的影响。作为研究的第一阶段，限于篇幅，重点研究在船型和航速保持不变的情况下，海上波浪(波长、波幅等因素)和船体大幅运动，尤其是垂向运动对甲板上浪的影响。

对于迎浪航行的船舶，垂荡和纵摇等垂向运动对甲板上浪具有主要影响。为了研究在不同波长下，波幅对船体大幅运动的影响，及垂向运动对甲板上浪的作用，对上述两组海况进行了模拟。图 4、5给出了λ/L=1.2时，不同波幅下的船体运动响应与受力情况；图 6、7为λ/L=2.0时的船体运动响应与受力情况。

图 4~7中T_z表示船体垂荡运动的相对位移；R_y表示船体纵倾角度，艉倾角为正；C_z为无因次化垂荡力系数；C_y为无因次化纵摇力矩系数，

$ {T_z} = \frac{z}{d} $

(8)

$ {C_z} = \frac{{2\left( {{F_z} - {F_{z0}}} \right)}}{{\rho {U^2}S}} $

(9)

$ {C_y} = \frac{{2{M_y}}}{{\rho {U^2}SL}} $

(10)

式中:z为船舶垂荡运动的位移；F_z为船体沿z方向所受合力，F_z0为船舶在静水中航行时的受力；ρ为流体密度，其值为998.4 kg/m³；U为航速；d、L、S见表 1。

从图 4和6可看出，不同海况下船体运动响应差别明显。图 4中垂荡与纵摇幅值的相位差约为0.14 T_e，明显小于1/4周期，而图 6所示相位差约为0.20 T_e，接近1/4周期。即当船体上升达到极值时，λ/L=1.2波况下的船舶处于艉倾状态，而λ/L=2.0的处于水平状态。垂荡运动与纵摇运动耦合的相位差，使得船体运动响应不同，从而影响甲板上浪现象的发生。

3.4 甲板上浪过程及原理分析

朱仁传等^{[2, 7]}研究表明：甲板上浪是由海上大幅波浪和船体的垂向运动(垂荡和纵摇运动)组合作用的结果。为了分析甲板上浪与船体运动响应和波浪载荷的关系，图 8给出了λ/L=2.0, H/λ=0.020中浪海况下的一个周期内典型节点的船体运动姿态，并结合船体运动响应曲线(图 4)和波浪载荷时历曲线(图 5)进行了分析。

图 8(a)为周期开始时刻，此时船艏位于约15%波峰处，波浪涌上甲板并大量分布在甲板上。从图 4和图 5中可以看出此时船舶整体上升，流体作用力大于静平衡时受力，因此处于加速上升状态；同时船体绕Y轴的力矩为负值，纵摇角为正值，船体呈现艉倾姿态且艉倾角逐渐增大，但艉倾速度逐渐减小。在随后的1/4周期内，船体继续上升，艉倾角增大到幅值后开始逐渐减小，达到图 8(b)所示状态，船艏处在波谷附近，此时涌上的波浪少量残留于甲板边缘。船舶接近达到最大上升高度，垂荡方向受力处于谷值附近；同时纵摇角仍为正值，受到负力矩影响，船体呈现艉倾姿态且艉倾角逐渐减小，纵摇速度逐渐增大。然后，船体在达到上升幅值后开始下沉，艉倾角持续减小，最终减小到负值，即图 8(c)的状态，船艏处在约15%波谷处，甲板上已无水体积存。船舶整体有下沉趋势，流体作用力小于静平衡时受力，因此处于加速下沉状态；负倾角受正向力矩的影响，呈现艏倾姿态，艏倾角不断增加，但艏倾速度逐渐减小。随后船体一直处于下沉过程，艏倾角增大到极值后开始逐渐减小，如图 8(d)所示，船艏处于波峰附近，波浪开始涌上甲板。船体的下沉接近极值，垂荡方向受力在峰值附近；在艏倾姿态下，持续受到正向力矩的影响，艏倾角不断减小且纵摇速度增加。然后船体由下沉极值位置逐渐开始上升，同时由艏倾逐渐过渡到艉倾姿态，继续开始下一个上浪周期。通过以上船体各运动状态的分析可得到：图 8的船体运动和甲板上浪过程，与图 4、5的数据能够一一对应；图 4、5、8能够较完整地描述甲板上浪的整个过程。

3.5 波幅对甲板上浪的影响

通过3.4节的分析，可以较为清楚地了解甲板上浪发生的过程及原理，本节将根据图 4、5以及图 8和图 9研究和分析海上波浪(波幅)对甲板上浪的影响。

为了研究波幅对甲板上浪的影响，保持入射波长不变(λ/L=1.2)，改变波幅模拟上浪过程。图 9给出了波陡为H/λ=0.007轻浪海况时的船体运动姿态(图 8为H/λ=0.020的中浪海况)。图 9没有出现甲板上浪现象，而图 8所示海况中出现了严重的甲板上浪现象。对比图 4、图 5所示的轻浪和中浪海况可知：两组海况的运动响应与水动力特性有非常相似的规律。但中浪下的垂向运动幅值均大于轻浪海况，这是图 8出现严重甲板上浪，而图 9没有甲板上浪的主要原因。根据图 4、图 5和图 8、图 9的对应展示可以得到：对于相同波长下的同一船型，在波幅达到一定高度时才会发生甲板上浪，且波幅越大，甲板上浪现象越严重，船体的运动响应幅值越大，受到的流体力和力矩也相应增大。此外，轻浪和中浪海况下的船体运动过程存在一定的相位差。

3.6 特征波长对甲板上浪的影响

在设置波长时已对多个波长进行了研究和分析，本节将重点展示两组波长海况的模拟结果，通过船体姿态和运动响应更详细地分析海上波浪(波长)对甲板上浪的影响。

图 10中的波浪幅值与图 8的相同，均为0.03 m。图 8(波陡H/λ=0.020)出现了严重的甲板上浪现象，而图 10(波陡H/λ=0.012)只是出现了轻微的甲板上浪现象，表明在研究的两组海况下，大波陡更易发生波浪越顶和甲板上浪。Stansberg等^[17]对无航速FPSO的甲板上浪进行试验研究时，也给出了大波陡是产生波浪越顶和甲板上浪先决条件的结论。同时也印证了前文得到的结论：相对波长为1.2附近的波浪更易发生甲板上浪现象。

图 11展示的是两组海况下的船体运动响应。由于波长不同，波浪传播至船体的时间相差明显，因此图中的运动响应幅值会存在较大的相位差。从中可以看出，船舶在相对波长为2.0的波浪中航行时，运动响应幅值更大，而不是甲板上浪更明显的相对波长为1.2的波浪。纵摇角的艏倾幅值要略微大于艉倾幅值，其可能原因是甲板上浪现象导致。船舶艏倾后，波浪涌上甲板，产生使船体继续艏倾的作用力矩，由于短波海况的甲板上浪更严重，其对应的最大艏、艉倾角差值明显大于模拟的长波海况，该现象也印证了上述分析。此外，从纵摇运动来看，一个周期内长、短波的艏倾运动周期均略长于艉倾运动周期。

图 12为两组海况下船体的受力情况。从中可以看出，短波海况的垂荡方向受力幅值大于长波海况，且受力更不规则，尤其是船舶在上升过程中，其受到波浪影响大于长波海况；短波海况的纵摇力矩幅值也大于长波海况。此外，船体艏倾时受到波浪的砰击作用，会使得船舶达到艏倾力矩幅值的时间滞后，由于短波海况的甲板上浪现象更明显，其滞后效果也更明显，其载荷变化规律呈现出不规则的正弦曲线。

3.7 甲板及上层建筑的冲击载荷

甲板上浪导致的问题一直是船舶结构强度和航行性能的重要影响因素之一。限于篇幅，本节重点研究甲板上浪现象对甲板及上层建筑造成的冲击载荷。选取甲板及上层建筑表面沿船中线均匀分布的点作为记录点来记录各时段的动压力变化。记录点1~6为甲板中线上的点，间隔0.05 m；记录点6~16为上层建筑前表面上的点，间隔0.005 m。根据选择的海况，本次模拟并未记录上层建筑上表面的载荷情况；因上层建筑前表面点数量较多，文中重点展示具有代表性的记录点6~13的数据。通过这些数据进行甲板上浪过程的分析，记录点在船上的具体位置如图 13所示。

图 14为在中浪及大浪海况下，甲板上各记录点的载荷情况。在中浪海况下，甲板上浪现象造成甲板上的记录点的载荷峰值由船艏至船艉依次减小，且在上层建筑底部会造成两次冲击。冲击载荷峰值与上浪的剧烈程度有关，在上浪过程平缓的中浪海况下，第一次冲击的峰值与记录点2处的载荷峰值基本一致；在上浪过程剧烈的大浪海况下，第一次冲击的峰值明显小于刚上浪时波浪对甲板的冲击载荷峰值。

从图 15(a)可看到船体发生了剧烈的甲板上浪现象，波浪大面积涌上甲板，因此在图 14中可以看到此时处于甲板中前段的记录点2、3处的载荷均处于峰值附近，峰值约为H/λ =0.020中浪海况的四倍，而靠近上层建筑的记录点4~6处的载荷变化较小。随后船体由艏倾姿态运动到艉倾姿态，1/4周期后达到图 15(b)所示时刻，上浪过程经过发展，涌上甲板的水体被上层建筑阻碍，在记录点6处产生第一次冲击，但峰值明显小于船艏记录点2，记录点4~6处的载荷处于峰值附近，而记录点2、3处的载荷明显下降。图中所示的上浪情况剧烈，水体大面积覆盖在甲板上，对应图中记录点的峰值区间延长，且对上层建筑产生第一次冲击的时间提前。甲板上的水体继续发展到上层建筑前表面即图 15(c)所示状态，记录点2~5处已基本没有载荷。最后作用在上层建筑前表面的水体部分落回甲板，形成第二次冲击，但其峰值远小于其他记录点峰值，约为第一次冲击峰值的10%~30%。

对比图 14中两种海况的甲板载荷情况可以看出，大浪海况下的甲板上浪更加剧烈，波浪冲击到上层建筑的时间更短，冲击载荷成数倍增加；两种海况均会在上层建筑底部形成二次冲击，大浪海况下水体回落至甲板的时间滞后。

图 16重点展示了上层建筑前表面的部分记录点处载荷，其载荷变化情况与图 15所示运动姿态吻合。图示两种海况下，底部记录点6处的载荷峰值最大，向上记录点载荷峰值依次减小；中浪海况下的峰值减小梯度远大于大浪海况，且各记录点峰值存在相位差，大浪海况的载荷峰值无明显相位差；大浪海况下波浪作用在上层建筑的时间远长于中浪海况。

4 结论

1) 对相同波幅同一船型，相对波长为1.2，遭遇周期在0.95 s附近的波浪更易出现甲板上浪现象，且与相对波长为2.0的波浪相比，船体受力更大，而船体的运动响应幅值较小；

2) 相同波长或波幅条件下，波陡越大越易产生甲板上浪现象且现象更明显；

3) 一个完整周期内，艏倾角幅值大于艉倾角且艏倾运动周期略长于艉倾运动周期，此外甲板上浪现象还会导致船体受力出现不规则的正弦规律。

4) 甲板与上层建筑的交界处附近会产生两次冲击，第一次冲击的载荷峰值与上浪情况有关，上浪过程平缓时，其峰值接近甲板前端的最大载荷峰值，随着浪级增加，该处的峰值逐渐小于最大载荷峰值，且二次冲击载荷的峰值远小于第一次冲击。

文中研究重点在于探索甲板上浪的发生机理以及其对甲板和上层建筑的冲击载荷。船舶大幅运动和甲板上浪的模拟，可为甲板上浪的预报和船型的优化设计提供一定的技术支持。下一步将以真实船型为研究对象，对冲击载荷响应进行更为系统的研究。