含分布浮体的复杂构型深海缆线的动响应及其时空演化1)

郭双喜 陈伟民 ,2) 严定邦 宋吉祥 † 沈义俊

* (中国科学院力学研究所流固耦合系统力学重点实验室,北京 100190)

† (中国科学院大学工学院,北京 100049)

** (湖南应用技术学院,湖南常德 415100)

†† (中国水利水电科学研究院流域水循环模拟与调控国家重点实验室,北京 100038)

*** (海南大学南海海洋资源利用国家重点实验室,海口 570228)

引言

随着海洋油气、采矿工业向深海的发展,各种深水立管和输送管线系统相继问世,由于其结构复杂和较高的经济成本,用于油气、矿产输送的立管和电力输送缆线已成为整个生产系统的重要组成部分[1-3].而随着水深增加,已有深水立管面临新的挑战,例如SCR (steel catenary riser)上端点张力迅速增大,限制了其在深水海域的应用[4-6],TTR (top tension riser)由于顺应性较差,上部船体运动会使立管中产生很大张力,而HTR (hybrid tension riser)的接头部分有较为难处理的问题.为满足深水海域作业要求,一些新构型的立管形式被提出,包括油气行业的顺应式立管[5-6]和深海采矿系统的马鞍形立管[7]等,为了适应超深水油气工业和采矿的需求,近年来一些学者提出了Stepped 型和Double-stepped构型立管[6].为了减小结构上的张力,这些新构型深水管线多采用了分布式浮力模块安装形式,浮力模块的存在使深海柔性结构成为非均匀结构,使其在流体作用下的运动响应和应力分布比均匀结构更加复杂[8-11].为了保障立管系统的结构安全、有效运行和足够使用寿命,需要综合考察立管在环境荷载和移动边界条件激励下的张力、应力和位移等结构响应,提供结构详细设计和分析的依据.

针对SCR,TTR,lazy-wave riser 等构型的立管设计、响应分析已有丰富研究结果[8,12-18].Yin 等[14]使用数值模拟和模型实验两种方法研究了船体运动引起的TTR 立管的动响应,Cheng 等[15]提出了一种研究懒波型立管在船体运动和波流载荷作用下的动态特性的时域数值方法,并研究了浮体分布长度对立管动态特性的影响.唐达生等[16]建立实验室平台系统,测试了扬矿管道的扭转和弯曲振动特性,研究表明管道振动产生的动态应力远高于静态力,振幅放大系数约为6.75.Chatjigeorgiou[17]研究了水下可伸缩悬链型立管在顶部激励下的非线性动力学行为,结果表明只有平面外响应受内部流动的影响.Guo 等[8]研究了顶部船舶运动的幅值/频率以及浮力数沿结构长度的分布对懒波塔式立管和混合塔式立管动力响应的影响,Li 等[18]用有限元数值方法结合水动力模型研究了TTR 在线剪切流场中的多频涡激振动.

随着水深和立管跨度的增大,来自结构本身和周围流体的惯性和阻尼效应变得更加明显,这种动态行为会对结构的张力产生较明显影响,甚至当松弛发生时会导致顶部张力的急剧增加.Vassalos等[19]推导了运动缆索在水平激励下的非线性方程,Huang 等[20]提出了一种预测船用缆索在交替张紧-松弛条件下的拉断载荷的数值方法.考虑到一阶几何非线性,Mansour 等[21]得到了SCR 曲率的精确表达式,并定义了悬垂特征对SCR 运动的影响;Guo 等[22]研究了水下悬浮隧道运动的幅值和频率对系泊缆线位移和动张力的影响,在数值模型中考虑了几何非线性和非线性流体力.Zhang 等[23]发现紧-松状态会显著放大最大张力幅度,最大张力可能会增大5 倍.基于实验室实验,Hsu 等[24]发现由于锚链动力学的参与,冲击张力的最大值可能是原始张力的1.6 倍.近年来,由于有限元法善于处理复杂结构、边界条件,有良好的结构形状和边界条件适应性,开始被用于深水柔性结构的响应分析.Li 等[25]使用有限元法研究了SCR 系泊线的动态响应,发现考虑结构动态特性后,其回复刚度出现了明显的迟滞行为.吴天昊等[26]研究了平台运动和内流共同作用下立管的响应,研究发现平台运动引起的立管涡激振动具有不稳定性.

