基于计算流体力学的船舶风载荷计算方法分析

王 澄，张高峰，袁昌耀，吕晓勇，杨秀萍，王收军

（1.天津理工大学 机械工程学院天津市先进机电系统设计与智能控制重点实验室，天津 300384；2.南京晨光集团有限责任公司，江苏 南京 210006）

在全球经济一体化持续发展过程中，各国之间的贸易往来需求变得越来越大，这促使船舶运输往大型化、高效化发展.不管是从船舶能效还是航行安全性考虑，对影响船舶推进效率、稳定性、操纵性能的风载荷展开研究，都具有实际的工程意义.

目前，关于船舶风载荷的研究方法有：理论分析、风洞实验、计算流体力学（CFD）法.其中，CFD数值模拟相比于其他方法，具有发展成熟、费用少、周期短、通用性好等优点，也常被用于模拟计算结构的外流场情况.谢华平等[1]采用RNG k-ε湍流模型，研究不同工况下挡风系数、长宽比等因素对角钢格构塔的风载阻力系数影响；朱小海等[2]对港口起重机进行数值模拟，将流场各物理量计算结果与行业设计规范对比，验证数值模拟方法的可靠性；林一等[3]基于CFD方法对自升式平台所受风载荷进行数值模拟计算，结果与风洞试验相近，相比之下，规范计算结果则较为保守；倪歆韵等[4]采用SST k-ω湍流模型，评估两种特殊作业船所受环境载荷，模拟结果与试验结果比对吻合较好；贾宝柱等[5]采用Standard kε湍流模型，计算超大型油轮风载荷系数并与Isherwood公式估算结果比对，分析明显偏差存在的原因.

一般的干散货船都是单甲板船，所运输的货物品种单一，并且不需要被包装成捆，不怕挤压，装卸便宜.本文以干散货船为例，基于CFD方法对船舶周围风场进行数值模拟研究，根据干散货船型线图建立三维模型，建立内部、外部计算域并采用混合网格划分，计算平均风和梯度风两种风况下，船舶处于不同风舷角下的风载荷系数，并将数值计算结果与Isherwood公式计算结果进行比对，分析数值计算方法预估船舶风载荷的可行性.

1 风载荷系数与经验公式

1.1 船舶风载荷系数

船舶受风载荷的影响主要表现在纵向力、横向力和艏摇力矩，如图1，为便于在不同类型、尺寸的船舶间比较，常采用无量纲系数法来处理这三个载荷，这些无量纲系数的数值大小反映风场对船舶的影响程度，表达式为[6]：

式中：Cx为纵向风载荷系数；Cy为横向风载荷系数；Cn为艏摇力矩系数；Fx为风载荷纵向分力，Fy为风载荷横向分力，Mz为风载荷艏摇力矩；ρ为空气密度；V为计算风速；AT为水线以上船舶正投影面积；AL为水线以上船舶侧投影面积；LOA为船舶总长.

图1 风载荷示意图Fig.1 Schematic diagram of wind load

1.2 经验公式

目前求解风载荷系数的经验公式有：Isherwood方法、Hughes方法、Taylor方法等，从中选择应用广泛、精度较高和适用性好的Isherwood法计算风载荷系数，作为数值计算结果的参照[7]，表达式为

式中：A0～A6、B0～B6、C0～C5均为Isherwood常数；AT、AL、LOA同式（1）、（2）、（3）；B为船舶型宽；S为水线以上部分侧投影面积周长；C为侧投影面积中心与船首间距；M为侧投影中支柱或桅杆数目；Ass为上层建筑侧投影面积.

2 数值模型

2.1 模型建立

干散货船主要特征参数尺寸如表1.利用Solid Works绘制船舶水线以上部分结构特征，以1∶100比例建立三维模型，如图2所示.由于风载荷作用在船舶外表面，内部结构对其没有影响，所以不对船舶内部建模.同时，忽略船舶甲板上对风载荷计算影响微小的细节特征，比如镂空栏杆、缆绳钢索等，以便之后划分船舶外流场计算域.

表1 船舶特征参数尺寸Tab.1 Dimensionsof ship characteristic parameters

2.2 流场多计算域与网格划分

通过在ICEM中建立合理的外流场计算域来模拟船舶周围风场，数值模拟中的计算域存在边界，然而现实中的风场不存在边界，为了尽可能减小风场模型产生的误差，工程计算时数值风洞的阻塞比须低于3%[5]，阻塞比ε计算公式为.

图2 船舶水线以上特征模型Fig.2 Characteristic model of shipsabovewaterline

式中：A为船舶正投影于流场入口上的面积；AN为流场入口面积.

根据船舶模型尺寸，计算域的长取10倍船长，宽取30倍船宽，高取15倍船高，船体距离流场入口为3倍船长，距离流场出口为6倍船长，这保证了计算域中湍流可以完全发展.外计算域如图3所示，同时，为了便于计算不同风舷角时的风载荷，在船舶周围建立直径为2 m的内计算域，内计算域嵌套在外计算域中，内计算域如图4所示.在已建立好的流场计算域中，90°风舷角时，阻塞比最大，为0.7%，满足工程要求.

图3 外计算域示意图Fig.3 Schematic diagram of external calculation domain

图4 内计算域示意图Fig.4 Schematic diagram of the inner computing domain

计算域网格划分采用混合网格，外计算域使用六面体结构网格，如图5所示，由于船舶特征细节多为不规则表面，因此，内计算域采用四面体非结构网格，如图6所示，在船舶表面设置边界层网格，内外计算域之间设置interface面来实现不同网格类型间的数据交换，最终生成的计算网格数量为557万，网格整体质量优良.

