考虑大变形的大柔性立管涡激振动模型

吴学敏，黄维平

(中国海洋大学 山东省海洋工程重点实验室，山东青岛 266100)

涡激振动(Vortex Induced Vibration，VIV)可导致深海大柔性立管发生疲劳破坏。深水立管在深海复杂环境中，除与海底、平台接触外无任何支撑。研究发现，深海立管在漩涡脱落影响下发生多种振动现象:多模态振动、高阶模态振动、非对称弯曲大变形等［1－6］。正确预测深水立管涡激振动特征对海洋立管设计意义重大。

由于深水立管自身特性及所处复杂深海环境等特点，研究深水立管涡激振动特征用试验方法不易操作。据获得作用于结构的流体力方法不同，预测深水立管涡激振动主要有①半经验方法，该法中的尾流振子模型较常见，典型的为基于非线性van der Pol方程的尾流振子模型，该模型因具有自激、自限性特征，故能捕捉到涡激振动现象，如:锁定现象。但在采用该模型时，需作假设:如当锁定现象发生时，假设漩涡脱落沿立管长度全相关。②预测方法(Computational Fluid Dynamics，CFD)，该方法被公认为较理想的预测细长立管涡激振动方法，但由于三维模拟计算量大等问题，此类模拟较少。而大多数CFD数值模拟方法中，均选取立管若干截面，在各个截面上进行二维流体数值模拟计算。此二维模拟忽略展向尾流的影响。Yamamoto等［7］用CFD方法模拟细长圆柱体涡激振动，结果与试验结果吻合较好;秦延龙等［8］给出较全面的CFD方法研究深水立管涡激振动综述。

目前对海洋立管的动力响应分析模型中，大多采用线性化弹性细长杆模型，忽略了深水立管的几何非线性效应。由于海洋立管涡激振动不仅会出现多模态参与振动现象，且会出现弯曲大变形现象。因此深水立管的几何非线性效应不可忽略。对此本文提出考虑深海立管大变形弯曲振动模型，并用有限元方法分析均匀流中某大长径比深海顶张力式立管的涡激振动特性。

有研究表明，深海立管的顺流向涡激振动振幅与横向振动振幅具有相同数量级［9］。当立管顺流向与横流向固有频率比与顺流向与横流向激发频率比相等时，立管的涡激振动出现双向共振特征，立管呈“8”字形状［10－11］。即顺流向振动对深水立管疲劳损伤同等重要。本文在所提新模型基础上对深海立管进行顺流向、横流向自由度涡激振动大变形振动响应特征分析。

1 数学模型

在笛卡尔坐标系下建立运动方程［12－15］。假设立管不可伸长，即变形前后弧长s不变。图1为深水顶张式立管微元段力学模型。由动量、动量矩守恒可得长ds微段立管单元平衡方程:

图1 微元段Fig.1 Microelement

式中:F为截面内力;q为立管单位长度分布外力;ρ为立管单位长度质量;m为单位长度分布外力矩;ρ=ρtAt+ρiAi为截面力矩，为弯矩与扭矩之和，表示为:

式中:B=EI为截面弯曲刚度;H为截面扭矩。

由式(3)得:

将式(4)代入式(2)并令H=0、m=0，即忽略单元扭矩、外力矩，得:

式中:λ为Lagrange乘子;T(s，t)=r';F'为立管张力;κ2= －r'r‴，κ(s，t)为局部曲率。

将式(5)代入式(1)得大挠度细长梁运动方程为:

将式(7)写成顶张力式立管大变形弯曲振动微分方程为:

2 数值计算

为验证模型式(8)、式(9)，本文开发出深水立管涡激振动分析程序。由于深水立管涡激振动的非线性特征，程序采用Newmark－β隐式积分法求解。Newmark－β法为直接积分法，β=0.25时无条件稳定。

设大变形深海立管顶点为坐标原点，u轴为垂直方向坐标轴，v轴为水平方向坐标轴。将弯曲振动方程式(8)、(9)离散后得深海立管大变形弯曲振动有限元方程为:

式中:［Mv］为立管横向质量矩阵:

［Mu］为立管竖向惯性系数矩阵:

［Kv］为立管横向刚度矩阵:

［Ku］为立管竖向刚度矩阵:

［Cv］为立管横向阻尼矩阵:

［Cv］= ［Mv］+ ［Kv］

［Cu］为立管竖向阻尼矩阵:

［Cu］= α［Mu］+ β［Ku］

其中:α，β 为瑞雷祖尼系数，单位 s－1，s，两者取值为［16］:取体系基频，取体系第三阶自然频率;{Fu}为立管竖向荷载向量;{Fv}为立管横向荷载向量;{}为立管横向加速度向量;{}为立管横向速度向量;{v}为立管横向位移向量;{}为立管竖向加速度向量;}为立管竖向速度向量;{u}为立管竖向位移向量;［N］为立管单元插值函数矩阵;l为单元长度;n为单元数量。