然而,已有的研究大多关注TTR,SCR 和懒波型立管的响应,针对超深水Stepped 型和Doublestepped 型立管的研究大多关注构型的静力问题[27].Mao 等[28]建立了实际隔水管构型的分析模型,研究了钻井隔水管的力学行为,发现构型对隔水管的力学行为,特别是对结构弯矩的影响很大,这主要是通过影响结构弯曲刚度、横截面张力来实现的.基于虚功原理和变分原理,Lou 等[2]推导了顺应式立管的静态几何非线性方程,并对立管的静态特性进行了综合分析.对于深水复杂构型立管而言,其长度方向分布有浮力模块,需要承受结构重力、不均匀浮力、波浪和水流等载荷,其运动响应和结构响应时空演化变得更加复杂,对结构安全带来严峻挑战.此外,顶端船体运动和底端采矿车作业产生的移动边界条件,并且动态效应变得更加明显,会给立管变形和动力响应带来附加影响.为保证立管的安全运行和使用寿命,需要进一步发展针对复杂结构形式和考虑动态效应的动力学分析模型,深刻认识立管结构动响应在结构上的空间分布和时间演化规律和机制.

本文考虑深海采矿系统作业条件下受到的波浪力、船体位置变化等载荷工况,以Double-stepped缆线为研究对象,基于其流固耦合载荷特性和模型表征,建立了考虑非均匀深水管线复杂结构形式和动态效应的动力学控制方程.基于有限元法和二次开发建立了缆线的动响应计算模型,分析了环境载荷下缆线的响应和传播规律;研究了海面船体位移引起的缆线空间位置、张力变化对结构响应和波传播特征的影响,并给出位移时空演化规律,最后基于WKB (Wentzel-Kramers-Brillouin)理论分析了变张力结构的响应幅值和波长演化规律.

1 分析模型

1.1 缆线流固耦合控制方程和数值模拟

Double-stepped (图1)这种构型缆线主要针对超深水油气、采矿工业背景提出,其特点是通过精细设计浮体的展向分布,整体构型具有多个拱弯段、垂弯段,浮体位置、浮力值的配比优化设计更加复杂,由于其构型特点,Double-stepped 缆线具有更优异的运动变形顺应性和对边界激励的缓冲效应,对运动边界有更好鲁棒性.由于复杂构型缆线分布有浮体装置,与以往均匀结构不同,首先需要对这种含浮体的非均匀结构的变形重构特性和流固耦合载荷进行表征.

取dS微段(如图1 所示),由y向力平衡条件建立方程

式中,V为剪力,m为单位长度质量,S为结构的长度,g为重力加速度,fbuoyancy为缆线受到的浮力,fbuoyball为浮体上的浮力和阻力.S*为浮体安装位置,当S≤S*时,δ (S-S*)=1,表示浮体引起的浮力和阻力仅对浮体安装位置以下的结构起作用;当S＞S*时,δ (S-S*)=0,即浮体引起的浮力和阻力对浮体安装位置以上的结构无影响.流体荷载采用经典的莫里森公式计算

式中,CD为拖曳力系数,U为结构与水的相对速度,D为结构的外直径,阻力计算时缆线的阻力系数取1.2,浮体结构的阻力系数取0.6.将上式代入式(1),简化可得

由力矩平衡条件对单元的微段右截面和s轴的交点取弯矩,可得

式中,ρs为结构的密度,ρb为浮力球的密度,Db为浮力球外径.对式(6)沿弧长s积分,并代入边界条件,化简可得含浮体结构的控制方程

