Chinese Journal of Ship Research
恶劣海况下船舶砰击颤振响应特性数值计算与试验研究
本文网址：http://www.ship-research.com/cn/article/doi/10.19693/j.issn.1673-3185.03249 期刊网址：www.ship-research.com
引用格式：张涵韬,王一雯, 孔祥韶, 等.恶劣海况下船舶砰击颤振响应特性数值计算与试验研究[J]. 中国舰船研究, 2024, 19(2): 148–158.
ZHANG H T, WANG Y W, KONG X S, et al. CFD-FEM numerical simulation analysis of ship whipping response under extreme sea conditions[J]. Chinese Journal of Ship Research, 2024, 19(2): 148–158 (in Chinese).张涵韬1,2，王一雯*1，孔祥韶1，郑成1，吴卫国1
1武汉理工大学绿色智能江海直达船舶与邮轮游艇研究中心，湖北武汉 430063 2武汉理工大学船海与能源动力工程学院，湖北武汉 430063
摘 要:［目的］针对恶劣海况下船舶所受砰击颤振响应现象，探究船舶非线性波浪载荷与瞬态高幅值砰击载荷的耦合作用。［方法］采用计算流体动力学与有限元方法（ CFD-FEM）相结合的双向流固耦合方法对S175集装箱船进行数值仿真计算，并与试验结果及切片理论计算结果进行对比验证；采用分段变截面弹性龙骨梁模型开展船舶的砰击颤振特性模型试验，基于 CFD-FEM双向流固耦合方法开展船艏砰击载荷及高频非线性砰击颤振响应特性分析，并与模型试验结果进行对比验证。［结果］结果显示，波浪砰击载荷对船艏颤振响应的影响不可忽视，6级海况下由砰击颤振诱发的二阶高频成分分量占低频波浪弯矩的59.86%。［结论］采用基于 CFD-FEM的双向流固耦合方法可准确计算船首砰击颤振响应；在高海况下船舶所受非线性波浪载荷及结构动态响应易受船首瞬态砰击载荷的影响，在船舶结构设计与安全评估中需考虑高频砰击颤振的情况。关键词：极端海况；砰击颤振响应；计算流体动力学；有限元方法；双向流固耦合方法中图分类号: U661.44文献标志码:A DOI：10.19693/j.issn.1673-3185.03249
CFD-FEM numerical simulation analysis of ship whipping response under extreme sea conditions
ZHANG Hantao1,2, WANG Yiwen*1, KONG Xiangshao1, ZHENG Cheng1, WU Weiguo1
1 Green & Smart River-Sea-Going Ship Cruise and Yacht Research Center, Wuhan University of Technology, Wuhan 430063, China
2 School of Naval Architecture, Ocean and Energy Power Engineering, Wuhan University of Technology, Wuhan 430063, China
Abstract: ［Objective］This paper investigates the coupling effect of ship nonlinear wave load and transient high-amplitude slamming load in view of the whipping response of ships under extreme sea conditions. ［Methods］A two-way computational fluid dynamics-finite element method (CFD-FEM) fluid-structure interaction method is used to carry out the numerical simulation of an S175 container ship, and the results are compared with the test results and 2D linear strip theory. The testing of a segmented ship model with a variable cross-section elastic backbone is carried out, and the analysis of the bow's slamming load and high-frequency nonlinear whipping characteristics is undertaken based on the two-way CFD-FEM fluid-structure interaction method. The results are then verified through comparison with the model test. ［ Results］ The influence of wave slamming load on the whipping response of the bow cannot be ignored. The second order highfrequency component induced by whipping accounts for 59.86% of the low-frequency wave bending moment in sea state 6. ［Conclusions］The proposed two-way CFD-FEM fluid-structure interaction method can calculate the whipping response accurately. The nonlinear wave-induced load and dynamic structural response of the ship under extreme sea conditions are affected by the transient slamming load of the bow. In conclusion, the high-frequency whipping response should be taken into consideration in the safety design and assessment stages of ship structures.
Key words: extreme sea conditions；whipping response；computational fluid dynamics (CFD)；finite element method (FEM)；two-way fluid-structure interaction method
收稿日期: 2023–01–19 修回日期: 2023–03–21 网络首发时间: 2023–05–25 16:27
基金项目: 基础加强计划重点基础研究资助项目（2020-JCJQ-ZD-225·11）
作者简介: 张涵韬，男，1998年生，硕士。研究方向：船舶结构安全可靠性。E-mail：zhanghantao@whut.edu.cn王一雯，女，1990 年生，博士，助理研究员。研究方向：船舶结构安全可靠性。E-mail：yiwenwang90@whut.edu.cn *通信作者:王一雯
0 引 言
在恶劣的海况下，船舶不仅会发生剧烈的船波相对运动，还会引起船体结构的高频动态响应。非线性波浪载荷会导致船体结构产生较大的动态弯矩，从而对结构安全产生较大影响。此外，高海况波浪往往会导致船体局部区域遭受波浪砰击载荷作用，从而诱发砰击高频颤振响应。波浪载荷试验是分析船舶水弹性效应的重要途径。按试验方式划分，主要包括实船试验和模型试验2 种。Lee 等[1]采用应力监测系统记录了8063TEU集装箱船在两年运营中的船体垂向弯矩，对载荷进行了长期预报并与计算值比较。然而，因实船试验所需成本较高、周期较长，且试验工况的选取也具有较大的局限性，所以绝大多数的波浪载荷试验都是以水池模型试验的方式开展。试验采用的模型主要分为2类：整体弹性模型和分[2]段弹性模型。林吉如等 于上世纪 90年代初针对弹性船模开展试验研究，其采用 TE702 材料制作了 2艘船模，测量航模在规则波与不规则波中的弯矩、甲板应变、加速度及运动。整体弹性模型试验虽然能较好地模拟实船的局部刚度及船体扭转振动模态特征，但对模型设计及制造工艺的要求较高。而分段试验模型则具有制造工艺简单、剖面载荷测量较为准确的优点，从而被广泛采用。Zhu 等 采用矩形龙骨梁和17分段船模研[3]究了 13 000 TEU 超大型集装箱船在规则波与不规则波中的垂向、横向和扭转振动响应，其中分段船模在垂向与扭转刚度方面实现了与原型的相似，但梁的测量因其位于甲板处并未满足中和轴[4]及扭心位置的相似条件。王一雯等 采用非线性水弹性时域计算方法及U型龙骨梁分段模型试验研究了宽扁型江海直达船的非线性波浪载荷高频响应特性，探究了宽扁型江海直达船在不同装载条件下的高、低频垂向载荷响应以及扭转高频响应特性，结果显示，在遭遇频率为5.2 rad/s 的入射波浪下，其6阶倍频与船体湿模态较为接近，诱[5]发了明显的高频振动现象。汪雪良等 采用2根非均匀截面梁的10 分段 VLCC船模进行试验，并与二维切片理论和三维水弹性理论计算结果进行对比验证，结果显示，在中、大波高下的波浪载荷的非线性现象会加剧，高频响应成分可达垂向总载荷的 20%～ 30%。Fonseca 等 通过水弹性模型[6]试验分析了高阶弯矩成分与低频波浪弯矩成分，认为 2阶弯矩成分的最大值约为低阶波浪弯矩的 50%。Drummen 等 计算比较了某集装箱船的[7]垂向弯矩的高阶响应传递函数与试验测量值，发现高阶垂向弯矩成分的计算值通常小于测量值。由于势流理论仅适用于小波陡、船体微幅运动与结构变形的情况，在恶劣的海况下，船舶高速航行时，船体运动响应及遭受的波浪载荷呈明显的非线性特征。基于势流理论的边界积分法不能有效模拟由流体分离、涡流、波浪破碎或甲板上浪效应产生的黏性效应，低估了流体阻尼效应且不能考虑由黏性引起的涡流耗散等，故需要求
解雷诺时均方程（ RANS）以解决上述问题。随着计算机的发展，采用计算流体动力学（ CFD）与有限元方法（FEM）相结合的流固耦合计算方法能较好地求解水弹性问题。根据计算得到的弹性结构的运动和变形是否反馈给CFD求解器，可将流固耦合过程分为单向耦合与双向耦合2类。有学者采用 CFD-FEM法对一艘驳船在波浪中的动态响应进行了数值计算，并在后续研究中以 S175集装箱船为对象对单、双向耦合方式进行对比，发现采用这2种方式均能充分捕捉流场的非线性因素，但对于高频非线性结果，双向耦合较单向耦合的计算结果更接近于试验值[8-9]。从 Jiao 等 开展的数值计算结果来看，采用 CFDFEM [10]双向流固耦合方法能较好地预测S175 集装箱船在波浪中的高频波激振动和颤振响应。综上，在恶劣海况下船舶易遭遇非线性波浪载荷与瞬态高幅值砰击载荷的耦合作用，高频砰击颤振响应也是不可忽略的因素，多种非线性波浪载荷的耦合作用将使船体结构中的总应力水平达到较高量级，影响船体结构强度。因此，在船体结构强度评估中需考虑船体的水弹性效应对非线性波浪载荷和砰击载荷等复杂载荷的影响。本文拟针对 S175 集装箱船，采用 CFD-FEM方法开展双向流固耦合数值仿真研究，并与文献 [11-13]中的数值模拟及模型试验结果进行对比分析，验证本文所采用方法的可靠性与有效性。然后，在此基础上开展分段变截面龙骨梁砰击颤振模型试验，并采用 CFD-FEM 双向流固耦合方法对考虑波浪下船首瞬态砰击载荷的颤振响应进行数值计算，探究非线性波浪载荷特性及船体结构动态响应的影响规律，为船舶结构的安全可靠性及结构轻量化设计提供支撑与建议。1 非线性水弹性数值计算方法及有效性验证1.1 有限体积法
STAR-CCM+采用有限体积法（ finite volume
况，则需采用隐式格式的双向流固耦合方法，CFD和 FEM模型在各时间步中多次交换数据信息，使得采用 CFD-FEM协同仿真方法可以模拟强非线性瞬态流固耦合问题，并能考虑结构动态响应及船首砰击等瞬态非线性载荷对流场的影响。船舶砰击颤振响应作为典型的强非线性瞬态流固耦合问题，对于其所遭遇的非线性波浪载荷与结构动态响应之间的强耦合作用，需采用隐式格式的双向流固耦合方法开展非线性水弹性时域分析。1.4 计算结果与比较对波长船长比 λ/Lpp = 1、航速 V = 22 kn、波高H = 4.8 m工况下的 S175 集装箱船模（图 2）进行了基于 CFD-FEM双向流固耦合的非线性水弹性时域数值仿真计算。箱船进行了系列波浪载荷模型试验与数值仿真研究。本文对文献 [11-14] 中使用的数值模拟和试验方法予以了总结，见表1。相比文献 [11-14]，本文考虑了流体黏性效应的影响，采用 CFD-FEM双向流固耦合方法对 S175集装箱船的运动特性和垂向波浪弯矩进行了计算，并与文献中的理论分析、模型试验以及数值计算结果进行了对比，如图3所示。从图 3可以看出，二维切片理论方法受物面非线性的影响，其垂荡和纵摇结果与 CFD-FEM方法计算结果间的差异分别达15% 和−25%。文献 [12] 开展的 30缩尺比船模试验其垂荡和纵摇运动数值较 CFD-FEM的结果分别偏小 13% 与12%。文献 [13]在试验中未测量船舶重心垂荡运动，其纵摇运动相比 CFD-FEM 方法偏大了 7%。文献 [12]认为由于势流理论忽略了流体的黏性效应，使得在诱发高频共振响应的遭遇频率附近，其垂荡运动响应被明显低估。文献[11] 中采用 CFD-FEM方法但未考虑黏性效应，与考虑了黏性的本文方法相比其垂荡响应略高，而垂向弯矩的波频成分略偏低，纵摇运动响应则较为接近。在该工况（ λ/Lpp = 1， V = 3.36 m/s， H = 0.12 m）下，船模稳态运动及船中非线性波浪载荷时域分析如图 4所示。图中，Z为垂荡幅值， θ为纵摇角为船中垂向弯矩。从图中可以看出，船舶在遭遇稳态规则波时，船体在规则波中呈现稳态运动且仅有低频成分，船中垂向弯矩出现明显的高频分量，船波相对运动诱发了船体结构的高阶波激振动与砰击颤振响应，在两者的耦合作用下，使得垂向弯矩的波峰与波谷呈现明显的不对称性。将船中垂向弯矩时间历程曲线通过傅里叶变换，得到频域结果如图5所示。由频域分析可知，
垂向弯矩响应包含有遭遇频率的低频成分，此外，还包括与波浪遭遇频率成2～8倍大小的高阶倍频成分。在该工况下， 8近船体的一阶固有频率，诱导其高频弯矩值急剧增大，导致船舶发生剧烈的非线性波激振动现象。虽然低频波频分量幅值占船中垂向弯矩的主要部分，但其各阶高频分量幅值的占比也不可忽略。采用带通滤波技术将船中垂向弯矩分解为波浪弯矩和高频弯矩两部分，如图 6所示。由图可见，低频波浪弯矩的时间历程响应呈现典型的简谐特性；高频弯矩部分因各阶高频波激振动分量和砰击颤振分量叠加而呈现明显的不对称性。由此可见，计及非线性因素影响的 CFD-FEM 双向流固耦合方法能较好地捕捉结构高频响应特性。2 变截面龙骨梁水弹性模型试验及数值计算为了探究高海况下船舶非线性波浪载荷与砰击颤振高频响应的流固耦合效应，开展变截面龙骨梁水弹性模型试验，针对非线性波浪载荷、船首瞬时砰击载荷以及砰击颤振响应展开分析，明确高、低频载荷的幅频特性以及砰击颤振的响应特性。试验拖曳水池长 131.32 m，宽 10.8 m，水深2 m；位于水池末端的20个液压造波板相互独立运动，可模拟不同周期与波高的规则波及不规则波。综合考虑拖曳水池尺度、造波能力、有效测试距离等因素的影响，设计加工缩尺比为1︰25 的分段船模。为能准确体现其高阶振动模态，本船模的分段数量为8段，其分段处分别位于 St3，St5， St7.5， St10， St12.5， St15 和 St17 站号处，各个分段之间通过铝制变截面龙骨梁连接。由于目标船的刚度分布沿船长方向变化较大，为使各剖面满足结构刚度相似的关系，采用了变截面形式的铝制龙骨梁，通过切削加工龙骨梁外径，使其满足模型与实船结构动力相似的条件，变截面梁尺寸如表 2所示。船模分段之间的外壳通过硅胶连接，可在保持船模外型线光顺的前提下，也能使各分段之间避免在船波相对运动的过程中相互接触。分段船模由玻璃纤维增强复合材料（ FRP）制作。所示。2.2 计算工况及参数实船及船模的主要参数如表3所示。模拟船
舶在复杂高海况下的作业，船速V = 15 kn，其波长船长比 λ/Lpp = 0.29，0.42，0.57，0.65，0.74，0.85，0.93， 1.0，1.13，1.35，1.6，1.8，波高 H = 2.5 m。2.3 流体计算模型根据 ITTC关于计算域的建议及要求[15]，船模
的坐标系原点位于船底中部、中横剖面和中纵剖面的交点处，其中，x轴与船长方向平行，y轴与船宽方向平行， z轴与波浪高度方向平行。在船舶中纵剖面采用对称平面边界以将计算域减半来减少网格数量，提高计算效率。计算域的范围为：
−2Lpp< x <1.5Lpp， 0< y <1.5Lpp， −2Lpp< z ITTC[15] 中提到，为了产生稳定不衰减的波，每个完整周期的波浪在纵向方向和垂向方向至少分别由 40 个和 16个网格单元构成。本文研究中，每个完整周期的波浪在纵向方向与垂向方向分别由 150 个单元和 20个单元构成，以精准捕捉自由液面处的船波相对运动。船体周围网格及自由表面细化网格如图9所示。2.4 结构计算模型 采用基于 FEM 的 Abaqus 软件开展结构动态响应分析，以及隐式动力学的方法求解，并使用 Hilber−Hughes−Taylor 隐式算子（ Newmark-β 的一种扩展） [16]进行时间积分。有限元模型与试验模型相似，同样采用龙骨梁模型，包含船体表面及变截面龙骨梁。龙骨梁用B31梁单元模拟，梁单 元纵向位于船舶重心处，并在梁单元上施加质量点，以满足各个分段的重量分布与试验模型相近。船体表面采用面网格单元 SFM4D4 进行模拟，梁单元与面网格单元通过运动耦合约束连接，三维壳单元将外部流体载荷传递到结构单元进行结构响应分析，结构状态也将反馈给表面壳体单元并传递至流体进行下一个时间步的水动力分析。全船结构有限元模型如图10所示。约束船体中纵剖面上的单元的横荡、横摇与艏摇运动，实现半船模型的对称约束，并考虑重力加速度的影响。本次计算考虑了阻尼的影响，并在B31梁单元材料属性中增加了瑞利（Rayleigh）阻尼系数，使得其更加接近于船模试验系统的阻尼。 2.5 求解器设置在 STAR-CCM+中采取隐式非稳态求解方法进行计算，为了满足计算要求且兼顾计算效率，时间步为 0.001 s。在各时间步内进行4换，并在每次数据交换的过程中进行3次内部迭代。在 Abaqus 软件中，隐式求解采用的时间步为 0.01 ms，选用较小的初始增量步有利于 Abaqus软件在耦合求解最初阶段位移的收敛[17]，也有利于整个耦合计算的收敛。Abaqus 软件在连续若干个分析步的收敛之后，增量步尺寸会自动增大至一定大小进行计算。3 结果分析与比较3.1 运动结果分析比较在规则波迎浪条件下，采用CFD-FEM 方法对船舶重心处的运动响应进行模拟，并与试验测量结果进行比较分析。该工况下的船舶垂荡与纵摇运动响应传递函数（ response amplitude operators， RAO）如图 11所示。图中， Z/ζ为垂荡传递函数， θ/cζ 为纵摇传递函数，其中c 为波数，ζ为波幅。在垂荡运动RAO对比中，船舶垂荡运动的无量纲化幅值是随遭遇频率的变大而降低的，试验结果与仿真结果吻合较好，但对于波长较大的工况，数值计算结果略低估了垂荡响应。在纵摇运动 RAO对比中，随着遭遇频率的变大，纵摇运动的无量纲化幅值变小并逐渐趋于平缓，采用CFDFEM图 11 模型试验和 CFD-FEM 方法的船舶运动响应RAO 对比Fig. 11 Comparison of ship motion RAOs results between model test and CFD-FEM method双向流固耦合方法计算的数值较为接近试3.2 弯矩结果分析比较对分段弹性龙骨梁模型试验中变截面龙骨梁St3， St5， St7.5， St10， St12.5， St15， St17 截面处的垂向弯矩结果进行傅里叶变换，其一阶遭遇频率对应的无量纲化垂向弯矩如图12所示。图中， g为重力加速度， B为船舶型宽。如图 12 所示，当遭遇频率 ωe = 4.987 rad/s 时，各截面的低频波浪弯矩达到极值。随着遭遇频率ωe 从 3.40 rad/s 降低至 3.18 rad/s，各剖面的垂向弯矩幅值显著降低，船中处的低频波浪弯矩成分降低至峰值的31.95%。此外，在不同遭遇波频下，沿船长方向的垂向弯矩在船中处达峰值，且越靠近船中位置处其所受到的低频波浪弯矩越大。针对不同波浪遭遇频率下的无因次化船中垂向弯矩响应，采用 CFD-FEM双向流固耦合方法所得结果与试验结果的对比如图13所示。由图可见，两者吻合较好，当ωe = 4.99 rad/s 时，其垂向弯矩达到峰值，当 ωe < 4.99 rad/s 时，仿真数值略微高于试验数值，特别是在ωe = 3.40 rad/s（波长船长比 λ/Lpp = 1.8）时，数值计算结果相比试验结果偏高 32.6%，这是因为当入射波浪的遭遇频率较小时，入射波波长达 11.36 m，较长的波浪对流体 边界造成了反射现象，导致计算区域内的稳态波浪个数较少，这对船舶垂荡运动及垂向弯矩计算精度均有影响。后续，可考虑通过增加流体域尺寸来提高计算资源，以满足计算量的需求，进而提高低频下非线性波浪载荷的计算精度。3.3 船体运动以及载荷的时域分析对波高 H =5m，波长船长比 λ/Lpp = 1.35，航速V = 15 kn工况下的波高、垂荡运动、纵摇运动以及船中垂向弯矩响应采用 CFD-FEM 双向流固耦合方法进行计算分析，其时域结果如图 14 所示。结果显示，各物理量水弹性时域仿真计算值与试验测量值的一致性较好。3.4 砰击载荷分析 在波高H =5m、波长船长比 λ/Lpp = 0.85 的波浪条件下，以 15 kn航速迎浪航行时船首出现了明显的入水砰击现象，船首受到瞬态高幅值的砰击载荷作用，波浪因受到挤压沿船首表面攀升出现了射流飞溅现象，如图15 所示。通过模型试验以及数值计算，可以得到船首瞬时砰击压力时间历程。当测点位于水下时，还会受到静水压力的作用，进行对比时，将这部分静水压力略去只保留水动压力。采用 CFD-FEM双向流固耦合方法计算船首压力，测点P1～P4 的位置如图 16 所示，其压力时程如图 17所示。由图可见，靠近船底处的测点P4受到的砰击作用严重，其峰值超过了 2.5 kPa。测点 P1～ P4的压力峰值较平钝，并随测点高度的变高而减小，这是因为船体表面向内倾斜，其斜升角为负值。采用试验与 CFD-FEM双向流固耦合方法得到的测点 P4的砰击压力如图 18所示。由图可 图 18 模型试验和 CFD-FEM 方法在测点P4处的砰击压力对比Fig. 18 Comparison of slamming pressure between model test and CFD-FEM method at point P4见，两者的砰击周期一致且砰击压力幅值相近。3.5 砰击颤振结构响应分析在波浪中航行时，船首底部在出水与入水的过程中会产生剧烈的瞬态砰击载荷，导致船体产生高频砰击颤振。图19所示为 St3 和 St7.5 截面处的垂向弯矩试验时程曲线。由图可见，在波峰处均存在高频非线性简谐特性，且相比船中附近的非线性特性，船首附近的高频成分更加明显。对其进行傅里叶变换，得到两者的频谱如图20 所示。图 20 St3 与 St7.5截面处垂向弯矩频域结果比较 Fig. 20 Comparison of frequency domain results of VBM at St3 and St7.5由频域分析结果可见，St3截面处的垂向弯矩因受船首砰击颤振的影响，除了遭遇频率波频ωe成分以外，还诱发了与遭遇频率成倍数关系的各阶倍频 nωe（n =2，3，…，7）成分，其分量数值随阶数的变大而减小。采用 CFD-FEM 方法可以准确捕捉 7阶倍频分量，其低频波浪分量与试验测量值几乎一致；在高频分量中，其数值略微小于