0 引 言

立管在涡激振动的作用下发生大变形，导致各立管间距变小，在立管遭受涡激振动疲劳损伤的同时大大提高了立管发生碰撞断裂的机率。在流场干涉效应下立管群的运动响应特性与独立立管有所区别^[1～7]，如即使未发生锁振立管也可能出现大振幅的驰振现象。而目前为止，学者们对于涡激振动抑制措施的研究基本都以单一结构物为设计对象^[8～15]，关于多根柔性圆柱体涡激振动抑制措施的研究依然很少^[16]。因此，管群相互干涉下的涡激振动特点及抑制措施对管群之间相互作用的影响是今后研究工作的一个重要发展方向。本文采用商业软件 Ansys15.0 对均匀流下 2 根串联排列、轴间距为 5D 的相同深海立管模型在流场干涉效应下的涡激振动特性进行流固耦合研究，分析侧板对系统中不同位置处立管涡激振动的抑制效果。

1 数值计算方法 1.1 流体基本控制方程

当外力在不考虑能量损耗的情况下，牛顿流体的运动可由质量守恒与动量守恒加以描述。

1）连续方程

连续方程是质量守恒定律在流体力学中的表现形式，对于不可压缩流体该方程可以表述为单位时间内流入该微元的净质量等于单位时间流出该微元的净质量。

$\frac{\partial \rho }{\partial t}+\frac{\partial }{\partial {{x}_{j}}}(\rho {{u}_{j}})=0$

(1)

2）动量方程

动量方程是动量守恒定律在流体力学中的表现形式，其本质是满足牛顿第二定律，该方程可以表述为对于一给定的流体微元，其具有的动量对时间的导数等于作用在该微元上的各种力之和。

$\frac{\partial }{\partial t}(\rho ui)+\frac{\partial }{\partial xj}(\rho uiuj)=-\frac{\partial p}{\partial xi}+\frac{\partial }{\partial xj}(2\mu Sij)$

(2)

式中 x_i，x_j，u_i，u_j 分别为笛卡尔坐标系下的位置和速度分量；ρ，t，p 分别为密度、时间与压力；μ 为分子粘性；S_ij 为应变率张量：

${{S}_{ij}}=\frac{1}{2}(\frac{\partial ui}{\partial xj}+\frac{\partial uj}{\partial xi})$

(3)

1.2 结构动力学控制方程

在流场作用下，顶张力柔性立管的结构运动微分方程为：

$\begin{array}{*{35}{l}} \frac{{{\partial }^{2}}}{\partial {{z}^{2}}}\left( EI\frac{{{\partial }^{2}}v}{\partial {{z}^{2}}} \right)-\frac{\partial }{\partial z}\left( T\frac{\partial v}{\partial z} \right)+m\frac{{{\partial }^{2}}v}{\partial {{t}^{2}}}+c\frac{\partial v}{\partial t}={{F}_{(t)}}, \\ m\frac{{{\partial }^{2}}u}{\partial {{t}^{2}}}-EA\frac{{{\partial }^{2}}u}{\partial {{z}^{2}}}=0 \\ \end{array}$

(4)

式中，EI，EA，m 分别为抗弯刚度、抗拉刚度和单位长度的质量；T和F_（t） 为张力和流体动力载荷；v和u 为横向位移和轴向位移。

对控制方程采用基于三维单元的有限元方法进行离散得到有限元方程：

$M\ddot{d}+C\dot{d}+Kd={{F}_{(t)}}$

(5)

式中：d，$\dot d$，$\ddot d$ 分别为单元节点处位移、速度和加速度矢量；M，C 和 K 分别为质量矩阵、阻尼矩阵和刚度矩阵。

1.3 计算模型