求解上述方程得到结构的变形后,将沿水流方向的流体力分量在结构的长度上积分便可得到施加在结构上的总阻力

进一步引入重构数R,重构数定义为相同几何形状的柔性结构阻力与刚性结构的比值,强调了柔性对阻力的影响[24]

由式(10)可得到

通过式(9)～式(11)可以得到柔性结构上的载荷形式

式中,ν为Vogel 指数[29-30].得到Vogel 指数后,即可用式(13)表征作用在柔性结构上的与结构变形和运动相关的流固耦合载荷.

基于有限元法对缆线响应进行数值模拟,为了模拟柔性缆线的变形重构和流固耦合载荷作用,采用改进的梁单元离散结构,梁单元(图1)的位移向量

其中,u和v分别为节点的轴向位移和横向位移,θ为转角,i为梁单元的节点编号.可以将梁的位移分解为轴向位移和横向位移分别进行考虑

单元的质量矩阵和刚度矩阵

常规梁单元的质量矩阵由式(16)给出,而含分布附体的梁单元质量矩阵中需考虑浮体的作用,对应的刚度矩阵变为

式中,Ae为横截面面积,ρ为缆线的密度,M0为浮体的水下重力.将式(15)和式(16)中的矩阵组装便可得到单元的刚度矩阵Ke和质量矩阵Me,进一步的便可得到整个缆线结构的刚度矩阵和质量矩阵.由此可得出缆线的动响应控制方程

式中,M,C和K分别为结构的质量矩阵、阻尼矩阵和刚度矩阵,U和P为结构的位移矢量和外载荷矢量.流体阻力通过式(2)计算,结构受到的附加质量力由Morison 公式给出

其中,CA为附加质量系数,文中取1.0;af和as分别为流体和结构的加速度.在计算模型中与流体相关的载荷项通过计算流体加速度求得后作为外载荷加载,与结构加速度相关的项则可通过改变结构的线密度加载.考虑到问题的强非线性和流固耦合特征,采用Newmark 法求解缆线的动响应

式中,Δt为时间增量,β和γ为计算常数.由于深水缆线结构柔性很大,运动过程中会出现大位移、大变形响应,为避免刚体位移引起的矩阵奇异,本文通过有限元二次开发采用多次加载-卸载方法,先预加载张力再卸载得到缆线的构型.

1.2 模型实验

为验证上述分析模型,设计搭建了测试柔性缆线构型和响应的流固耦合实验平台(如图2 所示).实验平台主要由主体设备、测量系统组成,其中提供流场环境的主体设备包括XYZ三轴滑台牵引系统、高透光性水箱、整体支持框架等.通过滑台牵引水下柔性试件在X,Y,Z三个方向运动,可模拟不同方向水流与结构的相对运动,在运动过程中采集流固耦合载荷、缆线形状、运动位置等信息.

图2 柔性结构流固耦合响应实验平台Fig.2 Experimental platform for fluid-solid coupling response of flexible cable

平台主要参数:高透光性亚克力水箱1.8 m ×0.6 m × 1.1 m,XYZ滑台的X方向行程2.0 m,最大速度2.0 m/s,Y方向0.7 m 行程和0.4 m/s 速度,Z方向0.3 m 行程和0.4 m/s 速度.为了稳固支撑滑台系统且不影响水箱,专门设计了实验台框架,使用4080型标准铝材,框架长2.0 m、宽0.7 m、高1.2 m;为了避免系统共振,设计了实验台框架的固有频率至少高于实验试件自振频率三倍.测量系统包含测试件拖曳力的力传感系统和测试件变形系统.拖曳力测量系统由测力传感器、数据采集仪和电脑客户端软件实现;测量试件变形系统由高速摄像仪,视频编辑软件和图片处理软件三部分构成.

