基于LS-DYNA的带球封头耐压结构深海碰撞仿真

郭桐桐，张爱锋，俞白兮

（中国船舶科学研究中心，江苏无锡214082）

0 引 言

海洋不仅是人类交通运输的重要通道也是自然资源的宝库，研究表明海洋蕴藏着大量的石油、天然气、煤、铁等矿产资源。具有一定自航能力的深海潜器在海洋资源开发和能源勘探中有着广泛的应用。类似于深海空间站这类大型深海潜器具有下潜深度大、海底情况复杂、工作环境恶劣的作业特点，在近海海底自航时存在着与海山、礁石等海底固定物发生碰撞的可能。而这类深海潜器正常航行时自身处于高应力状态，在静水压力和碰撞载荷的联合作用下极易发生破坏，这种载荷环境是最危险的工况，也是潜器发生水下事故的主要原因[1]。

据统计[2]，超过40%的潜艇水下事故是由于碰撞产生的，主要形式包括潜艇之间的相互碰撞、触礁和搁浅。2005 年美国的“旧金山”号核潜艇在关岛海域全速航行时与海底山脉发生碰撞，造成了艇体严重受损近百人死伤的严重后果；2009 年英国“前卫”号核潜艇和法国“凯旋”号在大西洋发生碰撞事故。因此，对潜器深海碰撞过程进行动力响应分析对提高其自航安全性有着重要意义[3]。

自上世纪五十年代以来，国内外学者对水面船舶的碰撞问题进行了系统的研究。Minorsky[4]等人采用经验法、简化解析法和试验方法在船舶碰撞方面做出了开创性的工作。对于潜器的水下碰撞问题，目前公开的资料较少，相关的研究也不多。通过理论方法研究难度较大，而试验方法成本又太高。因此，有限元方法是分析耐压壳体水下碰撞问题目前最经济、方便的方法。于昊、李陈锋[5]等人通过MSC.Dytran 软件对潜器舷侧碰撞进行了分析，得到了碰撞载荷和结构变形的情况。钛合金材料具有比强度高、耐高温、耐腐蚀、耐冲击、抗疲劳性能好等特点，常用于深海大型结构件。但较易发生屈曲破坏，对于其动态力学性能的了解也不够充分。陈敏、陈伟[6]等人通过大量试验对TC4 钛合金的动态力学性能进行测试，拟合出本构模型和失效模型。

本文基于显式动力学程序LS-DYNA 对某深海潜器钛合金耐压壳体与海底固定物发生碰撞进行数值仿真，从结构损伤变形、碰撞力和能量转换等方面对碰撞动力响应进行研究，分析了碰撞失效模式和剩余承载能力，重点讨论了撞击速度、静水压力、撞头形状等因素对碰撞结果的影响。

1 碰撞计算模型和参数选取

1.1 碰撞积分方法（显式积分法）

碰撞问题是一种非线性瞬态响应过程，位于碰撞区域的结构起初处于弹性状态，随后迅速进入塑性流动状态，并可能发生屈服或者屈曲等破坏。显式非线性有限元法是解决此类碰撞问题的主要方法。非线性有限元程序采用的计算方法为显式积分方法，LS-DYNA 解决碰撞问题即采用显式积分法。

1.2 接触算法的选用

ANSYS-DYNA 中提供了多种可供显式分析使用的接触算法，分为单面接触、点面接触和面面接触。本模型为耐压壳体与岩土结构的碰撞，接触区域较大且形状不规则，因此选用面面接触。面面接触常用于一个物体的面穿过另一个物体的面，是最为常用的接触算法。

由于耐压壳体采用壳单元建模，选用面面接触中的自动接触算法（ASTS）。自动接触算法可以考虑单元厚度的影响，允许接触出现在壳单元的两侧，在计算接触力时更加准确。

1.3 材料的动态力学性能

本文模拟某大型深海潜器的水下碰撞，耐压壳体采用钛合金材料（TC4）。海底的岩土种类较多，本文选取了一种较为典型的花岗岩材料作为海山的主要成分进行数值模拟。

1.3.1 TC4的动态力学性能

