球鼻船艏斜撞船舶舷侧结构的变形机理研究

王泽平， 胡志强， 刘 昆， 陈 刚

(1.上海交通大学 海洋工程国家重点实验室，上海 200240；2.哈尔滨工程大学 船舶工程学院，哈尔滨 150001；3.英国纽卡斯尔大学 工程学院，英国纽卡斯尔 NE1 7RU；4.江苏科技大学 船舶与海洋工程学院，镇江 212003；5.中国船舶及海洋工程设计研究院，上海 200011)

随着海上航行船舶数量的不断增加，船舶发生碰撞事故的概率也逐渐增加。尽管很多先进的避撞设备已经应用到船舶航行中，但船舶发生碰撞事故的风险仍然存在；船舶碰撞事故会导致环境污染、经济损失甚至人员伤亡。有统计数据表明，船舶发生斜撞事故的概率要高于船舶正撞事故[1]。然而目前对于船舶碰撞事故的研究主要集中在船舶正撞场景，对于船舶斜撞场景开展的研究相对较少，因此需要针对船舶斜撞场景开展相关结构耐撞性研究。

船舶碰撞的研究方法可以分为四种，包括经验公式法，试验法，有限元数值仿真法和解析计算法。Minorsky[2]作为经验公式法的开创者，通过对已发生的船舶碰撞事故结构损伤进行统计分析，进而提出内能吸收与结构损伤体积的经验公式。试验法分析船舶碰撞问题的准确性较好，不过研究成本较高。Pedersen等[3]通过试验法研究船舶结构的耐撞性；Amdahl等[5]通过开展模型试验研究船舶碰撞与搁浅时的结构变形机理；Liu等[6-7]通过落锤冲击试验研究了应变率的影响。随着计算机技术的发展，有限元数值仿真法被广泛应用于船舶碰撞问题的研究中。王自力等[8]通过数值仿真法研究了船舶碰撞动力学问题。刘敬喜等[9]通过数值仿真法研究了双壳油船舷侧结构的耐撞性。刘昆等[10]采用全耦合技术分析了船体结构的碰撞性能。此外，解析计算法是一种成本低且可以提供理想精度的方法，其通过提出结构的理论变形模型，运用塑性力学等理论，推导得到解析计算公式。Alexander[11]创造性地将解析计算法应用到薄板结构在轴向加载情况下的分析计算中。高振国等[12]通过解析计算法分析了强肋框结构在面内加载时的变形机理。Hong等[13]研究了船舶碰撞与搁浅过程中桁材结构的变形机理，并提出了解析公式。

近年来，随着船舶吨位的增加，船舶球鼻艏的安装也更加普遍，例如补给船，穿梭油轮等。船舶球鼻艏结构通常比较坚硬，发生船舶碰撞事故时更容易引起被撞船舷侧结构的损伤，进而导致舱室进水、漏油等事故的发生。为了减少因为球鼻艏撞击导致的损失，研究舷侧结构在球鼻艏撞击场景下的结构响应具有重要意义。然而大多数已有的研究都集中在球鼻艏正撞舷侧结构的场景，比如，Lee等[14-15]研究了舷侧外板在球鼻艏正撞场景下的结构响应。对于球鼻艏斜撞场景下舷侧结构耐撞性的研究较少。在斜撞场景下，舷侧结构的变形情况与撞击角度相关性很大，舷侧结构的主要构件，例如舷侧外板、横向肋板等构件在斜撞场景下的变形模式与正撞场景是不同的，因此需要针对舷侧结构在球鼻艏斜撞场景下的变形机理开展研究。

1 解析计算基本理论

解析计算方法可以快速评估被撞船的耐撞性，通过应用塑性力学理论，求解舷侧结构主要构件在球鼻艏斜撞场景下的抵抗力。根据上限定理，外力对结构的做功功率等于结构内部构件的能量耗散率[16]

(1)

(2)

(3)

(4)

N0=σ0t

(5)

式中：σ0为流动应力，t为舷侧外板的厚度。

2 舷侧外板的变形机理研究

2.1 舷侧外板几何变形模型

椭球形撞头斜撞一块边长分别为2a和2b的矩形板，如图1所示。该椭球形撞头可以用如下方程的形式表示

(6)

图1 椭球形撞头斜撞矩形板示意图Fig.1 A rectangular plate struck by an ellipsoid indenter

假设截面为椭圆的撞头撞击矩形板，如图2所示。该椭圆截面是椭球形撞头在xoy平面内的截面，那么椭圆形截面可以用如下方程的形式表示

