0 引　言

水面舰船在作战时会受到水中兵器的攻击，水中兵器在水下爆炸产生冲击波载荷和气泡载荷。冲击波载荷一般会导致舰船发生局部破坏，由于舰船已经有足够的强度，对于舰船来说局部的损伤已经不是致命的伤害。但是，气泡载荷作用下船体总纵强度将受到严重的影响，极端情况下会导致舰船在船中位置发生断裂，舰船沉没。此外，当结构与气泡一定距离时，气泡坍塌过程中会产生射流，对船体造成严重损伤。为了提高舰船生命力，有必要对水下爆炸气泡载荷对舰船结构损伤特性开展研究。

当前，水下爆炸载荷对结构影响以及气泡动态特性研究已经成为研究热点。张阿漫^[1]使用边界元方法对水下爆炸气泡的三维动态特性做了系统详细的研究，验证了边界元方法对于气泡动态特性研究的有效性，并且研究不同参数对水下爆炸气泡动态特性的影响，之后对近自由面、近壁面以及与弹塑性结构相互作用进行了研究。崔杰^[2]建立了一种更加接近实际情况以及模拟效果更好的气泡生成方法，对水下爆炸产生的气泡载荷以及其对结构的毁伤进行了试验研究。牟金磊等^[3]对爆炸气泡在不同边界条件下的动态特性进行了数值模拟和研究，并对存在不同边界时气泡运动受到的影响进行了仿真研究。李海涛等^[4]研究了箱型梁在近场爆炸作用下整体损伤特性，并结合试验结果和数值模拟结果对箱型梁的整体运动响应进行了分析。Ching等^[5]利用Abaqus软件的欧拉方法对水下爆炸产生的气泡运动的过程进行了模拟，并对一些特征进行了描述与讨论。李健等^[6]验证数值方法精度的前提下，对自由面附近的气泡与结构作用的整个过程进行了模拟，并对气泡脉动、射流等的产生影响的因素进行了研究。金辉等^[7]利用等效原理研究了射流载荷对船体的破坏效应，结果表明射流会使舰船受到严重的二次毁伤。张弩等^[8]建立气泡模型和弹塑性船体梁模型，对气泡载荷与舰船作用时，总体强度的响应进行了研究，对船体梁发生弹塑性损伤进行了机理和特性分析。Chen等^[9]对钢制船模进行了水下爆炸试验，监测冲击波和水下爆炸气泡载荷对船模作用时船模加速度、位移和应变等参数，并且对由于气泡运动产生的船模的动态响应进行了详细的研究，对船模的在气泡脉动作用下的模态进行了分析。Zong^[10]对漂浮于水面的自由-自由梁在水下爆炸气泡载荷作用下的动态响应进行研究，船体梁采用完全塑性材料。Zhang等^[11]对船体梁模型进行水下爆炸试验，使用高速摄像机拍摄船体梁在水下爆炸冲击波作用下的位移，通过动态图的形式展示了船体梁在冲击波作用下的中拱弯曲，并基于实验数据对船体梁的动态响应进行了分析，并对由于爆炸产生的局部变形情况进行了详细的分析。王恒等^[12]使用MSC.Dytran对近距离的非接触爆炸载荷下的加筋板的动态特性进行了研究，研究表明在计算模型合理的情况下，能够采用有限元软件对水下爆炸全过程进行模拟，文中还对多种载荷联合对加筋板的毁伤效应进行了研究。张弩等^[13]基于DAA方法，采用自行编制的结合BEM和FEM的计算程序，对船体结构在冲击波和气泡载荷作用下的动态响应进行研究。

针对舰船结构已经开展了很多试验、数值仿真研究，但是水下爆炸作用下船体梁结构响应研究还是不够清晰。采用实船进行水下爆炸冲击试验耗费巨大的人力、物力。因此，对于舰船结构损伤机理及特性的研究，一般采用简化船体梁模型。本文将针对简易船体梁模型在水下爆炸载荷作用下的动态特性开展研究。

1 数值方法 1.1 流体控制方程

对于确定的流体，其质量在运动过程中不生不灭，遵循连续性方程^[14]表示为：

$\frac{{\partial}\rho}{{\partial}t}+\nabla\cdot\left(\rho\boldsymbol{V}\right)=0。$

(1)

式中：$\rho$为流体密度，kg；$\boldsymbol{V}$为流体运动速度矢量，m/s²。

总动量随时间的变化率应等于作用在流体上面的所有体积力和表面力的和。对于可压缩粘性牛顿流体且其动力粘度为常数时，其动量方程^[14]为：

$\rho\frac{{\mathrm{d}}V}{\mathrm{d}t}=\rho f-\nabla p+\mu\nabla^2\boldsymbol{V}+\frac{\mu}{3}\nabla\left(\nabla\cdot\boldsymbol{V}\right)。$

(2)

式中：$p$为压强，Pa；$\mu$为动力粘度，N∙s/m²。

一般情况下，运动流体的能量方程可以表示为^[14]：

