0 引　言

大型海洋立管是连接海洋平台和水下生产系统的关键结构，整个立管系统中间没有任何支撑，对其结构动力学性能要求更高，技术含量也更高^[1]。当旋涡脱落频率与立管固有频率接近时，就会发生涡激共振现象，使得振幅在一定的流速范围显著加大，从而导致立管疲劳破坏，产生巨大的经济损失^{[2 − 4]}。目前，工程上主要通过安装一些扰流装置达到破坏旋涡结构或者改变涡脱模式来减弱立管的涡激振动^{[5 − 7]}。然而，安装了扰流装置后，立管受到的阻力将大大增加，此外扰流装置还有可能引起其它形式的振动^[8]。尽管增加扰流装置有着一些缺点，但是改变立管截面的形状能够非常有效的抑制涡激振动的产生。目前工程上依然是主要依靠扰流装置来抑制涡激振动。螺旋列板涡激振动抑制装置能够抑制各个来流方向的涡激振动，因此在海洋立管的涡激振动中使用最为广泛。

对于螺旋列板涡激振动抑制装置的设计，多数依赖于计算流体力学（CFD）方法和水池拖拽试验方法^{[9 − 11]}。在流体的作用下，流场中的结构会发生弹性变形，这种结构的弹性变形又对流场分布产生影响，从而使流体和结构形成一个相互联系、相互作用的复杂系统。因此，考虑涡激振动的立管的升力系数也会与静止情况下有较大差异。通过对比静态立管的升力系数难以准确证实螺旋列板的涡激振动抑制效果。目前，有学者通过二维弹性支撑柱体的模型来计算螺旋列板的抑制效果。但这种方法忽略了螺旋列板引起流体的三维效应，无法考虑螺距等参数的影响。同时这类方法并未给出从三维模型简化为二维模型的合理解释。随着计算机科学的高速发展，采用全三维立管的计算流体力学/有限元（CFD/FEM）双向流固耦合模拟成为可能，可以进行加装螺旋列板后全尺寸立管的涡激振动计算。但该方法对计算资源要求非常高，且计算非常耗时，并不适用于工程上螺旋列板结构的优化设计。在海洋中，沿着立管轴向的流速分布不同。在立管的涡激共振中，往往是低流速主导大振幅低频率振动，高流速激发立管小振幅高频率振动。螺旋列板抑制涡激振动的主要原理是破坏旋涡结构，大幅度降低涡脱频率，避开结构固有频率，从而避免发生涡激共振。因此，设计出的螺旋列板能够抑制低流速下的大振幅，那对于高流速下的涡激振动的振幅抑制也必然有效。另外，当涡激共振发生的情况下，流速激发的是立管某一阶的模态。立管的振动响应可以通过模态叠加的方法计算，对于涡激共振情况，除了所激发的立管模态，其它阶的能量非常小，可以忽略。因此，涡激共振所在的那一阶模态可以直接作为系统建模的动力学参数。

在亚临界雷诺数区间内，斯托怒哈尔数${S_t} = 0.2$，不同流速对应的旋涡脱落频率固定不变。本文首先计算全三维立管的振动模态，考察立管来流向最小速度所能激发的立管模态。以该模态参数为弹性支撑螺旋列板抑制装置结构的刚度参数，选取2个螺距长度的立管（加装螺旋列板）为研究对象，建立系统的流固耦合动力学模型，预测加装螺旋列板后的立管的横流向振幅和拖拽力系数。

1 螺旋列板结构设计及参数

图1为本文设计的一个单元的螺旋列板涡激振动抑制装置。其结构基本参数如表1所示。该装置主要由卡条、鳍、凹槽等结构组成。单个螺距长度的螺旋列板涡激振动抑制装置由6个单元的该装置组成。可以通过使用气动工具将绑带绑在凹槽部位，从而逐个的将该装置固紧在立管上，沿着立管进行安装。通过装置上的鳍破坏旋涡的结构，大幅度降低涡脱频率，避开立管结构固有频率，从而达到抑制涡激振动的目的。鳍存在合理的斜度以及装置上的凹槽区，降低了拖拽力的增加。

2 抑制装置涡激振动动力学建模

基于CFD商业软件Ansys Fluent和结构动力学基本理论以及嵌套网格技术，建立弹性支撑刚性立管段（加装螺旋列板抑制装置）的流固耦合动力学模型。流场的求解基于非定常不可压缩流体RANS方程为：

$ \frac{{\partial {{\bar u}_i}}}{{\partial {x_i}}} = 0，$

(1)

$ \frac{{\partial {\rho _f}{{\bar u}_i}}}{{\partial t}} + \frac{{\partial {\rho _f}\overline {{u_i}{u_j}} }}{{\partial {x_j}}} = - \frac{{\partial \bar p}}{{\partial {x_i}}} + \mu {\nabla ^2}{\bar u_i} - \frac{{\partial {\rho _f}\overline {{{u'}_i}{{u'}_j}} }}{{\partial {x_j}}}，$