(7)

图2 截面为椭球形的撞头斜撞矩形板示意图Fig.2 A rectangular plate struck by an indenter (b1+b2=2b)

矩形板的变形模式如图3所示。在矩形板受压变形的过程中，左侧部分矩形板的应变可以用如下方程表示

图3 矩形板的变形模式示意图Fig.3 Deformation model of the rectangular plate

因此，变形矩形板左侧部分的应变率是

(8)

式中:l1是矩形板左侧部分未发生弯曲部位的瞬时长度;α1是矩形板左侧部分的瞬时旋转角度。

(9)

同样，右侧部分矩形板的应变可以表示为

(10)

变形矩形板右侧部分的应变率是

(11)

将式(9)和(11)代入式(3)中，可以得到膜拉伸能量耗散率为

(12)

N0_sideplate=σ0tp

(13)

式中：N0_sideplate是塑性膜拉伸力；tp是矩形板的厚度。

椭球形撞头的撞击深度可以用下式表示

(14)

其中，b1,b2的长度可以分别表示为

b1=x1tanα1+l1cosα1

(15)

b2=x2tanα2+l2cosα2

(16)

式中，(x1,y1)和(x2,y2)是矩形板未弯曲部分与弯曲部分的交点。x1,x2分别为变形矩形板未弯曲部分与弯曲部分的交点在x轴方向上的坐标值，如图3所示。x1和x2可以表示为

(17)

(18)

由此可以得到椭球形撞头的撞击深度

(19)

椭球形撞头的撞击深度对瞬时旋转角度求导，可以得到撞击速率为

(20)

由于矩形板的弯曲能量耗散占总能量耗散的比例很小，所以通常忽略弯曲能量耗散。因此，矩形板在椭球形撞头撞击时的瞬时阻力可以表示为

(21)

同样，假设截面为椭圆的撞头撞击矩形板，该椭圆截面是椭球形撞头在xoz平面内的截面，那么椭圆形截面可以用如下方程的形式表示

(22)

矩形板受到撞击时的瞬时阻力计算方法与上文提出的方法相同，可以表示为如下形式

(23)

Haris等的研究发现椭球形撞头撞击矩形板的瞬时阻力与上文提出的矩形板受到两个方向上的椭圆形截面撞头撞击时的瞬时阻力关系可以表示为

(24)

在斜撞过程中，椭球形撞头和矩形板之间的摩擦力可以表示为

FS=μ·Fp·sinβ

(25)

由于摩擦产生的能量耗散可以表示为

(26)

2.2 舷侧外板破裂的预报

在碰撞过程中，如果外板发生破裂，那么外板的抵抗力会出现明显下降，因此能够较为准确地预测出外板的破裂十分有必要。通常临界断裂应变被当做板受到面外载荷作用时发生膜拉伸的最大应变，即认为当板发生膜拉伸的应变值达到临界应变值时，板出现破裂。这一方法被广泛使用，因此本文也采用这种方法预测外板的破裂。低碳钢的塑性拉伸应变值通常在0.2～0.35，在船舶的初步设计阶段，设计师可以根据具体情况确定这一数值的大小。将式(8)代入式(19)中，可以得到舷侧外板最大膜拉伸应变和球鼻艏撞深之间的关系

(27)

式中，εm为外板的临界断裂应变值。根据式(27)可以得到舷侧外板发生破裂时球鼻艏的撞击深度。

2.3 舷侧外板破裂后的解析计算公式

随着碰撞过程的进行，球鼻艏撞深增加，舷侧外板将会发生破裂。Wang等分析了板的变形过程，提出了预测板破裂后的阻力计算公式

(28)

式中：n为板破裂后的裂纹条数，n>2；l为裂缝的长度；φ为球鼻艏顶角的一半；μ为摩擦因数。

3 桁材变形模型

为了研究桁材结构在不同斜撞角度下的变形机理，首先对桁材的斜撞过程进行了数值模拟。选取三种典型的斜撞角度，数值模拟了球鼻艏结构在30°、45°和60°三种场景下斜撞桁材的过程。桁材在30°、45°和60°场景下的有限元模型、变形模型以及横截面的变形过程如图4～6所示。

(a) 有限元模型

(b) 变形模型

(c) 横截面变形过程图4 桁材在数值模拟中的变形(β=30°)Fig.4 Deformation of the web girder in numerical simulation(β=30°)

(a) 有限元模型

(b) 变形模型

