超大型浮式储油船的水弹性响应预报

翟钢军， 程 勇， 马 哲

(大连理工大学 深海工程研究中心,辽宁 大连 116024)

随着我国工业化进程的加快和国民经济的迅速发展，我国的石油消费量呈迅速上升的势头，而我国石油资源供应能力有限，因此借鉴国际先进经验，迅速建立国家石油战略储备，是最现实、最可靠的手段。目前石油储备主要是针对原油储备，其方式有5种[1]：地上浮顶油罐储备方式、地中油罐储备方式、地下岩洞储备方式、地下盐穴储备方式和海上石油储备方式。结合我国实际情况，采用以海上战略石油储备为主的多样化储备方式更加适合我国人口密度、地理环境和工业分布状况的特点。

图1 海上石油储备基地布置图

海上石油储备基地如图1所示，主要由数艘储油船并排靠泊组成，其中系泊系统由若干系缆桩和橡胶护舷构成。储油船属于超大型浮式结构物(VLFS)，即水平尺度与高度的比值非常大，是一个极为扁平的柔性结构物[2]。由于储油船的弹性变形与刚体位移几乎属于同一量级，计算结构在流体作用下的水动力特性时，不能再把船体当作单一的刚体进行计算，必须分析储油船在流体荷载作用下的底部变形，以及船体变形对流场的影响，这是一个典型的流固耦合问题，需用水弹性理论进一步分析。

水弹性力学成立于70年代，由Bishop等[3]提出，它需要考虑惯性力以及流体荷载和结构弹性力之间的耦合作用，关系到水动力学和结构动力学两门学科以及它们之间的联系，属于典型的交叉学科。日本、美国、英国和韩国等相继对超大型浮式结构物的水弹性进行了研究，推行了一系列VLFS研究计划，如Wu等[4]通过特征函数展开法来求解弹性薄板的二维水弹性问题；Hong等[5]分析了连接OWC防波堤的VLFS的频域水弹性响应和漂移力；Murai等[6]提出了一种子结构法来加速水弹性计算且不丢失原有精度；Fu等[7]预测了两个柔性连接的VLFS水弹性响应；Watanabe 等[8]提出了一个圆柱形VLFS的标准水弹性方法；Pham 等[9]进一步分析了配有减振装置的圆柱形VLFS的水弹性问题。我国也有很多学者从事该方面研究，如吴有生[10]发展了适用非常广泛的三维水弹性理论；崔维成等[2,11-12]总结了VLFS的研究现状和一些关键问题，并运用模态展开法分别计算了梁模型和矩形平板模型的水弹性响应；滕斌等[13]学者介绍了计算水弹性的模态展开法和直接计算法两种方法，并计算了五种振型函数下结构的水弹性动力学；张淑华等[14]学者采用直接法计算了弹性平板模型的挠度幅值。

对于超大型浮式储油船的水动力分析，上海交通大学的李文龙等[16]基于三维势流理论，建立了储油船的刚体运动方程，并计算了相关的水动力参数。然而储油船属于超大型结构，其弹性变形甚至比刚性位移更为重要，因此必须给予考虑。本文采用模态分析法，将船的水弹性位移表示成一系列模态振型的叠加，应用边界元法求出绕射势和一系列模态上单位速度运动引起的辐射势。再通过储油船的弹性体运动方程求出各模态的振型幅值，最后回代到水弹性位移表达式中，求出储油船的水弹性响应。在储油船的水弹性力学分析中，本文假设储油船为一薄板模型，通过与滕斌等[13]的数值结果及Maeda等[15]的实验结果进行对比，验证了本文方法的可行性，进一步分析了不同水深及波要素对船体运动响应及水动力参数的影响。

1 储油船的水弹性分析

考虑波浪与储油船相互作用时，储油船长为L，宽为B，吃水为d，且开敞水域水深为h，选定坐标系如图2所示，x-y平面位于未干扰自由水面，根据右手定则z轴垂直向上。整个流域Ω由储油船底面Sb，侧面Sc，自由水面Sf，无穷远柱面S∞和海底面Sd面组成。