由于碰撞过程具有大变形和高应变率等特点，需要对材料的动态力学性能进行准确的描述。各国学者提出了许多本构模型。陈敏、陈伟等人通过对TC4 本构模型的研究发现：独立考虑应变、应变率、温度影响的J-C 本构模型的计算结果和试验数据较为吻合，故采用这一本构模型来模拟TC4的动态力学性能。相应的失效模型采用J-C 失效模型，该失效模型不仅考虑了应力三轴度、温度和应变率对材料失效的综合影响，而且还考虑了变形路径的影响[7]。

本文所用材料的J-C本构和失效模型的参数如表1所示。

表1 TC4材料模型参数[8]Tab.1 TC4 material model parameters

1.3.2 Gruneisen状态方程

水下碰撞是大变形、高应力的碰撞过程。钛合金材料的弹性模量和剪切模量在发生大变形时不再保持常数，而是随着压力变化。固体材料在高压状态下几乎没有抵抗变形的能力，因此只需考虑体积的变化。为了考虑这一变化，在本例的有限元计算中采用Gruneisen状态方程。

Gruneisen状态方程的参量使用文献[9]中的数据，如表2所示。

表2 Gruneisen状态方程参量Tab.2 Gruneisen EOS parameters

1.3.3 岩土的动态力学性能

不同于金属材料，岩土材料的力学性能需要考虑静水压力（岩石上覆静压力）的影响。目前已经提出了20 多种屈服强度准则，考虑到岩土受碰撞冲击的影响，本文采用Holmquist-Johnson-Cook（HJ-C）模型来模拟岩土的动态力学性能[10]。

H-J-C 模型为一种率相关的损伤型本构关系，能够较好地描述材料在高应力、大应变和高应变率情况下的动态响应，在岩石冲击数值仿真中有着广泛的应用，LSDYNA 也引入了该本构模型。H-J-C 模型由状态方程、屈服面方程和损伤演化方程三块组成。状态方程（图1）可以分为三个阶段：第一阶段为线弹性阶段；第二阶段为过渡阶段，表示材料内部的空隙被完全压缩进塑性阶段；第三阶段为完全密实阶段，表示材料内的空隙已经被完全压碎，用一个三次多项式来表示，即

H-J-C模型的屈服面方程为

图1 H-J-C状态方程Fig.1 H-J-C equation of state

式中，σ*和p*为特征化等效应力和特征化等效压力；ε̇*为特征化应变率；D为损伤度；A、B、N和C为材料的强度参数。

H-J-C 模型为损伤型材料本构模型，采用等效塑性应变和塑性体应变累计来描述损伤，损伤的演化方程为

本文所用的参数设置如表3所示，表中ρ0为密度，fc为抗压强度，pc和uc分别为压垮时的静水压力和体积应变，pl和ul分别为压实时的静水压力和体积应变。

1.4 碰撞计算模型

采用的模型为某深海潜器载人舱耐压壳体，为球柱组合结构，在壳体内部均布内肋骨，外部设有外肋骨。模型的主要设计参数如表4所示。

表4 载人舱模型主要参数Tab.4 Main parameters of the model

载人舱模型为两端带球封头的环肋圆柱壳结构，如图2所示。由于该模型的厚度半径比为0.025，属于薄壳结构，故壳体部分采用shell 163 单元建模。被撞厚板采用solid 164 单元建模。碰撞主要发生在艏部球形封头处，因而对艏部壳体网格进行局部细化。为了简化计算，模型采用四分之一建模。通过在LS-DYNA 程序中定义加载曲线的方式，对耐压壳体进行动力加载，使得水压在0.09 s内从0 逐渐加载至15 MPa 后保持不变，再对耐压壳体进行显式动力学计算。壳体的初始运动速度为2.5 kn，与一块5×5×0.2 m3的方形花岗岩厚板发生碰撞。碰撞前初始的间距为0.125 m，确保在碰撞前完成加载过程。

图2 耐压壳体模型Fig.2 Model of pressure hull

2 水下碰撞动力分析过程

2.1 碰撞区域损伤变形分析

图3 碰撞前应力分布情况Fig.3 Stress situation before collision

图4 碰撞后应力分布情况Fig.4 Stress situation after collision

图5 结构的塑性应变Fig.5 Plastic deformation of structure

域的结构应力情况与碰撞前相差不大。结构的塑性应变情况如图5所示。随着撞深的增加，塑性应变以撞击点为圆心逐渐增大，但主要集中在艏部碰撞区域内。这体现了碰撞过程具有局部性的特征，提升结构的耐撞性能可以从增强易发生碰撞区域（艏部或者舷侧）的结构来实现，整体结构加强没有必要。

