基于声固耦合方法舰船结构冲击响应数值预报研究

2020-05-21嵇春艳郭建廷

江苏科技大学学报(自然科学版) 2020年2期

嵇春艳,季斌,郭建廷,陈琅

(江苏科技大学船舶与海洋工程学院,镇江 212003)

大型水面舰船在历次海战中起着重要的作用[1-2]．大型水面舰船在战争中可能会遭受导弹、鱼雷、水雷等武器的攻击[3-4]．接触爆炸主要产生冲击波，造成舰船结构以及设备的破坏,而非接触爆炸一般不击穿船体结构,但会使得舰船上重要设备受到冲击破坏[5],以至于舰船失去生命力及战斗力．因此,舰船的抗冲击性能显得格外重要[6-7]．

影响舰船冲击响应的因素有很多,包含武器攻击的角度、爆炸的位置、冲击因子等．目前,大部分学者在单个因素下对舰船的冲击环境展开研究．文献[8]中仅围绕了气泡脉动载荷对冲击环境进行研究,并没有考虑冲击波载荷对全船冲击环境的影响．文献[9]中基于特征谱速度和冲击环境本征方程提出了一种舰船冲击环境的工程化预报方法,该方法考虑了外载荷参数、舰船尺度特征,但预报对象太过依赖设备,很难把握整个舰船的冲击环境．文献[10]中仅在相同的冲击因子下,研究了冲击环境沿全船的分布规律,并没考虑冲击因子的变化对冲击环境的影响．文中基于声固耦合方法对全船结构进行数值计算,将公开发表的实船试验数据与基于声固耦合方法的数值计算结果进行了对比分析,验证该方法的计算精度．同时,提出远场爆炸作用下船体冲击响应随纵向、垂向以及冲击因子变化时分布规律的预报公式,给出了水下非接触爆炸下舰船冲击响应的分布规律与外界参数之间的关系．

1 水下爆炸冲击环境理论

1.1 水下爆炸载荷概述

水下爆炸过程中,关于冲击波阶段和气泡膨胀与收缩阶段的压力模拟，在2002年已经得出了很好的估算公式[11]．该一系列公式建立在2个假设之上[7]:第一,水下爆炸过程中冲击波产生过后得到的均视流体为无粘无旋及不可压缩流体,不考虑外界压力变化对水密度的影响；第二,气泡脉动过程中将气泡的运动视为球形运动,忽略非球形变化对后续计算的影响．

1.2 基于声固耦合方法的结构动响应

大量文献表明声固耦合法计算水下爆炸载荷在水中的传播以及船体结构与周围流体相互作用具有较好的精度[12]．声固耦合法是对结构部分和流体部分按有限元法离散,并用声学单元来描述流场．冲击波在声学单元中传播,受到水下爆炸载荷作用的结构会发生材料和几何非线性响应．

流体有限元方程为:

(1)

通过虚功原理,可得与声场耦合的结构动力学方程为:

(2)

将式(1、2)合并的声固耦合分析的有限元方程为:

(3)

式中:p为流体单元任意一点的声压;u为流体与结构交界面上的结构位移;Ma、Ca、Ka分别为流场的总体质量阵、阻尼阵、刚度阵;Ms、Cs、Ks分别为结构的总体质量阵、阻尼阵、刚度阵;A为声场和结构的耦合阵,且AT=-As;Fs为结构的载荷阵．

根据Geers and Hunters模型以及水下爆炸载荷的基本特点,运用MATLAB软件自行编制三维气泡动力学程序得到冲击波和气泡脉动载荷[11]．在INP文件中定义冲击波和气泡脉动脉动载荷以声固耦合的方法加载到舰船结构上,进一步计算得到水下爆炸载荷作用下舰船的冲击环境．

2 声固耦合方法的验证

2.1 实船试验结果

选取文献[12]中6个典型位置,分别为船体首部、船体中部、船体尾部、3甲板首部、3甲板中部、3甲板尾部，实船爆炸试验谱速度值如表1.

表1 实船试验谱速度值Table 1 Spectral velocity value of real ship test m/s

2.2 数值计算结果