(c) 横截面变形过程图5 桁材在数值模拟中的变形(β=45°)Fig.5 Deformation of the web girder in numerical simulation(β=45°)

(a) 有限元模型

(b) 变形模型

(c) 横截面变形过程图6 桁材在数值模拟中的变形 (β=60°)Fig.6 Deformation of the web girder in numerical simulation

从图4～6中可以看出，在不同的斜撞角度下，桁材腹板的折叠高度比是相似的，因此根据数值模拟的结果提出了桁材的理论变形模型，如图7和图8所示。其中，b1为受压桁材结构的一侧长度，另一侧长度为b2=b-b1。

图7 桁材撞击后的理论变形模型Fig.7 The deformation model of the web girder after collision

图8 横截面受压折叠变形过程(黑点代表塑性铰所在位置)Fig.8 Folding deformation process of the cross section (the black dots represent the location of the plastic hinges)

根据所提出的桁材理论变形模型，通过推导可以得到适用于不同斜撞角度的桁材变形阻力计算公式。限于篇幅，详细推导过程请参阅Wang等[17]9-10。推导后得到的桁材结构在形成第1个褶皱时的瞬时阻力公式为

(29)

在形成第1个褶皱时的平均阻力公式为

(30)

根据上限定理，当平均阻力最小时，能量耗散也最小，即：

(31)

将式(30)代入式(31)中，可以得到：

(32)

同样，当桁材结构形成第n个褶皱时(n≥2)，瞬时阻力公式为

(33)

式中：δ′为球鼻艏在斜撞方向上的位移距离；t为腹板的厚度。

桁材上带板的阻力公式可以采用梁理论进行计算

(34)

式中：a是带板宽度的一半；b1+b2是带板的长度；tf是带板的厚度。

4 横向肋板变形模型

球鼻艏斜撞横向肋板时，假设横向肋板受到长度

为d的线荷载作用，切向碰撞力使横向肋板上边缘产生水平位移，垂向碰撞力使横向肋板产生二维折叠变形，所提出的横向肋板理论变形模型，分别如图9和图10所示。

图9 横向肋板碰撞前后的横截面Fig.9 The cross section of transverse frame before and after collision

图10 横向肋板受压变形模式Fig.10 Deformation mode of the crushing transverse frame

根据所提出的横向肋板理论变形模型，通过推导可以得到适用于不同斜撞角度的切向力和垂向力的计算公式。限于篇幅，详细推导过程请参阅Wang等[17]10-12。推导后得到的切向力和垂向力的计算公式如下

(35)

(36)

式中：δ′为球鼻艏在斜撞方向上的位移距离；t为肋板的厚度。

5 数值仿真

本文通过开展数值仿真分析来验证所提出解析公式的合理性与准确性。

5.1 有限元模型

在数值仿真中，选取一艘双壳油轮作为被撞船，一艘60 000 t船舶的椭球形球鼻艏作为撞击船艏。撞击船艏和被撞船舷侧结构的有限元模型如图11所示，表1列出了被撞船的主尺度。

图11 油轮舷侧结构和撞击船艏的有限元模型Fig.11 Finite element model of the side structure of the tanker and the bow

表1 被撞船主尺度表Tab.1 Dimentions of the struck ship

使用有限元软件LS_DYNA进行数值仿真，使用Johnson-Cook GTN材料模型[18]。撞击船艏的材料属性定义为刚体，被撞船舷侧结构的材料参数如表2和表3所示。撞击速度定义为3 m/s，摩擦因数定义为0.3[19]。有限元模型中采用的是四节点四边形的Belytschko-Tsay单元。为了在保证数值模拟精度的同时有合理的计算时间，有限元模型的网格尺寸设置为200 mm。被撞船舷侧两端采用固定边界条件；本文研究的是船舶斜撞过程中相对危险的情况，即撞击船在碰撞过程中不发生偏转，被撞船舷侧结构在较大撞深情况下的损伤变形机理。在数值模拟中，被撞船和撞击船之间的接触设置为面面接触。

表2 GTN模型参数Tab.2 GTN parameters

表3 Johnson-Cook模型参数Tab.3 Johnson-Cook parameters

5.2 碰撞场景定义

船舶碰撞时斜撞角度的定义如图12所示，图13为三种典型的碰撞位置，两种碰撞场景的定义分别如表4和图14所示。工况1中的碰撞位置为位置1，工况2中的碰撞位置为位置2，随着碰撞的进行，撞击船艏会经过位置3。

图12 斜撞角度β的定义Fig.12 Definition of the collision angle β