图5 外计算域剖视图Fig.5 Sectional view of the outer computing domain

图6 内计算域剖视图Fig.6 Sectional view of theinner computing domain

2.3 垂直风速剖面

大气边界层内风的水平运动受地表摩擦力影响，风速随高度增大而增大，对于近地风场风速的表示有指数风速剖面和对数风速剖面，我国风工程结构抗风设计中常用指数风速剖面，故本文采用指数风速剖面，风速呈梯度递增，也称作梯度风，表达式为

式中：VH为距地面高度H处的平均风速；Vref为参考高度Href处的平均风速，参考高度Href取10 m，Vref取10 m/s[8]，考虑到模型按1∶100比例建立，Href代入计算时取0.1 m；α为地面粗糙度指数，取0.12.

2.4 湍流模型与边界条件

湍流模型选择基于RANS的Realizable k-ε模型，相对于其它湍流模型，该模型引入旋转和曲率相关的内容，收敛性好并且收敛时间短.文献[6]中选取了多种基于RANS的湍流模型，对船舶进行数值计算并比对结果，得出从模拟精度和计算时程两方面综合考虑，Realizable k-ε模型最合适.

流场入口设为速度进口边界，平均风速取10 m/s[9]，梯度风速为指数律形式，利用Fluent用户自定义函数UDF功能实现.雷诺数计算公式：

式中：ρ为空气密度，取1.225 kg/m3；L为特征长度，取2 m；μ为空气动力黏度，取1.8×10-5N·s/m2.湍流强度I=0.16（Re）-1/8=2.7%，湍流尺度l=0.007L=0.14 m[10].流场出口设为完全发展出流边界.左右壁面设为对称边界，这样可以减少计算域边界带来的计算误差.其它面采用默认的无滑移壁面.求解器设置为稳态计算的压力基求解，SIMPLEC算法，湍流动能和湍流耗散率方程采用二阶迎风格式.在监视器中设置，每一步迭代时计算船舶的纵向、横向阻力系数和力矩系数.

2.5 数值模拟结果

现实风场中，自然风作用于船身时的风向角是随机的不确定的，为了分析不同风舷角下流场特性，在0°～180°风舷角内，每隔10°对船舶进行一次数值计算，由于船舶整体具有对称性，不再模拟计算180°～360°风舷角内的船舶风载荷.从两种风场中分别选取0°、40°、90°、150°共4个风向角下的计算结果，对比分析不同风速下流场内的压强、速度分布情况.

图7—图10为风舷角不断增大的过程，纵向对比中可见船舶背风面的流场分布不再均匀，流线紊乱并且开始相互缠绕，尤其是在船体艉部上层建筑周围的流场中，漩涡现象更为明显；横向对比各个风舷角下的不同流场，可以发现相对于平均风风场，梯度风风场的压强范围和流速都略大（流速大小由流线颜色深浅判断）.

图7 0°风舷角流线和压强分布Fig.7 Streamline and pressure distribution at 0°

图8 40°风舷角流线和压强分布Fig.8 Streamline and pressure distribution at 40°

图9 90°风舷角流线和压强分布Fig.9 Streamlineand pressure distribution at 90°

图10 150°风舷角流线和压强分布Fig.10 Streamlineand pressure distribution at 150°

3 风载荷系数对比分析

图11—图13为两种风场数值模拟所得的风载荷系数与基于Isherwood公式计算结果的曲线对比，在部分风舷角下三条曲线出现轻微偏差，但整体上三条曲线走势一致，在量值大小上契合度较好；纵向风载荷系数曲线中，极值点出现在40°和140°风舷角位置，此时船舶受纵向风载荷影响最大；横向风载荷系数曲线中，极值点位于90°风舷角，此时两种风场结果偏差最为明显，由图9中船体艉部流线变化可看出，船舶背风面出现钝体绕流现象，两种风场中风速的不同，致使形成的气流漩涡尺度大小不同，能量消耗不同，使得两种风场的计算结果出现明显偏差；艏摇力矩系数曲线中，40°和150°附件三条曲线相差较大，主要是由于大角度绕流的复杂性，公式计算时无法考虑流场分布的复杂性.而风速不同时，产生的绕流运动发展程度不同，对船舶的影响也将不同.

图11 纵向风载荷系数Fig.11 Longitudinal wind load factor

图12 横向风载荷系数Fig.12 Lateral wind load factor

图13 艏摇力矩系数Fig.13 Coefficient of sway moment

4 结论

本文对船舶风载荷进行数值模拟计算，包括不同风舷角和不同风速下的模拟，将数值计算结果和Isherwood公式计算结果对比分析，所得结论如下.

1）计算船舶风载荷时，两种方法所得结果虽然在部分风舷角下存在一定偏差，但在整体走势和量值上吻合较好，表明数值模拟方法应用于船舶风载荷计算的可行性，通过数值模拟方法计算船舶风载荷，可以直观了解流场中风速、压强分布，进一步分析复杂风况对船舶的影响，提高预估船舶风载荷的精度，这对于船舶提高能效和安全性具有实际的工程意义.

2）在纵向和横向风载荷系数结果曲线中，相对于平均风风场，梯度风风场与Isherwood公式的契合度更高，采用垂直风速剖面数值模拟船舶周围流场，更符合实际风况.

3）40°、140°风舷角时，船舶受横向风载荷和艏摇力矩的影响最大，遭遇恶劣天气时，建议尽可能避免船舶风舷角处于这两个角度左右.