2.2 碰撞力分析

图6为碰撞力时历曲线，显示了碰撞力随时间变化的曲线。碰撞力的每一次卸载过程都代表了结构失效或者破坏，在0.114 s时碰撞力达到的峰值为2.43×106N；碰撞持续时间为0.058 s，冲量为9.4×104N∙m，碰撞力均值为1.62×106N。定义峰值与均值的比值为曲线的斜度K，则K=1.5。而正弦函数的斜度K约为1.5，因此碰撞力时历曲线可以用正弦函数来表示，其中fp为碰撞力峰值，Δt为碰撞持续时间：

图6 碰撞力时历曲线Fig.6 Time history curve of collision force

2.3 能量分析

图7 和图8 分别为碰撞过程中耐压壳体的动能和内能时历曲线。由图7 可知耐压壳体的动能先减小至0，随着反弹又有增大，动能为0时变形能最大。由于静水压力加载过程为线性加载，耐压壳体碰撞前的内能随着时间逐步增加，内能以弹性变性能的形式储存在壳体上。碰撞过程损失的动能为29.4 kJ，而耐压壳体的内能增加量为70 kJ，后者远大于前者。分析可知，碰撞过程中一直承受着15 MPa 的静水压力，由于静水压力作功，系统的机械能并不守恒（整体增加）。在整个过程中，由于被撞物体刚度较大，吸收内能很少，静水压力在碰撞区域内作功（占大部分）和耐压壳体的动能减小量转化为耐压壳体的内能，以结构变形能和断裂能等形式储存在耐压壳体上。由此可见，水下碰撞是一种极其危险的碰撞，会对深海潜器自航带来巨大的威胁。

图7 动能时历曲线Fig.7 Time history curve of kinetic energy

图8 内能时历曲线Fig.8 Time history curve of internal energy

2.4 剩余承载能力分析

2.5 失效模式分析

由以上分析可知，载人舱耐压壳体在1 500 m 水深以2.5 kn 的速度发生碰撞时会产生大量的局部塑性变形，但整体不发生破坏。进一步计算发现，随着碰撞速度的增加，达到某一个临界速度时，耐压壳体结构在静水压力和动态集中力的联合作用下发生失稳，结构整体破坏。图9显示了速度为1.3 m/s和1.4 m/s 时的碰撞结果。以1.3 m/s 碰撞时，在艏部碰撞部位产生了大量的塑性变形，但结构整体保持稳定。速度达到1.4 m/s 时，结构突然失稳，艏部大量凹陷，此时可认为耐压壳体在碰撞力和静水压力共同作用下发生了动屈曲，屈曲临界速度为1.4 m/s。图10是速度为1.3 m/s和1.4 m/s时艏部节点位移的时历曲线，反映了结构失稳时迅速压溃。

图9 速度为1.3 m/s和1.4 m/s时艏部变形情况Fig.9 Bow deformation at 1.3 m/s and 1.4 m/s

图10 速度为1.3 m/s和1.4 m/s时艏部节点位移的时历曲线Fig.10 Time history curve of bow node displacement at 1.3 m/s and 1.4 m/s

3 撞击参数对水下碰撞响应的影响

耐压壳体水下碰撞过程受碰撞速度、碰撞角度、静水压力、被撞结构形状等撞击参数的影响。本文主要研究环肋球柱组合结构与海山、礁石等深海固定物的碰撞情况，因而选取碰撞速度和静水压力作为主要的研究参数。海山凸出部分的形状选取三种比较典型的形状对比研究。撞击角度取为90°，即只考虑正撞的情况。

3.1 静水压力对碰撞的影响

深海潜器工作的水深不同，承受的静水压力也不同。为研究静水压力对水下碰撞过程的影响，结构模型与材料、碰撞速度等参数不变，静水压力改为0、5 MPa 和10 MPa，分别对应水面、水深500 m 和水深1 000 m处，与水深1 500 m的碰撞结果进行对比。