储油船在流场中受流体作用过程中，包括结构和流场两部分计算，船受到流体荷载作用要产生变形，属于结构计算，而船体的存在和变形改变了流场，属于流场计算，反之，流场的变化也对船体的弹性变形产生影响，总之两者相互联系，相互耦合。其中对流域的求解可采用边界元法，对结构的弹性运动方程求解可采用直接法或模态分析法(本文采用该方法)，下面将对这两个方面进行详细论述。

图2 坐标示意图

1.1 流场控制方程和边界条件

假定流体为无旋、不可压缩的理想流体，则流体运动可以通过速度势Φ来表述。当考虑入射波浪为一频率为ω的规则波时，所有含时间变量t的物理量均可分离出时间分量，表达如下：

Φ=Re{φ(x,y,z)e-iωt}

(1a)

W(x,y,t)=Re{w(x,y)e-iωt}

(1b)

P(x,y,t)=Re{p(x,y)e-iωt}

(1c)

式中：W(x,y,t)为船的垂向运动响应；P(x,y,t)为船底部所受的流体压力。

进一步假定波浪的波陡和船体运动量足够小，满足线性势流理论。因此，复速度势φ满足以下边界值定解问题：

2φ=0 在整个流域

(2)

(3)

(4)

(5)

(6)

ω2=gktanh(kh)

(7)

式中：n为单位法向量，指出流域为正；r为参考点到坐标原点的水平距离；复值变量w(x,y)为储油船底部的垂直位移，包括刚体位移和弹性变形两部分。

1.2 储油船运动方程

对于超大型浮式储油船，水平尺度远大于垂直尺度，因此可将其简化为弹性薄板，则问题转化为求解弹性薄板在波浪作用下的水弹性响应。典型的薄板振动微分方程为：

D4W(x,y,t)=Q(x,y,t)

(8)

式中：D为薄板的弯曲刚度；Q为薄板单位面积上的横向荷载，包括惯性力和流体作用力两部分，表达如下：

(9)

式中：ms为薄板单位面积质量，W(x,y,t)和P(x,y,t)表达如上。将式(1b),式(1c)代入到式(8)中可得储油船的弹性运动方程为：

D4w(x,y)-ω2msw(x,y)=p(x,y)

(10)

对于自由漂浮的储油船，其端部边界上的弯矩和剪力为零：

(11)

(12)

式中：v是泊松比。

求出弹性位移后，储油船内部弯矩可求得为：

(13a)

(13b)

2 储油船的运动响应及速度势的模态展开形式

由储油船的弹性运动方程(10)可以看出，结构运动响应w(x,y)和流体运动是耦合在一起的，并且由物面条件式(5)可知流体运动可由结构运动响应得到。因此，将问题分解为速度势的水动力学问题和薄板振动的结构动力学问题，船的运动响应和速度势都进行模态展开。

2.1 振型函数的一般形式

对于储油船的振型函数可以表示成x和y两方向振型函数的乘积形式，因此以一维梁模型为例。振型函数的选取分为干模态法和湿模态法，本文采取前者方法。干模态法中，以在空气中自由振动的振动模态为自然模态。在水弹性问题中，自然模态包括两种表达形式，分别由Maeda等[15,17]推导得到。对于一些难以求出自然模态的复杂结构来说，Eatock Taylor等[18-19]推出了一种非常简单的振型函数表达式：

(14)

式中：f1(x)表示刚体振型模态，本文采取此振型函数进行储油船的水弹性分析。

根据振型函数的表达形式，可得储油船的运动响应函数w(x,y)表达式：

(15)

式中：ζj为j模态上的振型幅值。

考虑到线性势流理论的假定，复速度势φ包括入射波产生的入射势、结构存在形成的绕射势以及物体运动产生的辐射势之和，表示如下：

φ(x,y,z)=φI(x,y,z)+
φD(x,y,z)+φR(x,y,z)

(16)

(17)

(18)

2.2 速度势及水动压力分析

关于绕射势和辐射势边界值问题的数值求解，通常采用格林函数法。将绕射势φD和辐射势φj(j=1,….N)分别与格林函数代入到格林第二公式中，得边界积分方程：

(19)

其中(x0,y0,z0)和(x,y,z)分别为源点和场点坐标，G为满足除物面条件外其他所有边界条件的有限水深格林函数：