图13 三种典型碰撞位置Fig.13 Three typical collision positions

表4 不同碰撞场景中的碰撞位置Tab.4 Definition of collision positions in different collision scenarios

(a) 工况1

(b) 工况2图14 两种斜撞场景的定义Fig.14 Definition of two oblique collision cases

5.3 碰撞过程中参与构件的确定

在5.2节中确定了碰撞场景，可以结合撞击球鼻艏结构的轮廓图和被撞船舷侧的结构图确定参与斜撞过程中的船体构件。图15(a)是工况1，撞击过程通过将撞击球鼻艏结构的轮廓图和被撞船舷侧的结构图结合起来进行模拟。在碰撞开始时，与球鼻艏接触的构件有舷侧外板和桁材，所以舷侧结构的抵抗力由这两个构件产生，舷侧外板的变形区域为CDHI；随着碰撞的进行，横向肋板也参与到舷侧结构抵抗力的计算中，舷侧外板的变形区域为CEGI，如图15(b)所示。图15(c)和(d)是工况2，在碰撞开始时，与球鼻艏接触的构件有舷侧外板和横向肋板，舷侧结构的抵抗力由这两个构件产生，舷侧外板的变形区域为AHJK，如图15(c)所示；随着碰撞的进行，舷侧外板的变形区域为FGJK，如图15(d)所示。

(a) 工况1

(b) 工况1

(c) 工况2

(d) 工况2图15 碰撞过程中变形区域及参与构件的确定Fig.15 Determination of deformation area and participating components during collision

6 数值验证与讨论

使用有限元软件LS_DYNA对两种斜撞场景进行数值仿真，将得到碰撞过程中球鼻艏撞击深度与碰撞力和能量耗散之间的关系曲线；解析方法中的能量耗散是通过对膜拉伸能量耗散率积分得到，将数值仿真结果与解析计算结果进行对比，如图16～图19所示。

图16 碰撞场景1中碰撞力-撞深曲线对比结果Fig.16 The comparison results of resistance-depth curves in case 1

图17 碰撞场景1中能量-撞深曲线对比结果Fig.17 The comparison results of energy-depth curves in case 1

图18 碰撞场景2中碰撞力-撞深曲线对比结果Fig.18 The comparison results of resistance-depth curves in case 2

图19 碰撞场景2中能量-撞深曲线对比结果Fig.19 The comparison results of energy-depth curves in case 2

从图16～图19的对比结果中可以看出，数值仿真方法和解析方法得出的结果吻合良好。随着碰撞过程的进行，碰撞力和能量耗散曲线呈上升趋势。当碰撞场景1和2中球鼻艏的撞击深度达到3.2 m时，舷侧外板发生破裂；舷侧外板发生破裂后抵抗力会下降，能量耗散曲线的增长速度也出现下降，解析方法和数值仿真方法都能捕捉到这一特点。在舷侧外板破裂初期，解析方法和数值仿真方法得到的碰撞力-撞深曲线存在一些偏差，这是因为在舷侧外板破裂初期，膜拉伸能量耗散仍然是舷侧外板最主要的能量耗散形式，舷侧外板的变形模式不会在很短的时间内转化为撕裂变形模式，这是解析方法和数值仿真方法得到的结果之间存在偏差的主要原因。总体来说，解析方法和数值仿真方法得到的结果对比良好，这也验证了本文所提出的解析公式的合理性。

7 结 论

本文研究了船舶球鼻艏结构斜撞船舶舷侧时，舷侧结构主要构件的变形机理，提出了舷侧外板、桁材和横向肋板三种构件的理论变形模型，运用塑性力学理论推导得到了三种构件在不同斜撞角度下变形过程中的解析计算公式。定义了两种典型的斜撞场景，通过有限元软件LS_DYNA进行了数值仿真，数值仿真结果验证了解析计算结果的准确性，解析计算公式可以很好地捕捉到舷侧外板发生破裂后碰撞力和能量耗散增长速度下降这一特点。

球鼻艏也是船舶常见的设计形式，船舶斜撞事故发生的概率要高于船舶正撞事故，而目前针对船舶碰撞的研究主要集中在船舶正撞场景，对于船舶斜撞场景的研究是不充分的，因此，研究舷侧结构在船舶球鼻艏结构斜撞场景下的变形机理是很有意义的。本文对舷侧结构在球鼻艏斜撞场景下的变形机理进行了研究，所得到的解析公式对于船舶舷侧结构耐撞性的评估具有较好的工程应用价值。