利用所搭建实验平台,分别测试了含分布浮体缆线的静态构型和顶部边界激励下的缆线响应,实验缆线模型如图3(a)所示,缆线参数见表1.缆线上下两端点通过铰支约束,浮体长5 cm、直径3 cm.图3(b)为缆线静态构型,通过图像处理提取缆线坐标点得到构型;静态构型实验结果与数值计算结果的对比如图3(c)所示,可以看出计算得到的缆线位置与实测值吻合良好.通过滑台控制缆线上端点做周期为5 s、幅值10 cm 的水平方向(X方向)周期运动,实验记录缆线顶端动张力并与数值结果进行对比.从图3(d)所示的动张力对比可以看出实验值和计算值基本吻合,数值结果略低于实验值;图3(e)所示的计算得到的重构数随CY数变化与实验值也基本吻合.

图3 柔性缆线响应数值和实验结果对比Fig.3 Comparisons between numerical and experimental results of flexible cable response

表1 实验缆线参数Table 1 Cable parameters of experiment

2 环境载荷作用下的缆线结构响应

考察环境载荷作用下缆线响应和海面船体运动对缆线响应的影响,分别从结构张力、位移响应及其时空演化等方面,分析缆线响应特征和机制解释.缆线材料参数见表2,浮体总体分布位置示意见图4(a),通过表3 两种工况考虑浮力对缆线构型和响应的影响.缆线构型如图4(b)所示,可以看出浮力减小会使缆线高度明显减低,case 2 中缆线高度最大降低了约52 m,所以case 1 构型更加安全,因为缆线高度降低会导致下端点附近区域缆线张力下降,使得环境载荷下结构响应幅值增大,从而增加触底风险.

表2 Double-stepped 缆线结构参数Table 2 Structural parameters of double-stepped cable

图4 Double-stepped 缆线的浮体分布和构型Fig.4 Buoyancy distribution and configurations of double-stepped cable

表3 工况条件Table 3 Loading cases

2.1 规则波下的响应

规则波周期为10 s,有效波高3.0 m,作用在缆线上的流体载荷随着水深呈对数衰减,作用水深200 m[31];缆线上、下端点铰支约束,船体运动的影响将在下一节讨论.图5(a) 给出了缆线均方根(RMS)位移响应对比,响应沿展向向下传播过程中幅值呈衰减趋势,最大响应幅值约为5D(D为缆线直径).由于Double-stepped 缆线轴向张力沿长度的非单调变化(图5(b)),响应在传播过程中因为结构张力降低会出现局部幅值增大的现象,在300 m 附近局部峰值约为2.73D;事实上,垂弯曲和浮体分布区域都会出现局部大位移.另外,对比两种工况的响应,响应沿结构传播到底部时case 2 幅值会变得高于case 1,由于F2浮力减小缆线底部响应增大,从而增大了触底风险,也进一步说明上节case 1 构型更加安全.结构张力对响应波传播影响机制将在下一节讨论.

图5 缆线RMS 位移和张力对比Fig.5 Comparisons of RMS displacement and tension

结构位移时空演化如图6 所示,由于分布浮体的存在使得结构参数轴向变化且不连续,复杂构型缆线的响应时空演化变得更加复杂,呈现驻波、行波混合效应.瞬态阶段为行波主导,即响应从顶部激励区域沿结构向下传播,然后(t=90 s)进入稳态响应阶段,这时呈现行波、驻波两种波形.响应传播过程中尽管水平位移逐渐减小,但垂直位移幅值会出现局部峰值(图6(b)、图6(d)).进一步分析响应时空演化(图6(c)、图6(d)),稳态响应阶段仍存在明显的波传播特征(向下传播),而且响应传播过程中波长也是变化的,例如,F2浮体(Y=900 m)附近响应波长变短、幅值增大.

图6 case 1 缆线位移的时空演化结果Fig.6 Spatial-temporal evolutions of cable displacement in case 1

2.2 极端波浪下的响应