(20)

v=ω2/g为深水波数，

k0为色散方程ω2=gktanh(kh)的根，

J0为第一类贝塞尔函数。

积分方程(19)经过物面单元离散，可表示为离散线性方程组形式：

[A]{φj}={B}

(21)

由于储油船关于x-z和y-z两平面对称，此时系数矩阵及向量可表示为:

(22)

图3 对称区域简化示意图

各参数定义如图3所示，对于矩阵A中各元素，满足以下对称关系：

A11=A22=A33=A44

(23)

A12=A21=A34=A43

(24)

A13=A31=A24=A42

(25)

A14=A41=A23=A32

(26)

为减小矩阵维数，可采用以下转换:

(27)

其中矩阵R满足以下关系：

(28)

这样将式(27)和(28)代入线性方程组(21)可得：

(29)

根据矩阵[A]和[R]的特征，矩阵[A]是一个块对角矩阵，则方程(29)为:

(30)

进一步将线性方程组(30)转化为四个区域的求解形式，各区域表达式如下：

区域1：

(31)

(32)

(33)

(34)

区域2：

(35)

(36)

(37)

(38)

区域3：

(39)

(40)

(41)

(42)

区域4：

(43)

(44)

(45)

(46)

求得速度势后，储油船底板所受的压力分布p(x,y)为：

p(x,y)=iωρ[φI+φD+φR]-ρgw(x,y)

(47)

2.3 储油船的水弹性运动方程求解

由边界元法求出速度势及底板压力分布后，将(15)和(47)式代入薄板的振动微分方程(10)可得：

(48)

对上述方程采用Galerkin积分法，即方程两边同时乘以fi(i=1,2,…,N)，并沿储油船底部积分可得：

(49)

3 数值结果分析

根据上述数值方法，本文采用Fortran语言进行程序开发，并将计算结果与已知实验和数值结果进行比较，验证了本文程序的可行性。进一步利用所开发的程序，对储油船在不同水深及波长下的水弹性响应及变形、弯矩、广义波浪激励力和广义水动力系数进行了分析。

3.1 程序的验证

为了验证本文数值方法及程序的可行性，选用Maeda等[15]实验模型作为算例，并将计算结果与实验结果以及文献[13]的数值结果进行比较，如图4所示。实验模型的主要参数如表1所示。为了保证数值收敛和准确性，通常1个波浪长内至少要包含14个单元，因此模型沿板长方向划分60个单元，板宽方向16个单元，板高方向1个单元。

表1 实验模型主要参数

图4为波长分别2 m和4 m时，薄板左侧即迎浪端的水弹性响应位移最大时的结果。由图可见，本文计算结果与实验及文献[5]结果基本吻合，只是在两端位移数值上，本文计算结果稍大。且由图可观察到，波长较小时结构变形程度要大于波长较大时的变形程度。

图4 波长为2 m和 4 m时薄板水弹性位移分布

3.2 广义波浪激励力和广义水动力系数的计算

本文所采用的储油船数值模型按原型比尺1∶100取参数值：长L=4 m，宽B=0.9 m，吃水d=0.25 m，水深h=0.4 m、0.7 m、1 m，泊松比v=0.3，弯曲刚度D=50.832 Nm。型长方向划分100个单元，型宽方向划分20个单元，水深方向划分5个单元。对不同水深条件下，运动方程(49)中储油船所受的广义波浪激励力Fi随相对波长变化进行了分析(取前6阶模态)，如图5所示。由图可知，随着水深的增加，同阶模态激励力变化趋势基本一致，且1、3阶模态激励力基本不变，其他阶模态激励力缓慢降低。

图5 水深0.4 m,0.7 m和1.0 m时广义波浪激励力Fi计算结果图

图6 水深0.4 m时储油船广义附加质量和辐射阻尼计算结果图