引入大变形，系统质量矩阵、阻尼矩阵、刚度矩阵与节点曲率相关。因此在每个时间步，有限元方程的质量矩阵、刚度矩阵、阻尼矩阵均未变化的。为提高计算精度，对每个时间步进行迭代，时间步长dt=0.01 s。

3 算例分析

用以上程序对1 500 m水深顶张式立管(TTR)进行VIV大变形涡激振动响应特征分析。该TTR为双套管(双屏)生产立管，外套管外径324 mm，内套管外径222 mm，油管外径114 mm。基于弯曲刚度等效原则将立管等效为单层管分析模型，等效后内径为292 mm，顶张力系数取1.4。计算模型边界条件一端固定、一端简支。

3.1 计算结果验证

由于Shear7采用基于模态叠加频域分析方法，计算结果为响应的均方根值。因此，本文将计算结果统计分析，给出响应均方根值。图2、图3为固定－简支边界条件的均方根位移与均方根应力沿立管轴向分布曲线，虚线为本文分析模型计算结果，实线为商用软件Shear7分析结果。可以看出，本文所提分析模型与Shear7模态叠加法计算结果吻合较好。

图2 均方根位移曲线(固定－简支)Fig.2 RMS Displacement along TTR

图3 均方根应力曲线(固定－简支)Fig.3 RMS Stress along TTR

3.2 动力特性分析

图4为 Re=1.2×105时水深 Z=300 m、750 m、1350 m处三水平横截面立管节点横向振动时程曲线。并经傅里叶变换得图5三水平截面处横向振动响应频谱。由图5(a)看出，立管在Z=300 m处有六个峰值频率:fv=0.029 9 Hz，0.069 9 Hz，0.119 9 Hz，0.159 9 Hz，0.199 9 Hz，0.259 9 Hz;在 Z=750 m 处有五个峰值频率:fv=0.029 9 Hz，0.069 9 Hz，0.119 9 Hz，0.159 9 Hz，0.259 9 Hz;在 Z=1 350 m 处也有五个峰值频率:fv=0.029 9 Hz，0.069 9 Hz，0.119 9 Hz，0.159 9 Hz，0.259 9 Hz。可见，立管水平位置不同，其控制模态亦不同，均发生多模态参与现象，参与的模态数目也不相同，且参与模态间频率间隔较小。研究发现此时所对应的频率峰值均相同。

图6为Re=1.2×105时水深 Z=300 m、Z=750 m、1 350 m处三水平横截面立管表面某点顺流向振动响应频谱图。由图6看出:立管在Z=300 m处有六个峰值频率:fv=0.029 9 Hz，0.069 9 Hz，0.119 9 Hz，0.159 9 Hz，0.199 9 Hz，0.259 9 Hz;Z=750 m 处有四个峰值频率:fv=0.029 9 Hz，0.119 9 Hz，0.159 9 Hz，0.259 9 Hz;Z=1 350 m处有五个峰值频率:fv=0.029 9 Hz，0.069 9 Hz，0.119 9 Hz，0.159 9 Hz，0.259 9 Hz。可见，同样在立管三个不同水平位置处，立管顺流向振动的控制模态也不同，但对应的频率峰值相同。均发生多模态参与振动现象。

图4 Re=1.2×105时立管横向振动时程曲线Fig.4 Time-history curve of cross－ flow vibration when Re=1.2 ×105

图5 Re=1.2×105时立管横向振动频谱图Fig.5 Spectrogram of cross-flow vibration，when Re=1.2 ×105

图6 Re=1.2×105时立管顺流向振动频谱图Fig.6 Spectrogram of in-line vibration，when Re=1.2 ×105

图7 Re=6 300时立管横向振动频谱图Fig.7 Spectrogram of cross-flow vibration，when Re=6 300

图8 当Re=6 300时立管顺流向振动频谱图Fig.8 Spectrogram of in-line vibration，when Re=6 300

在Re=1.2×105时，立管的顺流向、横向振动均发生多阶模态参与振动现象，频率峰值基本相同。

图7、图8分别为Re=6 300时水深Z=300 m、750 m、1350 m处三个水平横截面节点横向、顺流向涡激振动响应频谱图。由二图可见，低雷诺数时立管涡激振动以低阶单模态为主，漩涡脱落并未激发多个模态参与振动，此时，横向频率峰值为fv=0.024 9 Hz，顺流向频率峰值为fv=0.024 9 Hz，两向频率峰值相等。

图9为两种雷诺数时立管振动横向振幅达最大的变形图。Re=6 300时，立管横向最大振幅为0.195 5 m，发生在Z=825 m处;顺流向最大振幅0.005 8 m，发生在Z=795 m处，且立管振动以一阶单模态振动为主。当Re=1.2×105时，立管横向最大振幅为0.503 3 m，发生在Z=1 395 m处;顺流向最大振幅发生在Z=810 m处，且漩涡脱落激起高阶模态参与振动。立管涡激振动变形并非完全对称，且有大变形现象。随雷诺数的增大，立管高阶模态被激发。