(2)

$ - \rho \overline {{{u'}_i}{{u'}_j}} = {\mu _t}\left( {\frac{{\partial {u_i}}}{{\partial {x_j}}} + \frac{{\partial {u_j}}}{{\partial {x_i}}}} \right) - \frac{2}{3}\rho {k_t}{\delta _{ij}}。$

(3)

式中：$ {\rho _f} $为不可压缩流体的密度；${u_i}$为$i$方向上的瞬时速度分量；${u'_i}$为$i$方向上速度脉动量；${\bar u_i}$为速度的时间平均值；$ - \rho \overline {{{u'}_i}{{u'}_j}} $为雷诺应力张量；${x_i}$、$t$、$p$、$ \mu $ 分别为笛卡尔坐标系、时间、压力、运动粘度；$ {\mu _t} $为湍流黏度，下标“t”为湍流；${k_t}$为湍动能；$ {\delta _{ij}} $为“Kronecker delta”符号，当$i = j$时，${\delta _{ij}} = 1$，当$i \ne j$时，${\delta _{ij}} = 0$。湍流模型采用SST k–ω模型进行流场仿真，计算得到三维立管段（加装螺旋列板）表面的压力分布，进而可以得到三维立管段（加装螺旋列板）升力和阻力系数：

$ {F_D} = \frac{1}{2}{C_D}{\rho _f}{U^2}{{D}}，$

(4)

$ {F_L} = \frac{1}{2}{C_L}{\rho _f}{U^2}{{D}}，$

(5)

根据结构动力学理论，两自由度（2-DOF）弹性支撑的柱体（可加装螺旋列板抑制装置）运动的控制方程可以写为：

$ \begin{split} & m\ddot x + c\dot x + kx = {F_D}(t)，\\ & m\ddot y + c\dot y + ky = {F_L}(t)。\end{split} $

(6)

式中：$m$为圆柱体的质量；$c$为结构阻尼系数；$k$为结构刚度系数。

无量纲形式的式(6)又可以写为：

$ \begin{split} & \ddot x + 2\varsigma {\omega _0}\dot x + \omega _0^2x = {F_D}(t)/m，\\ & \ddot y + 2\varsigma {\omega _0}\dot y + \omega _0^2y = {F_L}(t)/m。\end{split} $

(7)

式中：柱体固有频率${\omega _0} = \sqrt {k/m} $；阻尼比$ \varsigma = {c /{2\sqrt {km} }} $。

综合考虑计算资源的情况下，流场域的尺寸大小取尾迹区域30D，立管段（加装螺旋列板抑制装置）前端和上下距离柱体都是10D。流场入口边界条件为速度入口，出口为压力出口，上下壁面为滑移壁面，柱体表面即动边界为无滑移壁面。

通过动网格技术来实现流场中刚性立管段（加装螺旋列板）边界的运动。基于嵌套网格技术的动网格方法适用于刚性边界运动问题。流场域网格划分采用嵌套网格，如图2所示，靠近立管表面为边界层网格（${Y^ + } < 1$）。采用基于嵌套网格技术的动网格方法，可以无需担心网格畸变以及负网格导致求解失败等问题。如图2(e)所示，嵌套网格的嵌套过程是求解器识别嵌套网格边界，对被组分网格遮蔽的背景网格部分进行“挖洞”，然后对嵌套区域边界单元进行插值，将背景区域的边界单元变量信息插值到嵌套区域的边界单元，最后进行流场计算。如图2(c)和图2(d)所示，对于加装螺旋列板后的立管划分结构化网格非常复杂，因此直接采用非结构化网格，为了减少计算量，将非结构化网格进一步转化为多面体网格。整个流场的计算网格如图图2(f)所示。对于流场的数值计算，时间项采用全隐式积分方法，对流项则采用二阶迎风离散格式。控制方程中速度分量与压力的耦合则采用COUPLED算法进行处理。

初始条件为$x(0) = \dot x(0) = y(0) = \dot y(0) = 0$。时间步长按照公式$t = D/{\text{(}}U\times{\text{(}}10\sim 20{\text{))}}$来计算，其中，D为立管（加装螺旋列板）的水力直径，U为来流速度。

3 网格无关性验证

如表2和表3所示，对2倍螺距长的裸立管和加装螺旋列板后的立管分别画出3套网格进行网格无关性验证。计算裸立管（不同流速）和加装螺旋列板（流速为0.5 m/s）后立管的升力系数、阻力系数，进行网格无关性验证。从表中数据可知，裸立管流场计算网格为663万时，水动力系数基本不变；加装螺旋列板后的立管流场网格数量为1286万时，水动力学系数基本不变。综合考虑计算资源，裸立管和加装螺旋列板后的流场网格都选取第二套网格。

4 计算结果及分析

仿真计算了来流速度为0.25、0.5、1、1.5、2 m/s下裸立管和加装螺旋列板抑制装置后的立管横流向振幅、拖拽力系数、来流向立管弯曲位移等结果。从图3可知，当来流速度为0.25时，横流向振幅最大，此时发生了涡激共振现象，拖拽力系数也同时被放大很多，达到了1.75。当流速大于0.5之后，流速处于超临界区间，该区间内旋涡脱落没有明确的频率，因此不会发生涡激共振。

从图3(a)可知，在流速为1、1.5、2 m/s时，裸立管的横流向振幅非常小，接近于0。对比加装了螺旋列板抑制装置之后的立管横流向振幅，可以看出，所有流速下，横向振幅都非常小，横流向振幅抑制效率超过98%。在流速为1.5、2 m/s时，相比裸立管的横流向振幅略微变大，这是由于超临界区间，振动形式以随机振动为主，加了螺旋列板之后在横流向反而引起了微弱的振动，但振幅都非常小，可以忽略其对立管的影响。

从图3(b)可知，加装本文设计的螺旋列板之后，拖拽力系数都小于1.6，满足工程要求。随着来流速度的增大，阻力系数略微变小，整体稳定在1.5～1.6之间。从图3(c)可知，在流体力的作用下立管在来流向发生一定的位移，代表整个1500 m长立管在该点处发生的弯曲位移。从图中可知，流速越大，立管弯曲位移越大。裸立管段在最大流速下，弯曲位移小于3 m，而加装了螺旋列板之后，立管段的来流向最大弯曲位移小于12 m。因此，对于加装螺旋列板引起的弯曲位移也会导致立管的承受巨大的弯曲应力，对立管结构产生不利的影响。在工程设计中，也需要合理考虑加装螺旋列板后引起立管大幅度弯曲位移带来的影响。

图4和图5为0.25 m/s时，裸立管段的涡激振动仿真结果、涡核云图。对应地，图6和图7为0.25 m/s时，加装螺旋列板之后的立管段涡激振动响应结果以及涡核云图。图8为以上2种工况对应的频谱分析图。图8(a)～图8(d)为裸立管段在0.25 m/s来流速度下的升力系数响应、阻力系数响应、横流向振动响应、来流向振动响应。可知，各响应的幅值都很大。对图4(a)中的升力系数做频谱分析，计算结果如图8(a)所示，说明裸立管在该工况下的旋涡脱落频率为0.075 Hz。该裸立管段对应的湿模态频率的计算公式为$ k = (m + {m_{{\text{water}}}})\omega _{{\text{water}}}^2 $，需要进一步考虑立管段排开的水质量$ {m_{{\text{water}}}} $。换算出的立管段湿模态频率$ {\omega _{{\text{water}}}} $=0.081 Hz，与此刻的旋涡脱落频率接近，因此在0.25 m/s时刻，裸立管段发生的是涡激共振现象。从图4(e)可知，运动轨迹成“8”字形，这是涡激共振典型的特性。从图4(f)可知，升力系数与振幅的响应曲线是相同相位的，因此升力系数相比静止刚性立管的升力系数放大了近10倍。在涡激共振发生的情况下，阻力系数相比静止立管的工况也被放大了。图5为立管段从0 s到120 s的涡脱模式和立管位置变化。可以清晰的看到立管在来流向有一定的位置偏移，在横流向有较大的振幅，涡脱模式以“2P”形式为主。

对比于图6中加装螺旋列板的情况，可知本文设计的螺旋列板能有效减小升力系数响应的幅值，其中阻力系数响应的振幅趋近于0，横流向的振幅为0.00145 m，来流向振幅几乎为0。计算结果显示，螺旋列板的涡激振动抑制效率超过98%，螺旋列板对立管涡激振动有很好的抑制作用，但是会一定程度的增加立管的弯曲应力。从图7可知，加装螺旋列板后立管段尾部的涡核被打碎，没有了明确的涡脱频率。从图中还可知立管相对于0 s时刻位置除了来流向有一定的位移，横流向几乎没有振动。通过对图6(a)中的升力系数响应进行频谱分析，得到图8(b)中的频谱分析结果。加装螺旋列板之后，改变了原有立管有序的涡脱频率，涡激力也失去了尖峰特征，能量非常小。计算旋涡脱落频率为0.020 Hz，距离立管段的固有湿模态频率0.081 Hz较远，因此不会发生涡激振动现象。

5 结　语

本文基于有限元模态计算、计算流体力学（CFD）方法、结构动力学理论以及嵌套网格技术，提出一种适用于螺旋列板涡激振动抑制装置的动力学设计方法。相比传统水动力计算方法更可靠，相比全尺寸涡激振动抑制装置流固耦合计算更节省计算资源以及计算时间。根据仿真结果，加装螺旋列板后涡核被打碎，涡脱频率远离结构固有湿模态频率，且涡激力的能量非常小，因此不会产生涡激振动。本文算例中螺旋列板横向振幅抑制效率达到98%以上，拖拽力系数小于1.6，可以满足工程需求。但加装螺旋列板又会增加整根立管的弯曲位移，增加立管的弯曲应力，在工程设计中也需要考虑弯曲应力带来的影响。