此外，对不同水深条件下，运动方程(49)中储油船的广义附加质量Maii和辐射阻尼Bbii随波长的变化也进行了分析，如图6～8所示。其中横纵坐标取无量纲形式，广义附加质量和辐射阻尼都取到前6阶模态。由图可知，不同水深下的水动力系数变化趋势也不相同：当水深较浅时广义附加质量表现为低频特征，广义辐射阻尼表现为波频特征；随着水深的增加，广义附加质量也逐渐表现为波频特征；当水深增加到一定程度时，水动力系数表现为波频和高频两种特征，因此水深的选取对储油船的安全稳定至关重要。从图中可发现前6阶模态水动力系数的累加基本可近似表达总的水动力系数值，且第1阶模态即刚性模态水动力系数所占比例最大，第2阶模态的水动力系数所占比例较小，第4阶模态的广义附加质量在无量纲波数0.5～1.25时高于第3阶模态值，其余模态水动力系数依模态阶数降低。

3.3 不同水深下储油船的水弹性响应及变形

进一步研究不同水深条件下，储油船在波长分别为0.8 m和1.0 m时的水弹性响应，结果如图9所示。由图可得，当水深从0.4 m～0.7 m变化时，储油船水弹性响应逐渐减小，但当水深增大到1.0 m时水弹性响应突然增大，并且水深为0.4 m～0.7 m时，结构水弹性响应对波长的敏感度小于水深1.0 m时结果。另外，由图可知储油船的迎浪端水弹性变形比较大，这一点在设计时必须引起注意。

本文计算了整个储油船底板的弹性变形，下面给出了随着水深变化，储油船迎浪端出现最大响应时的弹性变形图，如图10～13所示，其中左边为波长0.8 m时计算结果，右边是波长1.0 m计算结果。由图可发现，水深为0.4～0.7 m时，结构水弹性响应对波长的敏感度要小于水深为1.0 m时的结果。

图9 波长分别为0.8 m和1.0 m时储油船的水弹性响应

图10 水深0.4 m时储油船底板变形图

图11 水深0.55 m时储油船底板变形图

图12 水深0.7 m时储油船底板变形图

图13 水深1.0 m时储油船底板变形图

图14 波长为0.8 m和1.0 m时储油船的弯矩响应分布

3.4 储油船弯矩响应分布

根据方程(13)计算不同水深时整个储油船的弯矩分布，如图14所示，其中上图为波长0.8 m时的弯矩响应分布，下图为波长1.0 m时的弯矩响应分布。由图可知，当水深从0.4～0.7 m变化时，整个储油船的弯矩分布逐渐降低，到水深0.7 m时基本稳定，但当水深增加到1.0 m时，弯矩又突然增大，因此0.7 m水深可作为储油船锚泊地点的相对理想水深；对于相同水深时，波长较大时的弯矩分布峰值大于波长较小的情况。

4 结 论

本文采用模态分析法中的模态展开法，将储油船的水弹性响应及辐射势用一系列模态叠加表示，再应用基于三维势流理论的边界积分方程和薄板振动微分方程求出各模态幅值，最后回代入储油船的运动表达式得到水弹性。主要得出以下结论：

(1)通过与已知实验及数值结果比较，验证了本文方法能够比较准确地预报超大型浮体水弹性响应。

(2)随着水深的增加，储油船所受的同阶模态广义波浪激励力变化趋势基本一致，但广义水动力系数变化趋势并不相同。广义附加质量随着水深的增加由低频特征变为波频特征最后同时具有波频和高频特征，广义辐射阻尼由波频特征变为同时具有波频和高频特征。

(3)随着水深的增加，储油船的水弹性响应、整个结构变形和弯矩响应分布均随之减小，到一定水深时趋于稳定，但当水深继续增大时，水弹性响应和弯矩会增加，因此超大型浮体系泊水深的选取对于水弹性响应计算具有重要影响。

(4)整个储油船的变形程度随着入射波波长的增加而减小，这是由于船底部流体压力分布的变化程度随着波长的增加而减小造成的。

参 考 文 献

[1]李文龙，刘亚东，谭家华. 基于模糊数学理论的我国战略石油储备方式决策[J].上海交通大学学报, 2004, 38(11)：1892-1896.

LI Wen-long, LIU Ya-dong, TAN Jia-hua. Decision-making for the facilities of strategic pectroleum reserve in china based on fuzzy mathematics theory[J]. Journal of Shanghai Jiao Tong University,2004,38(11):1892-1896.

[2]崔维成. 超大型海洋浮式结构物水弹性响应预报的研究现状和发展方向[J]. 船舶力学，2002，6(1)：73-90.