Tezduyar 和 Shih^[17]提出，为避免选取矩形流域时下游壁面边界及侧壁面边界对圆柱周围流场产生影响，下游壁面边界距圆柱中心距离需大于 14.5 D，侧壁距圆柱中心需大于 8 D，三维圆管展向长度应大于 πD，其中 D 为圆柱直径。本文中选取计算流域如图 1 所示，上游立管中心距流域入口为 10 D，下游 30 D，距两侧壁面均为 10 D，展向长度取为立管长度。对于串列双立管，立管圆心间距取 5 D。入口条件为速度入口，U = 0.42 m/s，来流方向为 X 轴正向；出口条件为压力出口，忽略纵向速度梯度；上下壁面及圆柱表面满足无滑移固壁条件。柔性立管模型参考 Lehn（2003）模型试验^[18]，参数如表 1 所示。

1.4 流固耦合数值模拟过程

Ansys Workbench 平台流固耦合计算各模块连接如图 2 所示。首先利用 Design Modeler 模块建立流域及立管模型；其次在 Fluid Flow（Fluent）模块中划分流域网格并设置流场参数，网格划分如图 3 所示。流场采用标准 k-ω 模型，压力和速度耦合采用稳定性更高的 SIMPLEC 方法；压力项离散格式采用 2 阶精度，对动量方程的离散采用 2 阶迎风格式；由于动网格的使用时间离散方式只能采用 1 阶隐式；动网格技术采用网格光顺中的扩散光顺方法；计算中采用欠松弛技术以改善计算收敛条件，残差大小通常控制在 1.0E-05。然后在 Mechanical 模块中划分结构网格，见图 4，并将立管表面设为 Fluid Solid Interface，立管底部为简支端，上部顶张力为 817 N 且沿轴向立管可自由运动，XY 方向固定。最后利用 System Coupling 模块进行流固耦合及数据传递，3 个模块中的时间步长及计算总时长必须保持一致。

2 结果分析 2.1 结构动力学特性分析

图 5 为立管横向无量纲振幅均方根曲线。首先，由本文光滑单立管的数值结果与文献^[18]实验结果对比，2 组数据立管中下部模拟较好，在 Z/L = 0～0.3 及 0.9～1 范围内数值结果较实验值偏小，中上部偏大，由于试验与数值模拟的差别以及数值计算中网格划分方式、参数设置不同等原因，导致一定的误差，但误差尚可以接受；且 2 组数据都反映了立管在此来流情况下能够发生 2 阶振动，下部变形大于上部，很好地验证了本文流固耦合数值模拟的可行性。对于串列双立管，上半部分变形较单立管有所增加，下半部分有所减小；下游立管无因次均方根振幅略大于上游立管。崔洋洋^[19]对相同立管模型、间距为 4D 的光滑双立管进行数值计算，与单立管对比得到了相似的结论。

图 6 和图 7 分别为光滑双立管和附加侧板双立管横向（Cross Flow，CF）及顺向（In Line，IL）位移时间图。对于光滑立管，下游立管在上游立管的尾流作用下产生振动，并且比上游立管更快达到稳定振动状态；在横流向，同一时刻对应立管不同位置处的振动幅值有正有负，且正值与负值各自对应点关于 Z = 0.5 L 对称，说明立管最终锁定在二阶振动模态，此外，对于不同立管，同一节点振动方向相反，说明上下游立管在横向为反相振动；在顺流向，同一极值时刻对应各点的幅值均为最大值或最小值，说明在顺流向上下游立管均为 1 阶振动模态，且为同相振动。对于附加侧板的双立管，在横流向，同一极值时刻，同一立管各点对应幅值均为最大值或最小值，说明了双立管在横流向为 1 阶振动，同样两立管为反相振动，最大振幅点在 Z/L = 0.5 节点处；在顺流向，由图可以看出，上游立管在来流作用下发生了大变形，到达平衡位置后做缓慢振动，其振幅几乎为 0；下游立管在尾流作用下同样发生了大变形，然后在平衡位置振动，其振幅大于上游立管。