考虑极端海况的波浪周期13.3 s、有效波高14.0 m,case 1 缆线响应与常规波浪条件下的对比见图7,总体来看由于载荷增大,极端海况下的响应明显高于常规海况.值得注意的是位移响应波数会减小,1000～1500 m 范围内,波数从常规海况下的4 减少到了极端海况下的3,原因是激励周期增大.另外,可以看出响应经过浮体分布区域后(长度分别为900 m 和200 m)会有大幅的降低,均方根位移幅值从1050 m 位置的13D降低到到900 m 处的4D,从350 m 处的5.2D降低到了123 m 处的3.2D,幅值分别降低了69%和38%.表明浮体分布形成的垂弯区和拱弯区对响应有很好的缓冲效果,即提高了对环境激励的鲁棒性.

图7 两种海况下RMS 位移沿结构展向分布Fig.7 Distributions of RMS displacements for two loading cases

图8 为缆线位移时空演化云图,与上一节图6给出的结果相比,极限波浪载荷下位移沿着长度方向的传播特征,即行波效应更加明显.因为极端海况下缆线响应幅值更大,导致流体阻力增加(式(2)),因此在响应向下传播过程中流体阻力消耗能量更大,响应幅值衰减更加明显,从而呈现明显行波效应.

图8 极端海况缆线位移时空演化云图Fig.8 Spatial-temporal evolutions of cable displacement under extreme wave condition

2.3 顶部船体偏移幅值对响应的影响

实际作业过程中,由于船体和矿车作业的移动会使得缆线上下两端点相对位置变化,一方面会引起缆线空间位置、张力分布以及结构曲率的改变,另一方面也会导致不同的缆线响应和波传播特征.为考察缆线上下两端点相对位置变化对结构响应的影响,计算了船体偏移幅值不同时的缆线响应,观察不同状态下缆线固有频率、张力分布、均方根位移和响应时空演化.不同顶部船体与采矿车的相对运动距离下的缆线构型如图9(a)所示,随着顶部船体远离海底采矿车,缆线垂弯区的高度(y轴位置700 m 处)逐渐增加,而拱弯区域位置逐渐降低,缆线呈现逐渐张紧趋势.图9(b)为不同船体偏移幅值下缆线频率,轴向张力随着船体远离采矿车而逐渐增大,引起结构固有频率增加.

图9 缆线构型和固有频率随船体偏移幅值的变化Fig.9 Cable configurations and natural frequencies with different vessel offsets

常规波浪下典型缆线RMS 位移响应如图10(a)所示,响应传播过程中位移幅值总体呈降低趋势.图10(b)～图10(d)给出了不同区域的响应对比,可见在激励区域附近响应最大,随着船体偏移幅值增加响应逐渐降低(图10(a)和图10(b)),虽然A=150 m时RMS 位移的最大峰值比A=0 大,但对应的谷值较小;在此区域的响应平均值随着船体偏移幅值的增加而减小,表4 给出了响应峰谷值和平均值的统计结果.

图10 缆线均方根位移响应Fig.10 RMS displacement responses

图10 缆线均方根位移响应(续)Fig.10 RMS displacement responses (continued)

表4 1000～1500 m 范围内响应对比 (A/D)Table 4 Comparison of responses in 1000～1500 m region (A/D)

在缆线长度0～1000 m 范围内结构响应随着船体偏移幅值的增大明显减小(图10(c)和图10(d)),值得注意的是,在长度300～1000 m 范围内三种工况下均方根位移响应的波数是不同的,随着船体偏移幅值的增加波数从6 减小到4,这是由缆线上张力的变化引起的.另外,随着缆线张力的降低,响应在向下传播的过程中衰减幅度越小,造成底部产生较大的振幅,在下端点附近响应的幅值差超过了4 倍.图11(a)～图11(c)为船体在三个位置处缆线的位移时空演化,船体距矿车越远缆线1200～1500 m 范围内的响应波长越短,响应向下传播过程中的幅值也越小,这与均方根位移的结果一致.而300～900 m 范围内的三种工况的响应波数不同,在此范围内A=0 m 的工况下有5 个波峰,而其他两种条件下位移响应有4 个波峰.

图11 缆线位移时空云图Fig.11 Spatial-temporal evolutions of cable displacements