CUI Wei-cheng. Current status and future direction in predicting the hydroelastic response of very large floating structures[J]. Journal of Ship Mechanics,2002,6(1):73-90.

[3]Bishop R E D, Price W G, Wu Y S. A general linear hydroelasticity theory of floating structures moving in a seaway. Phil[J]. Trans. R. Soc. Of Londan (series)A. 1986, 316:375-426.

[4]Wu C， Watanabe E， Utsunomiya T. An eigenfunction expansion-matching method for analyzing the wave-induced responses of an elastic floating plate[J]. Applied Ocean Research, 1995, 17(5):301-310.

[5]Hong D C, Hong S Y. Hydroelastic responses and drift forces of a very-long floating structure equipped with a pin-connected oscillating-water-column breakwater system[J]. Ocean Engineering, 2007，34(5-6):696-708.

[6]Murai M, Kagemoto H, Fujino M. On the hydroelastic responses of a very large floating structure in waves[J]. Marine Science and Technology, 1999, 4(3):123-153.

[7]Fu S X, Moan T, Chen X J, et al. Hydroelastic analysis of flexible floating interconnected structures[J]. Ocean Engineering, 2007, 34(11-12):1516-1531.

[8]Watanabe E, Utsunomiya T, Wang C M, et al. Benchmark hydroelastic responses of a circular VLFS under wave action[J]. Engineering Structures, 2006, 28(3):423-430.

[9]Pham D C, Wang C M, Utsunomiya T. Hydroelastic analysis of pontoon-type circular VLFS with an attached submerged plate[J]. Applied ocean Research, 2008, 30(4):287-296.

[10]Wu Y S. Hydroelasticity of floating bodies[D]. Brunel University, UK,1984.

[11]金晶哲，崔维成，刘应中. 预报超大型浮体水弹性响应的模态函数展开方法和特征函数展开方法比较[J]. 船舶力学，2003，7(4)：86-98.

JIN Jin-zhe, CUI Wei-cheng, LIU Ying-zhong. Comparison of modal function expansion method with eigenfunction expansion method for prediction of hydroelastic responses of VLFS[J]. Journal of Ship Mechanics,2003,7(4):86-98.

[12]闫红梅，崔维成，刘应中.超大型浮体水弹性分析的平板林函数法[J].海洋工程，2003, 21(4)：8-21.

YAN Hong-mei, CUI Wei-cheng, LIU Ying-zhong. Hydroelastic analysis of VLFS using Plate Green Function Method[J]. The Ocean Engineering, 2003, 21(4):8-21.

[13]滕 斌，勾 莹.大型浮体水弹性作用的频域分析[J].工程学，2006, 23(s2)：36-47.

TENG Bin, GOU Ying. Hydroelastic analysis of very large floating structure in frequency domain[J]. Engineering mechanics,2006,23(增2):36-47.

[14]张淑华，韩满生.用直接法分析超大型浮体的水弹性响应[J].海洋工程，2004, 22(1)：9-18.

ZHANG Shu-hua, HAN Man-sheng. Analysis of hydroelastic response of VLFS by use of the direct method[J].The Ocean Engineering, 2004, 22(1):9-18.

[15]Maeda H, Masuda K, Miyajima S,et al. Hydrielastic reponsesof pontoon type very large floating offshore structure[J]. Jour of Soc Nav Archit Japan, 1995, 178：203-212.

[16]李文龙，谭家华.我国战略石油储备船系泊系统安全性评估[J].中国造船，2005,46(2):87-94.

LI Wen-long, TAN Jia-hua. Safety assessment of mooring system of oil storage vessels for strategic petroleum reserve in china [J].Shipbuding of China,2005, 46(2):87-94.

[17]Gran S. A course in ocean engineering[M]. Elsevier, Amsterdam,1992.

[18]Eatock Taylor R, Ohkusu M. Green functions for hydroelasticanalysis of vibrating free-free beams and plates[J]. Applied Ocean Research, 2000, 22：295-314.

[19]Eatock Taylor R. Wet or dry modes in linear hydroelasticity-why modes [C]. Hydroelasticity in Marine Technology, Oxford, UK, 2003：239-250.