$\rho\frac{\mathrm{d}}{\mathrm{d}t}\left(e+\frac{v^2}{2}\right)=\rho q+\nabla\cdot\left(k\nabla T\right)+\rho f\cdot\boldsymbol{V}+\nabla\cdot\left(\mathbf{P}\cdot\boldsymbol{V}\right)。$

(3)

式中：$e$为热力学内能，J；v²/2为动能，m²/s⁴；$q$为单位时间内控制体内单位质量流体获得的热量，J/(s∙kg)；T为温度，K；k为流体的导热系数。

1.2 材料状态方程

炸药采用JWL（Jones-Wilkins-Lee）状态方程，由美国Lawrence Livermore实验室的E.L.Lee等人在Jones和Wilkins的工作基础上于1956年提出^[15]，JWL方程一般压力形式为：

$p = A\left( {1 - \frac{\omega }{{{R_1}\overline V }}} \right){e^{ - {R_1}\overline V }} + B\left( {1 - \frac{\omega }{{{R_2}\overline V }}} \right){e^{ - {R_2}\overline V }} + \frac{{\omega E}}{{\overline V }} 。$

(4)

式中：$\overline V = {V \mathord{\left/ {\vphantom {V {{V_0}}}} \right. } {{V_0}}} = {{{\rho _0}} \mathord{\left/ {\vphantom {{{\rho _0}} \rho }} \right. } \rho }$；$E = {\rho _0}e$为单位初始体积内能；$A$、$B$、${R_1}$、${R_2}$、$\omega$等为由试验确定的常数，取$A$=373.77 GPa、$B$=3.7471 GPa、${R_1}$=11、${R_2}$=0.9、$\omega$=0.35；${\rho _0}$=1630 kg/m³，$e$=6.0 GJ/m^{3 [16]}。

水的参数由线性Us-Up Mie-Grüneisen形式状态方程^[17]表示：

$p=\frac{\rho_0c_0^2\eta}{\left(1-s\eta\right)^2}\left(1-\frac{\mathit{\Gamma}_0\eta}{2}\right)+\mathit{\Gamma}_0\rho_0E_m 。$

(5)

式中：$\rho$=1 000 kg/m³；${c_0}$=1 500 m/s；$s$=1.92；$\mathit{\Gamma}_0$=0.1；$\mu$=1e-6 N∙s/m²。

空气采用理想气体状态方程表征，根据国际标准大气约定进行属性定义：密度为1.225 kg/m³、气体常数R₀=287、大气压力p₀=101 325 Pa、定容比热容c_v=0.716 kJ/(kg·K)。

1.3 数值方法有效性

数值方法有效性的验证采用与小当量水下爆炸试验^[18]的得到的气泡脉动、气泡形态以及流场压力等数据进行对比。如图1所示，模拟工况为55 g　TNT炸药在水下3.5 m处，同时炸药距底部3.5 m，在水平方向距离药包中心0.7 m处设置压力测点P，数值模拟的气泡在不同时刻的形态如图2所示。

由图2可以看到气泡膨胀到最大之后收缩，在最后产生射流，并且经过对比气泡中心位置在气泡收缩阶段气泡产生向上的迁移运动，尤其是在气泡体积达到最小值时。由此可以得出此方法对气泡运动过程中形态的变化能够进行很好的模拟。图3为气泡在运动过程中，气泡表面最上部、最下部与气泡中心初始位置的相对距离，并且将得到的数据进行无量纲化与文献[18]中的试验结果进行对比。数值方法计算得到的无量纲气泡半径与试验吻合较好，表明本文使用的方法能够有效的模拟水下爆炸产生的气泡形态。

图4为在测点P处数值模拟得到的气泡射流载荷曲线与试验得到的气泡射流载荷曲线对比，射流的压力增长和衰减的趋势和试验射流载荷衰减趋势吻合较好，其中压强和时间无量纲化分别以标准大气压和气泡脉动周期为基准。

基于以上的对比分析，可以认为本文所使用的数值方法能够较好地模拟水下爆炸气泡动态特性和产生的载荷特性，可以应用于水下爆炸载荷对结构等作用的研究中。

2 水下爆炸载荷作用下船体梁结构响应分析 2.1 工况设置

本文对水下爆炸载荷对船体梁结构作用模拟使用TNT炸药，药量12 g，爆距0.25 m，有限元模型尺寸6 m×6 m×6 m，其中空气域高度1 m，水域高度5 m，使用Abaqus/Explicit求解器进行求解。模拟使用的简单船体梁模型（见图5），尺寸为1.5 m×0.2 m×0.1 m，且均匀分为5个舱段，板厚0.002 m，使用Q235钢材，弹性模量为2.1×10¹¹ Pa，泊松比0.3，在模拟中钢材考虑应变率效应。其中，P₁～P₃、P₅～P₆位于各个舱段底板中心位置，P₄、P₇位于边中心位置，P₈～P₁₀位于侧边上，且P₈、P₁₀位于两端，P₉位于侧边中心位置。

2.2 仿真结果分析及讨论