3 响应中的波传播和时空演化讨论

从缆线响应的位移及其时空演化结果可知,船体运动到不同位置时响应沿缆线长度传播过程中的波长和波数也发生响应改变.导致此种现象主要有两方面原因,一是船体位置变化导致缆线固有频率改变,二是结构张力本身沿长度变化.为了考察固有频率对波数和波长的影响,这里定义固有频率与激励频率的比值为频率数K,即K=fn/fw,式中fn为结构的第n阶固有频率,fw为激励频率.图12(a)给出了3.3 节三种船体位置下结构频率数,当频率数接近1 则该模态易被激发,可以判断出主要参与模态分别为9,10,11 阶.

图12 变参数缆线的频率数和波长分布Fig.12 Frequency numbers and wavelengths of cable with axiallyvarying tension

另外,结构张力不但影响结构响应幅值,也会改变响应传播波长,这里基于轴向变张力梁动力学方程从理论上解释波传播现象,含流固耦合载荷的变张力梁控制方程为

式中,EI为弯曲刚度,D和m分别为结构的直径和单位长度质量,t为时间,c为阻尼,Fb为浮体上的流体载荷和浮体浮力,V为流体速度,x*为浮体安装位置.假设解的形式为

式中,ω为结构振动频率,θ为相位角.小阻尼条件下阻尼对结构频率和波长影响较小,将上式代入式(20)并忽略阻尼和载荷项可得

假设方程的解为

其中

式中,k(x)和λ(x)分别为局部波数和局部波长.将式(23)和式(24)代入方程(22),得到

考虑张力沿轴向缓慢变化,进一步忽略高阶项,可以得到色散方程

将式(26)代入到式(25)中,通过化简可以得到

两边同乘模态函数a(x),并且考虑张力函数T(x) 沿轴向是线性变化的,式(27)可改写为

流体激励频率是实数,在式(28) 两边同时取实部可得

式(29)可以定性表征波数和波长沿缆线长度的变化.控制方程(20)的通解设为

其中,α1,α2,α3和α4都是常数.进一步 将式(30)代入简支边界条件,可以得到关于 α1,α2,α3和 α4的齐次线性方程组

将上式化简得到

其中,常数A,B,C分别为

由此可以得到简支缆线振型函数Y(x)

因为仅考虑模态形状沿轴向的变化,不考虑模态形状在各点的值,所以 α4可忽略.结合式(29)和式(34)可得波长随张力的变化,图12(b)给出了在长度300～1000 m 范围内缆线波长的变化,可以看出张力越大波长越长,这与图10(c)结果相符,其他区域波长变化受到边界条件的影响会更加复杂.需要注意的是,上述结果和讨论中出现的低张力区域的高位移幅值和短波长响应,会引起结构应力水平的显著提高,从而影响结构安全.

4 结论

考虑深海采矿中环境载荷和顶部船体运动工况,以Double-stepped 缆线为研究对象对其载荷特性进行表征,基于其流固耦合载荷特性和模型表征,建立了考虑非均匀深水管线的复杂结构形式和动态效应的动力学控制方程.基于有限元法和二次开发建立了缆线的动响应计算模型,分析了环境载荷下缆线的响应和传播规律;研究了海面船体位移引起的缆线空间位置、张力变化对结构响应和波传播特征的影响,并给出位移时空演化规律.最后,基于WKB 理论分析了变参数结构的响应幅值和波长演化规律和机理.研究表明:(1)响应沿着缆线长度向下传播过程中并非单调变化,在低张力区域会出现局部峰值,极端海况下缆线位移响应幅值可达16D;(2)由于分布浮体的存在使得结构参数轴向变化且不连续,复杂构型缆线的响应时空演化变得更加复杂,呈现驻波、行波混合效应;而且张力不但会影响位移幅值,还会引起响应传播过程中的波长改变;(3)低张力区域的高位移幅值和短波长效应,引起结构应力水平的提高需要在结构安全分析引起注意.