对各振动曲线进行傅里叶变换求得其振动主频，并和横向变形最大值、顺向位移最大值列于表 2，通过与同等条件下的单根柔性立管进行对比，分析 2 种情况下立管振动特性的不同。结果发现：对于光滑双立管，上游立管和下游立管的振动频率一致，说明在尾流干涉作用下双立管达到了同步振动，且顺向振动频率与横向振动约满足 2 倍关系；此外，顺向位移超过了横向最大变形，说明了柔性立管由于弯曲强度削弱导致顺流向更容易发生大变形。附加侧板后，双立管横向振动频率均由原来的 4.808 Hz 降至 2.345 Hz，使振动约减缓了 1 倍；上游立管在侧板的作用下横向振幅降幅显著，振幅最大值由 0.87 D 减小至 0.1 D；下游立管振幅反而略有增加，值为 0.89 D；由于侧板的存在，阻力增大，使上游立管顺向最大位移值约增大 1 倍，下游立管的顺向位移响应增幅很小。

图 8 和图 9 分别为光滑双立管和附加侧板立管横向振型图。光滑双立管在横流向为反相振动，锁定在 2 阶模态。附加侧板后锁定在 1 阶振动模态，且上游立管在横流向做微幅振荡，下游立管做大幅振动。

图 10 为光滑双立管和附加侧板立管顺向振型图。2 种工况下，立管在顺流向均为一阶振动，且上游立管顺向位移略大于下游立管，光滑双立管在运动过程中最小间距为 4.88 D，略低于静止时刻立管间距，附加侧板后，由于上游立管顺向位移变大，导致双立管间距进一步减小，约为 3.95 D，说明侧板的存在虽然降低了立管双向振动频率及上游立管横向振动幅值，但增加了立管顺流向碰撞的机率。图 11 为立管总变形图。

2.2 尾流场分析

对于光滑双立管，由涡量图 12 可以看出，当立管间距为 5 D 时，两立管之间有完整的旋涡脱落。在 Z/L = 0.1，0.3 和 Z/L = 0.5，0.7 截面处，双立管上下部分旋涡脱落方向相反，导致立管振动方向相反，立管间涡泄为 P + S 模式，下游立管后缘处为 2S 模式，且为双排涡，Z/L = 0.3 和 0.7 两个截面处旋涡横向间距较其他截面偏大，约为一倍立管直径的距离，与立管横向最大变形位置相对应。Z/L = 0.5 截面处两立管间旋涡为单排涡，下游立管尾流场由于立管中间位置变形较小呈现 2 条横向间距很小、几乎平行的尾涡。且由 Z/L = 0.3 和 Z/L = 0.7 截面图可以清晰地看出，下游立管近壁面处的尾涡脱落方向恰好与上游立管尾缘处旋涡脱落方向相反，说明同一高度处上下游立管横向振动方向相反。附加螺旋侧板后，两立管之间的流动情况较为杂乱，来流绕过上游立管壁面从侧板上分离。但下游立管在上游立管的尾流作用下其尾流场依然有相当量的旋涡脱落，导致下游立管发生了大振幅振动，对于不同截面，下游立管均偏向中心线的一侧，说明立管横流向振动为一阶模态。由附加螺旋侧板的立管三维涡量图 13 可以看出，下游立管的尾涡在轴向上呈螺旋状脱落。另外，通过不同截面处尾涡的差异及三维涡量图直观地表现了高速来流下立管尾流场旋涡发放的三维特性。

3 结 语

本文对 Re = 7 800 的均匀来流下长径比为 482 的串列双立管进行了数值模拟，研究其涡激振动特性及侧板对双立管涡激振动的抑制效果，主要结论如下：

1）当立管间距为 5 D 时，光滑双立管之间有完整的尾涡脱落，出现了多种脱落模式。上下游立管达到了同步振动，且顺向振动频率与横向振动均满足约 2 倍的关系；下游立管横向变形略大于上游立管；但顺向位移上游立管略大于下游立管，使立管顺向间距有所减小。

2）附加侧板后双立管振动频率约减缓了一倍，横向振动模态亦有所降低，上游立管横向抑制效果更为明显，下游立管两向位移响应变化不大，但由于上游立管顺向变形进一步增大导致立管最小间距较光滑立管有所减小，增大了碰撞机率。这说明对抑制措施在立管群中的应用需要做更加深入的研究，寻求合理的方式以提高抑制措施对各立管的抑制效率，或者寻求更加高效适用的方法削弱立管涡激振动损伤，提高立管安全使用寿命。