在ANSYS软件中建立某舰船有限元模型,同时在ABAQUS软件中以声固耦合方法进行水下非接触数值模拟．药包质量为162.5 kg TNT,爆距为35.67 m,攻角为30°,参数设置基本与实船试验工况一致．选取与实船试验相同的6个典型位置,具体谱速度值如表2.

表2 数值仿真谱速度值Table 2 Spectral velocity value of numerical simulation m/s

2.3 精度验证

将实船试验结果与数值仿真结果进行对比，如图1,图中的柱状图每个测点上方均表明了谱速度值误差,其中左侧柱状代表实船试验中测量值,右侧柱状代表数值计算值．从中不难发现,数值计算与实船试验之间误差均控制在35%以内,说明该数值方法具有良好的精度,用于工程计算是可信的．

图1 测点谱速度值对比Fig.1 Comparison of the spectral velocity values

3 典型舰船结构冲击响应计算

3.1 有限元模型

以某水面舰船为研究对象,运用有限元软件ANSYS进行1 ∶1比例三维实体建模,全船结构主要采用壳单元和梁单元,壳单元类型为S4R,梁单元B31,网格的平均尺寸为0.5 m．舰船主尺度如表3．

表3 舰船主尺度Table 3 Principal dimensions of ships m

在计算时,将增加MASS质量点单元来保证整船的质量分布及调整其重心位置．外部流场的单元类型为AC3D4,所建流场区域半径为船宽一半的6倍[8]．通过采用Tie的约束形式定义流场与船体结构面—面接触将流场与水域进行耦合;爆炸载荷则通过自编程序将其压力随时间变化的数据导入INP文件中,实现爆炸载荷在水域及船体结构中的传递．图2为舰船与流场耦合示意．

图2 全船有限元模型与流场耦合示意Fig.2 Sketch map of finite element model of ship coupling with fluid

3.2 材料参数的设置

常用材料本构方程主要包括Johnson-Cook(简称JC模型)和Cowper-Symonds(简称CS模型)两种模型[7]．而CS模型将材料动态屈服应力值与静态屈服应力值作为应变率的函数,参数较集中,方便研究者使用．所以文中采用Cowper-Symonds模型考虑应变率影响:

(4)

3.3 计算工况的设置

计算工况共8个,爆炸距离分别为1.16L、0.78L、0.58L、0.47L、0.39L、0.33L、0.29L、0.26L(L为船长),分别从船艏1/4L、船舯1/2L和船艉3/4L处设置药包,每个位置30°的攻角;8个工况中的药包质量均设定为1 000 kg，其中,攻角定义为爆点与参考点连线与水平面之间的夹角．图3为爆点设置示意,表4为工况设置．

图3 各爆点位置示意Fig.3 Sketch map of explosion location

表4 工况参数设置Table 4 Parameters of working conditions

数值计算结果采用无量纲的形式．无量纲船长和无量纲型深分别以水线间长及型深为标准量,无量纲谱速度以工况3中0.3L处对应的谱速度值为标准量进行无量纲化．

3.4 舰船结构冲击响应计算结果及分析

3.4.1 相同爆点下典型甲板的垂向冲击响应分布

水下爆炸冲击环境包括垂向和横向冲击环境,一般地垂向冲击环境较横向冲击环境危险,进而文中只研究水面舰船垂向冲击环境[12-13]．该舰船共4层甲板,除去01甲板外,1甲板、2甲板、3甲板在船长方向布置以及各甲板长度都较为接近,因此选择该舰船的1、2、3甲板作为研究对象,用以描述舰船的冲击环境在纵向的分布规律．

图4～6给出冲击因子为0.4时,L/4、L/2和3L/4爆炸位置下各甲板测点处无量纲速度沿纵向分布曲线．各甲板的测点均在中纵剖面上,其中横坐标各考核点到船艏距离的无量纲值．

从图4～6可知：各层甲板在爆点处的无量纲谱速度值达到最大,同时沿船长方向艏艉两端呈现整体的衰减趋势．由于爆点位置位于船体以下,相比之下,距离爆点较近的3甲板的谱速度整体上要大于其他两层甲板；2甲板的谱速度其次,距离爆点最远的3甲板谱速度最小．每张图上曲线的峰值基本上位于每个工况爆点处,符合船体冲击响应沿纵向的分布规律．