图9 横向振幅最大时立管振动变形图Fig.9 Deformation figures of riser，when the max CF amplitude

由数值模拟结果知，雷诺数不同，立管不各位置参与振动的模态也不同。雷诺数较低时，立管涡激振动以低阶单模态振动为主;随雷诺数的增大，立管涡激振动出现多模态参与振动现象，甚至立管高阶模态也被激发。涡激振动振幅也随雷诺数的增大而增大。雷诺数相同时，立管发生涡激振动的不同时刻，其最大振幅发生的位置亦不同，且立管涡激振动均发生非对称大变形现象。

表1 Re=6 300时立管涡激振动响应振幅Tab.1 Vortex-induced vibration amplitude of riser when Re=6 300

本文在Re=6 300、Re=1.2×105两种雷诺数下，将考虑大变形与不考虑大变形立管涡激振动的最大振幅进行比较，结果分别见表1、表2。由表1、表2数据知，大变形在低雷诺数时对振幅影响较小，在高模态数情况下对振幅影响增大。

表2 Re=1.2×105时立管涡激振动响应振幅Tab.2 Vortex-induced vibration amplitude of riser when Re=1.2 ×105

4 结论

本文在已有研究基础上，考虑深水立管因其大长径比特点，在深海易发生大变形现象，提出深水顶张力式立管大变形涡激振动分析模型，研究立管在不同雷诺数下的涡激振动特性。结论如下:

(1)深海立管多模态响应与非对称弯曲大变形等涡激振动现象在低雷诺数情况下未能体现;却体现在高模态数情况下。

(2)考虑大变形的几何非线性分析模型可用于深水、极深水立管系统涡激振动响应分析预测。该模型考虑深水立管在复杂深海环境下易发生大变形现象，将大变形几何非线性特征引入深水立管涡激振动分析中，为寻求更准确模型分析深水立管方法提供参考。

［1］Modarres-Sadeghi Y，Chasparis F，Triantafyllou M S，et al.Chaotic response is a generic feature of vortex-induced vibrations of flexible risers［J］. Journal of Sound and Vibrations，2011，330(11):2565 －2579.

［2］Wu X D，Ge F，Hong Y S.A view of recent studies on vortex-induced vibrations of long slender cylinders［J］.Journal of Fluid and Structures，2011，28:292 －308.

［3］Williamson C H K，Govardhan R.A brief review of recent results in vortex-induced vibrations［J］.Journal of Wind Engineering and Industrial Aerodynamics，2008，6(6－7):714－735.

［4］Willden R H J，Graham J M R.Multi-modal vortex-induced vibrations if a vertical riser pipe subject to a uniform current profile［J］.European Journal of Mechanics B/Fluids，2004，23(1):209－218.

［5］So R M C，Liu Y，Cui Z X，et a1.Three-dimensional wake effects on flow-induced forces［J］.Journal of Fluids and Structures，2005，20(2):373 －402.

［6］Trim A D，Braaten H， LieH， eta1. Experimental investigation of vortex-induced vibration of long marine risers［J］.Journal of Fluids and Structures，2005，21(3):335 －361.

［7］Yamamoto C T，Meneghin J R.Numerical simulations of vortex-induced vibration on flexible cylinders［J］.Journal of Fluids and Structures，2004，19(4):467－489.

［8］秦延龙，王世澎.海洋立管涡激振动计算方法进展［J］.中国海洋平台，2008，23(4):14 －17.QIN Yan-long， WANG Shi-peng.Computational methods progress of the vortex-induced vibration of ocean risers［J］.China Offsh Ore Platform，2008，23(4):14－17.

［9］Vandiver J K，Jaiswal V，Jhingran V.Insights on vortexinduced，traveling waves on long risers［J］.Journal of Fluids and Structures，2009，25(4):641 －653.

［10］Dahl J M，Hover F S，Triantafillou M S，et al.Dual resonance in vortex-induced vibrations at subcritical and supercriticalreynolds numbers［J］. JournalofFluid Mechanics，2010，643:395 －424.

［11］Dahl J M，Hover F S，Triantafillou M S.High harmonic forces and predicted vibrations from forced in-line and crossflow cylinder motions［C］.In:International Offshore and Polar Engineering Conference Proceedings，2008:481 －488.

［12］Nordgren R P.On computation of the motion of elastic rods［J］.ASME Journal of Applied Mechanics，1974，41(3):777－780.

［13］Garrett D L.Dynamic analysis of slender rods［J］.Journal of Energy Resources Technology，Transaction of ASME，1982，104(4):302－307.

［14］Pauling J R，Webster W C.A consistent large-amplitude analysis of the coupled response of a TLP and tendon system［J］. Proceedings ofOffshore Mechanics and Arctic Engineering(OMAE)，1986，3:126－133.

［15］Ma W，Webster W C.An analytical approach to cable dynamics:theory and user manual［R］.Sea Grant Project R/OE －26，1994.

［16］克拉夫R，彭津J，著.王光远，译.结构动力学(第二版)［M］.北京:高等教育出版社，2007:183－185.