图6为某些典型时刻计算结果，冲击波阶段采用流场中的压力云图展示，而气泡与结构的相互作用采用气泡与结构作用的形态图。

由图6(a)可以看到冲击波与结构作用发生反射和透射后继续传播。图6(b)~ 图6(f)气泡与船体梁模型相互作用形态图，图6(b)为气泡膨胀阶段，气泡在高压驱动下向外膨胀，同时船体梁模型受到气泡的作用产生向上的弯曲变形即“中拱”变形，这种变形是整体的变形并且在气泡作用下船体梁整体向上有一个刚体位移。随后，气泡继续膨胀（图6(c)），由于结构具有弹塑性，其在自身刚性作用下变形逐步恢复，变形不再明显，此时气泡内的压力继续减小。气泡达到最大半径以后（图6(d)），由于气泡内部压力小于外界压力，气泡开始坍塌，同时在流场中产生压力差，船体梁在压力差的作用下产生“中垂”变形，这是因为气泡坍塌时对结构的吸引，最初的中拱变形为气泡对结构的排斥造成。而后，气泡继续坍塌体积减小，同时结构也被气泡吸引产生更大的变形。最后，由于自由液面和结构的共同作用气泡被撕裂，上部产生向上的射流，下部分产生向下的射流。射流的方向由多种因素决定，在此不做深入讨论。可以看到，结构整体产生了很大的形变，当气泡能超过某个值以后，很可能船体梁被气泡吸引产生更大的形变，最后从船体梁的中部位置折断。

图7(a)为船体梁结构上的测点P₁、P₂、P₃位移时程曲线，图中的数据是相对于初始位置而言的，向上为正值，反之为负值。可以看到船体梁上的几个测点都有着相同的变化趋势，只是由于位置的不同变化幅度有一定的差异，在运动初期测点是向上移动的但在某一时刻测点开始向下移动，在另一临界时刻某些测点位置在初始位置之下。以上现象的原因就是由于冲击波载荷和气泡载荷作用不同效果造成的，由于P₁测点在船体梁底部船舯，当射流产生并且作用在底部时P₁有1个明显的向上的位移，当然其他测点也由于整体受到向上的力有1个整体的刚体位移，这与P₁测点的变形是不同原因导致的。图7(b)为沿长度方向某些时刻船体梁整体的位移变化，可以看到在气泡载荷阶段整个船体梁都发生了较大的中拱变形和中垂变形，这也说明气泡载荷对整个船体梁的总纵强度有非常大的影响。

图8为船体梁上测点P₁、P₁₀的应变时程曲线，可以看出在前期冲击波阶段P₁测点位置应变较大而P₁₀位置应变几乎没有变化，说明冲击波造成的损伤在接近爆源的位置较大，远离爆源位置的舷侧位置结构发生应变较小，主要原因是冲击波与底部是直接接触的，而侧面是传递过去的力。但在气泡载荷阶段由于船体梁发生整体变形，P₁₀测点的应变远远大于P₁测点位置的应变。

图9(a)为测点P₁位置的塑性应变时程曲线，可以看到应变先为正值，之后变为负值，这与之前讨论形态时观察到的现象吻合，由于变形方向的改变，测点位置所受到的力的性质也会发生改变，导致应变曲线发生改变。图9(b)为船体梁舷侧位置测点P₁₀的塑性应变曲线，可以看到P₁₀位置应变变化幅度远大于底部测点P₁，这说明在冲击波载荷、气泡载荷联合作用时，底部损伤不一定最大。对于实际船舶来说双层底的强度大于船舶舷侧结构的强度，因此舰船舷侧强度也需要引起重视。

图10为船体梁几个典型时刻发生变形时的Mise应力云图，从Mise应力云图可以清晰地看到船体梁模型不同部位产生的应力，不同时刻不同部位应力大小并不相同，但是显然船体梁中部以及中部侧面的应力很高，这和提取出的测点P₁、P₁₀等位置的应变以及塑像应变的现象相同，在其他位置的应力并不是很大。

图11为船体梁上测点的加速度对比图，船体梁在前期冲击波阶段产生较大的加速度，在气泡载荷作用阶段，很高的加速度没有延续，而是变得较小，由此可见联合作用时气泡对结构的加速度响应没有明显的影响，这是因为气泡载荷的变化周期远远超过冲击波，气泡载荷作用不如冲击波作用导致的加速度变化剧烈，但是这并不是说气泡对结构产生的毁伤小于冲击波。

3 结　语

船体梁在冲击波的作用下产生局部的变形，整体变形不大，值得注意的是整个船体梁的加速度比较大，这对于实际舰船来说，船上的设备、人员等在冲击波到达时需要承受较大的加速度，这需要得到特别的关注。

气泡载荷在冲击波载荷作用过后结构处于冲击波引起的动态响应时加载导致船体梁产生明显的中拱变形。随后，气泡收缩阶段，由于船体梁变形回弹以及水下爆炸气泡脉动引起流体流动对结构的吸引，使得结构整体发生位移，并且由于本文研究工况下气泡在船体梁中部正下方，导致船体梁中部发生的变形更加明显。