通过对以上曲线的分析,甲板谱速度值一般以爆点位置为分界点向船艏艉两端满足不同的数学衰减模型．如图7～9:对1、2、3甲板每个考核点的无量纲谱速度值进行拟合,拟合结果显示多项式数学模型拟合散点的效果较好,并且其他任何数学模型都可通过泰勒公式进行展开得到近似的多项式数学模型,具有一定的通用性,因此选用多项式模型对冲击谱速度进行预报是比较合适的．

各甲板在爆炸位置前后区域的数学模型均可用最小二乘法进行拟合．拟合曲线预报公式以三次多项式ax3+bx2+cx+d表示较为合适,具体参数如表5.

图4 船艏爆炸下各甲板谱速度分布曲线Fig.4 Distribution curves of non-dimensional spectral velocity in each deck for the L/4 condition

图5 船舯爆炸下各甲板谱速度分布曲线Fig.5 Distribution curves of non-dimensional spectral velocity in each deck for the L/2 condition

图6 船艉爆炸下各甲板谱速度分布曲线Fig.6 Distribution curves of non-dimensional spectral velocity in each deck for the 3L/4 condition

图7 1甲板纵向谱速度拟合分布曲线Fig.7 Longitudinal spectrum velocity fitting distribution curves in deck 1 under the third condition

图8 2甲板纵向谱速度拟合分布曲线Fig.8 Longitudinal spectrum velocity fitting distribution curves in deck 2 under the third condition

图9 3甲板纵向谱速度拟合分布曲线Fig.9 Longitudinal spectrum velocity fitting distribution curves in deck 3 under the third condition

表5 各甲板拟合曲线预报公式Table 5 Fast prediction formula of fitting curves for each deck

文中爆点位置位于船舯时,各典型甲板谱速度预报公式与文献[9]提出的通用型预报公式进行对比,为了验证表5所总结的多项式预报公式的精度．对比结果如表6,从中可知,表6所总结的多项式对谱速度的预报公式的预报结果与文献[9]中较为接近,误差控制在30%以内,具有一定的通用性．

表6 典型甲板预报结果对比Table 6 Fast prediction formula of fitting curves for each deck

3.4.2 不同冲击因子下垂向冲击响应分布

为了分析爆距的变化对舰船冲击环境规律的影响,以2甲板为研究对象,分别在船艏L/4、船舯L/2以及船艉3L/4处布置测点,测点均位于2甲板的中纵剖面．2甲板的0.25L、0.50L和0.75L处谱速度在不同爆炸位置随冲击因子变化曲线见图10～12,各自爆炸位置处的谱速度随冲击因子分布曲线如图13．

由图10～12可知,2甲板谱速度值随着冲击因子的增大而增加,而且在同一冲击因子下距爆点近的位置测点谱速度值要大于其他位置测点,这符合冲击环境的垂向冲击响应随爆炸距离减少而增大的基本规律．图13呈现的是0.25L、0.50L和0.75L位置处爆炸时谱速度随冲击因子的变化曲线,在不同冲击因子下,各位置谱速度值比较接近,但整体上船艏位置处最大,其次为船舯位置,船艉位置处谱速度值最小．可能由于在船艏位置处结构较轻、刚度变化大,而船艉位置处的支撑2甲板结构较多所致．

对不同爆炸位置处不同冲击因子的谱速度值进行拟合,如图14～16．图17为爆炸位置处在爆炸冲击时不同冲击因子的谱速度值拟合曲线．通过分析表明,随着冲击因子的增大,各爆点位置谱速度值基本呈线性的方式增长;在爆炸位置的曲线斜率较非爆炸位置的曲线斜率大,可见冲击波具有短时间内急剧增加,紧接着近似于指数规律进行衰减的特点,同时也说明了爆点位置一定时,典型纵向位置的变化将会影响冲击环境的谱速度值．在2甲板不同爆炸位置处谱速度随冲击因子变化的快速预报公式及其参数如表7．