图11所示为耐压壳体在不同水压下发生碰撞时碰撞区域的塑性变形情况。随着水压的增加，塑性变形的大小和范围也随之增大。图12 为相同速度不同静水压力下碰撞力的时历曲线。图中各条曲线的变化规律大致相同，波峰位置基本一致，碰撞力峰值与静水压力呈负相关。这是由于外界静水压力的存在降低了耐压壳体的碰撞刚度，且静水压力越大，刚度降低得越多。

图11 不同水深下碰撞区域的塑性变形情况Fig.11 Plastic deformation in collision zone at different depths of water

从损伤变形的时历变化中，可以发现随着静水压力的增大，结构损伤变形的大小和范围增大。这是由于静水压力的存在提高了结构的整体应力水平，同时静水压力在碰撞过程中作功扩大了损伤的范围。

图12 不同水深下碰撞力的时历曲线Fig.12 Time history curve of collision force at different depths of water

3.2 撞击速度对碰撞的影响

深海潜器以不同的速度自航时，发生碰撞的动力响应也不同。结构模型与材料、静水压力等参数不变，将速度改为1 kn、1.5 kn、2 kn和2.5 kn时的碰撞情况进行对比。

图13所示为耐压壳体在不同速度下发生水下碰撞的结构损伤情况。随着碰撞速度的增大，塑性变形的大小和范围也增大。图14为碰撞力在相同水深、不同速度下的时历曲线。由图可知碰撞力峰值与速度呈正相关，随着碰撞速度的增大，曲线震荡更加剧烈，这可能是结构快要发生失稳（动屈曲）所致。

图13 不同速度下碰撞区域的结构损伤情况Fig.13 Plastic deformation in collision zone at different velocities

图14 碰撞力在不同速度下的时历曲线Fig.14 Time history curve of collision force at different velocities

3.3 被撞结构形状的影响

由于海底环境复杂，海底固定物形状多种多样，与深海耐压装备发生水下碰撞的结构相撞也各有不同。因此，为了研究被撞结构形状对水下碰撞的影响，选取三种比较典型的结构形式：平板、四棱锥和半球形，如图15所示。

图15 被撞结构形状Fig.15 Structure shape of the impacted object

由表5对比可知，三种情况下均不发生整体破坏，按碰撞力大小依次为平板、半球和四棱锥,这是由于平板的整体碰撞刚度要远大于四棱锥的整体碰撞刚度，而半球又大于四棱锥。观察耐压壳体随时间变化的塑性变形情况可以发现，四棱锥的最大塑性变形最大，其次是半球型，最后是平板。对碰撞后的剩余承载能力进行分析，由剩余承载能力因子发现碰撞对结构的破坏程度按大小依次为四棱锥、半球形和平板。

表5 被撞结构形状的影响Tab.5 Effect of structure shapes on structural impact

4 结 论

（1）具有自航能力的深海潜器进行深海作业时，耐压壳体自身处于高应力状态，最大应力状态接近于材料屈服极限，艏部发生碰撞时会产生大量塑性变形。随着撞深的增加，最大应力部分向四周扩散，非碰撞区域的结构应力情况与碰撞前相差不大，塑性变形集中在碰撞区域，碰撞具有局部性的特征；若碰撞速度达到某一临界值时结构在动态集中力和静水压力的共同作用下会发生动屈曲，结构迅速崩溃。

（2）水下碰撞过程中时刻承受着较高的静水压力，因而在碰撞前很大的初始内能是以壳体内的弹性应变能的形式存在。碰撞过程中由于静水压力作功，系统的机械能整体增加。在整个过程中，静水压力在碰撞区域内作功（占大部分）和耐压壳体的动能减小量转化为耐压壳体的内能，以结构变形能的形式储存在耐压壳体上。

（3）静水压力越大，产生的塑性变形越大，结构损伤范围扩大，损伤变形提前，发生动屈曲的风险越大；碰撞力峰值与静水压力呈负相关，这是由于外界压力的存在降低了耐压壳体的碰撞刚度；碰撞速度越大，塑性变形的大小和范围越大，碰撞力峰值与速度呈正相关，随着碰撞速度增大，曲线震荡更加剧烈，这可能是结构快要发生失稳（动屈曲）所致。

（4）被撞结构形状对碰撞载荷有很大的影响，按大小依次为平板、半球和四棱锥。对碰撞后的剩余承载能力进行分析，由剩余承载能力因子发现碰撞对结构的破坏程度按大小依次为四棱锥、半球形和平板。