刚性联结串列双圆柱尾流致涡激振动研究

2022-02-15杜晓庆朱红玉吴葛菲林伟群

振动工程学报 2022年6期

杜晓庆，朱红玉，吴葛菲，林伟群，赵燕

（1.上海大学力学与工程科学学院土木工程系，上海 200444；2.上海大学风工程和气动控制研究中心，上海 200444；3.台州学院建筑工程学院，浙江台州 318000）

引言

大长细比柱群在实际工程中应用广泛，如大跨度桥梁的并列索、多分裂导线、海洋立管等［1-2］。当两根或多根圆柱并排布置时，下游圆柱受到上游圆柱尾流的影响会产生强烈的尾流激振现象。通常采用刚性联结器将多根圆柱联结，以抑制多圆柱结构的尾流激振［3-6］。目前针对无联结双圆柱的研究较多，而对刚性联结双圆柱的研究较少，刚性联结对双圆柱流致振动的抑振效果和控制机理尚未被澄清。

对于无联结串列双圆柱的流致振动，学者通过试验和数值模拟，探讨了圆柱间距、质量比、雷诺数、振动自由度、阻尼比等参数的影响。Du 等［7］通过风洞试验，研究了间距比P/D=4（D为圆柱直径，P为圆心间距）的串列和错列双圆柱的尾流致涡激振动，研究发现下游圆柱的尾流致涡激共振振幅均大于单圆柱涡激振动振幅。Du 等［8］采用风洞试验，在雷诺数Re=18000～168800、间距比P/D=4 条件下，研究了表面粗糙度和雷诺数对串列和错列双圆柱尾流致涡激振动的影响，研究发现：增加圆柱的表面粗糙度会明显改变尾流致涡激振动的动力特性，雷诺数对下游圆柱尾流致涡激振动的不稳定性具有强烈影响。Carmo 等［9］在雷诺数Re=150 和300、质量比m*=2 的条件下，研究了间距比（P/D=1.5～8）对串列双圆柱的涡激振动的影响，研究表明：对于间距比P/D=3～8 的工况，当结构发生尾流致涡激共振时，下游圆柱横流向最大振幅约为单圆柱涡激振动的1.5 倍；此外，下游圆柱在较高的折减速度下仍会保持较大的振幅。郭晓玲等［10］在雷诺数Re=150、间距比P/D=3～8、质量比m*=5.0～20.0 情况下，探讨了质量比、圆心间距及折减速度对下游圆柱尾流致涡激振动的影响。Nguyen 等［11］采用DES 方法在雷诺数为103～105，间距比P/D=4 和5 时，研究了串列双圆柱时的尾流致涡激振动，发现雷诺数对圆柱的振幅和振动频率均有影响。Borazjani 等［12］在Re=200，m*=2 和P/D=1.5 条件下，研究了振动自由度对串列双圆柱（上、下游圆柱均弹性支撑）涡激振动的影响，发现两自由度振动体系（横流向和顺流向）下游圆柱的横流向振幅明显小于单自由度振动体系（横流向振动）。

以往对刚性联结双圆柱尾流致涡激振动的研究较少。Zhao［13］在雷诺数Re=150、质量比m*=2 条件下，研究了串列和并列双圆柱在横流向单自由度振动体系的流致振动问题。研究表明：间距比对串列双圆柱流致振动有着重要的影响；此外，刚性联结并列双圆柱还会出现弛振现象。上述研究中，圆柱只有横流向振动的自由度，目前仅有Gao 等［6］考虑了横流向和顺流向两自由度的振动，在雷诺数Re＝150、低质量比m*=2.6 条件下，研究了来流攻角和间距比对刚性联结双圆柱流致振动的影响。

以往研究表明［14］，根据雷诺数不同，静止双圆柱的临界间距在3.5D～5.0D之间变化；且双圆柱的流致振动在该间距范围内的振动特性和流场特征较为复杂［15］；此外，该间距范围在实际工程中具有广泛应用，如中国西堠门大桥的吊索间距为3.4D～10D。

为了进一步探讨刚性联结对较大质量比双圆柱振动特性的影响规律及其流场机理，本文以圆心间距为4D的无联结及刚性联结串列双圆柱为研究对象，考虑横流向和顺流向两个自由度影响，在雷诺数Re=150，m*=20 时，采用数值模拟方法，研究了两类双圆柱发生尾流致涡激振动的动力响应特性和流场特征，探讨了刚性联结对双圆柱尾流激振的减振效果及其流场作用机理。

1 数值方法

1.1 计算模型与运动方程

图1给出了无联结及刚性联结串列双圆柱的计算模型示意图。本文将串列双圆柱流致振动简化为质量-弹簧-阻尼系统。两个圆柱的中心间距为4D（D为圆柱直径），上、下游圆柱的物理和几何参数分别用脚标1 和2 标记，无联结串列双圆柱的上、下游圆柱分别可做横流向和顺流向振动，刚性联结串列双圆柱可做横流向及顺流向振动。需要说明的是，本文中各工况下圆柱横流向与顺流向具有相同的自振频率。

图1 双圆柱计算模型Fig.1 Computational model of twin circular cylinders

对于无联结串列双圆柱，其运动方程为：

对于刚性联结串列双圆柱，其运动方程为：

式中xi，和分别为圆柱顺流向的瞬时位移、速度和加速度；yi，和分别为横流向的瞬时位移、速度和加速度；m为单位展向长度圆柱的质量；ζ为结构的阻尼比；ωn为圆柱的固有圆频率。FD（t）=0.5ρU2DCD和FL（t）= 0.5ρU2DCL分别为作用在单位长度圆柱上的顺流向及横流向气动力，其中，CD和CL分别为阻力和升力系数，ρ为空气密度。

1.2 流固耦合算法

本文通过动网格技术来实现数值模拟中圆柱与流场之间的流固耦合，过程如下：

（1）运用数值模拟方法求解流体控制方程，获得流场速度场、压力场及圆柱表面流体力；

（2）将流体力作用于两自由度振动的圆柱，以四阶Runge-Kutta 法求解圆柱的运动控制方程（式（1）～（6）），得到振动圆柱顺流向和横流向的动力响应；

（3）通过动网格技术，将振动圆柱的振动速度传递于网格系统，更新网格位置；

（4）返回第（1）步开始计算下一个时间步的响应，如此循环获得各时间步圆柱的动力响应，实现上述流固耦合算法。

1.3 网格划分及计算参数

双圆柱的计算域及网格划分如图2和3 所示。由图2可知，坐标原点O位于双圆柱的正中间，入口边界距坐标原点O为30D，上边界、下边界距坐标原点O为30D，出口边界距坐标原点O为35D。本文对计算域网格划分进行分块处理。所有工况圆柱近壁面采用结构化网格；圆柱限制在红色虚线圆形区域内运动，在此区域内采用三角形非结构化网格；红色虚线圆形区域外均采用结构化网格。

图2 计算域和边界条件示意图Fig.2 Schematic drawing of computational domain and boundary conditions

边界条件的设置如图2所示：入口边界设为速度入口（velocity-inlet）边界条件；出口边界设为自由出流（outflow）边界条件；两侧壁面采用对称（symmetry）边界条件；圆柱壁面采用无滑移壁面（wall）边界条件。本文在具体的数值计算中，压力和速度耦合采用SIMPLEC 法求解，动量方程采用二阶精度的离散格式。

图3 计算域网格方案Fig.3 Computational domain mesh scheme

本文的计算参数如下：串列双圆柱的圆心间距为4D，雷诺数Re=150，折减速度Vr=U/（fnD）=3～12（U为来流风速，fn为圆柱的自振频率），圆柱的质量比m*=m/（0.25ρπD2）=20（m为圆柱的单位长度质量，ρ为空气密度），参考文献［12，16］，为使圆柱较早达到振动稳定状态，节省计算时间，结构阻尼比ζ取0。

2 模型验证

为保证计算结果的可靠性，首先针对固定单圆柱，研究了周向网格数量、无量纲时间步长和阻塞率等参数对计算结果的影响，并与文献［17-18］结果进行对比，确定无量纲时间步长Δt*=0.005（Δt*=ΔtU/D）、阻塞率为1.67%、周向网格数为200，详细的结果验证工作可参考文献［19］。

进一步针对单圆柱的双自由度涡激振动问题，考虑了三种网格方案对计算结果的影响，网格方案如表1所示。其中，运动单圆柱的质量比m*=10，结构阻尼比ζ=0，折减速度约为6。将顺流向的平均振幅Xmean/D、顺流向的脉动振幅Xrms/D、横流向的振幅Ymax/D、平均阻力系数CD，mean、脉动阻力系数CD，rms和St=fsD/U（fs为圆柱的涡脱频率）与现有文献［16，20-21］的结果进行对比。从表1可以看出，Mesh 1～Mesh 3 网格下的计算结果随着网格数量的增加，计算结果趋于收敛，综合考虑计算资源及计算精度，本文网格模型选取Mesh 2。

表1 单圆柱模型网格方案和结果验证（Re = 100）Tab.1 Grid schemes and results verification of a vibrating single circular cylinder model for Re = 100

图4给出了圆柱涡激振动的验证结果。图4（a）为单圆柱横流向振幅Ay/D随折减速度的变化曲线，并与文献［22-24］的研究结果进行对比。此外，本文还选取刚性联结双圆柱单自由度涡激振动做进一步验证，图4（b）给出了刚性联结双圆柱单自由度涡激振动横流向振幅与文献［13］结果的对比。以上结果显示，本文计算结果与文献值吻合较好。

图4 单圆柱和刚性联结双圆柱的涡激振动结果验证Fig.4 Result verification of the vortex induced vibration of the single circular cylinder and rigidly coupled circular cylinders

3 结果分析

3.1 振动响应特性

图5给出了无联结、刚性联结串列双圆柱及单圆柱的无量纲振幅（Ay/D，Ax/D）随折减速度的变化曲线，其中Ax=（Xmax-Xmin）/2；Ay=（Ymax-Ymin）/2；Xmax，Xmin分别为顺流向位移的最大值和最小值；Ymax，Ymin分别为横流向位移的最大值和最小值。

图5 串列双圆柱及单圆柱振幅随折减速度的变化Fig.5 Variation of the amplitude for the double circular cylinders and the single circular cylinder with the reduced velocity

由图5可知，串列双圆柱和单圆柱均以横流向振动为主，顺流向振动的振幅较小。无联结时，上游圆柱与单圆柱类似，在整个折减速度范围内的顺流向振幅均很小；下游圆柱的顺流向振幅也较小，但当Vr=6.5 时，下游圆柱的顺流向振幅达到0.08D；上、下游圆柱的起振风速与单圆柱相同，均为Vr=4；上游圆柱的横流向最大振幅略大于单圆柱；而下游圆柱由于受到上游圆柱尾流的干扰，其横风向最大振幅达到1.1D左右（Vr=7），远大于单圆柱。

当上、下游圆柱之间进行刚性联结后，在Vr=6 时，双圆柱的顺流向振幅达到了0.23D，这与单圆柱和无联结串列双圆柱有很大不同；双圆柱的起振风速Vr= 5，大于单圆柱及无联结上、下游圆柱，在共振范围内双圆柱的振幅约为0.6D～0.8D，略大于单圆柱及无联结上游圆柱。与无联结下游圆柱相比，刚性联结串列双圆柱的振幅极值及出现共振的风速范围均有所减小，说明刚性联结对振动有一定的抑制作用。

图5还给出了无联结上、下游圆柱和刚性联结串列双圆柱典型工况下的振动轨迹。在锁振区域外，三者的振动轨迹均为“8”字形。在锁振区域内，对于无联结上、下游圆柱，振动轨迹出现了不规则情况；对于刚性联结串列双圆柱，在锁振区域内出现了类似于椭圆的闭合圆形振动轨迹。同时，在最大振幅时，三者的振动轨迹均呈“8”字形。

以往学者认为，当旋涡脱落频率与固有频率接近时（即fy / fn趋近于1.0），圆柱发生振动锁定［25］，参考文献［10］，本文以0.95 ≤fy / fn≤1.05 为锁振区域。图6给出了无联结串列双圆柱、刚性联结串列双圆柱及单圆柱横流向振动频率比随折减速度的变化曲线。由图6可知，在各折减速度下，无联结上、下游圆柱的振动主频均相同，锁振区域的折减速度范围均为Vr=4.5 ～ 7，其锁振区域的折减速度范围比单圆柱小，表明尾流致涡激振动的锁振区域较单圆柱涡激振动的锁振区域小。对于刚性联结串列双圆柱，其锁振区域的折减速度范围为Vr=5.5～7，明显小于单圆柱及无联结串列双圆柱，说明刚性联结对尾流致涡激振动有一定抑制作用。

图6 横流向振动频率比随折减速度的变化Fig.6 The vibration frequency ratio as a function of reduced velocity in the transverse direction

值得注意的是，这与Gao等［16］得到的结果有一定区别，其研究结果表明刚性联结串列双圆柱的锁振区域大于单圆柱的锁振区域。造成结果差异的原因可能是因为质量比的不同，Gao 等［6］研究中的质量比m*= 2.6，而本文的质量比较高。这与杨骁等［17］对单圆柱及无联结下游圆柱质量比效应的分析结果类似，即随着质量比的增加，圆柱的锁振区域会减小。

3.2 气动力特性

图7和8 分别给出了无联结及刚性联结串列双圆柱的平均阻力系数、脉动阻力系数和脉动升力系数随折减速度的变化曲线。

由图7可知，对于无联结串列双圆柱，由于“屏蔽”作用，下游圆柱的平均阻力系数在多数情况下小于上游圆柱，但在Vr= 7 时，出现了下游圆柱平均阻力大于上游圆柱的情况，故此时出现了下游圆柱的顺流向平衡位置偏移量大于上游圆柱的现象（限于篇幅未给出）；对于刚性联结串列双圆柱，下游圆柱的平均阻力系数在各个折减速度下都比上游圆柱要小。

图7 平均阻力系数随折减速度的变化Fig.7 The mean values of the drag coefficients as a function of reduced velocity

由图8（a）可见，刚性联结串列双圆柱与无联结串列双圆柱类似，上、下游圆柱的脉动阻力系数随折减速度先增加后减小，最后趋于稳定，这与顺流向振幅的变化趋势基本一致。由图8（b）可知，脉动升力系数的变化较为复杂。对于无联结串列双圆柱，多数情况下上游圆柱的脉动升力系数小于下游圆柱，当Vr= 5 和Vr=5.5 时，上游圆柱的脉动升力系数大于下游圆柱，结合图5可知，此时上游圆柱达到横流向振幅最大值。对于刚性联结串列双圆柱，下游圆柱的脉动升力系数均大于上游圆柱。刚性联结双圆柱的横流向振动与上、下游圆柱升力之和的脉动值有关，但当横流向振幅达到极值时（Vr= 6.5），升力相加的脉动值却低于0.5，这说明振幅不仅仅与脉动气动力的大小有关，还与气动力的频率相关。

图8 脉动气动力系数随折减速度的变化Fig.8 The r.m.s.values of the aerodynamic coefficients as a function of reduced velocity

3.3 位移、升力时程及其功率谱

为进一步分析横流向位移与气动升力的关系，图9给出了无联结下游圆柱及刚性联结串列双圆柱在一些典型折减速度下的位移、升力的时程曲线及其功率谱图。图中Ay和ACL分别为位移及升力的功率谱幅值，下标1 和2 分别代表上、下游圆柱。对于无联结下游圆柱，如图9（a）所示，在折减速度不断增大的过程中，位移与升力的时程曲线从同相变为反相。从功率谱图可知，在最大横流向振幅时，升力有两个峰值，其三倍固有频率处分量大于固有频率处的分量，但位移的主频仍然在固有频率处，这与锁振区域内的特性相符合。对于刚性联结串列双圆柱，如图9（b）所示，总的规律与无联结下游圆柱相似，但在最大横流向振幅时的功率谱图存在差异，刚性联结串列双圆柱升力的主频在固有频率处。

图9 圆柱的横流向振动和升力系数的时程曲线及功率谱图Fig.9 Time histories and amplitude spectra of oscillation in the transverse direction and lift coefficient of the cylindes

3.4 能量输入机制

图10和11 分别给出了在无联结上、下游圆柱及刚性联结串列双圆柱在各自振幅最大时的瞬时能量输入（其中，正值表示气动升力为结构输入能量，负值表示气动升力消耗结构的能量，能量输入用无量纲功率P*表示，P*=CL(t)⋅v(t)/U，v(t)表示圆柱瞬时的横流向速度，U为来流风速）。因为此时圆柱的振动以横流向为主，所以仅考虑气动升力系数及横流向振幅。

从图10的瞬时能量输入图可以看出：对于无联结串列双圆柱，上、下游圆柱的瞬时能量输入时程曲线均具有周期性，且上游圆柱的幅值很小。在一个振动周期内，无联结下游圆柱能量输入和耗散的时程曲线不对称，这可能是由上游圆柱的振动使流场不对称所引起的。在每个振动周期内，圆柱的气动升力输入和耗散的能量相等，所以其横流向振幅很稳定。

图10 无联结串列双圆柱的瞬时能量输入（Vr = 7）Fig.10 Transient energy input in transverse direction for the uncoupled circular cylinders in tandem arrangement（Vr = 7）

从图11的瞬时能量输入图可以看出：对于刚性联结串列双圆柱，上、下游圆柱的能量输入时程曲线都具有周期性。在一个振动周期内，上游圆柱的气动升力几乎都在耗散系统的能量，而下游圆柱的气动升力整体上为系统输入能量；在每个振动周期内，上游圆柱气动升力耗散的能量和下游圆柱气动升力输入的能量几乎相等，所以其横流向振幅很稳定。

图11 刚性联结串列双圆柱的瞬时能量输入（Vr = 6.5）Fig.11 Transient energy input in transverse direction for the rigidly coupled circular cylinders in tandem arrangement（Vr = 6.5）

3.5 流场特征

为了探讨无联结及刚性联结串列双圆柱尾流激振的流固耦合机制，进一步对无联结下游圆柱及刚性联结串列双圆柱在横流向振幅最大时所对应的折减速度下的流场形态进行分析，研究了单个振动周期内上游圆柱脱落的涡与下游圆柱之间的相互作用。

图12给出了Vr= 7 时无联结串列双圆柱的涡量图、横流向位移和能量输入的时程曲线。从图12（a）可以发现，此时下游圆柱运动到最高点，上游圆柱脱落的涡撞击在下游圆柱的下侧，并与下游圆柱下侧脱落的涡相互融合，使下游圆柱下侧的旋涡强度增强，促进了下游圆柱的振动。

图12 无联结串列双圆柱振幅最大时的涡量图（Vr ＝7）Fig.12 Vortex diagram of two uncoupled tandem circular cylinders at the maximum tranverse amplitude（Vr ＝7）

图13给出了Vr= 6.5 时刚性联结串列双圆柱的涡量图、横流向位移和能量输入的时程曲线。由图13（a）可知，当下游圆柱向上运动还未达到最高点时，上游圆柱下侧剪切层脱落的涡撞击到下游圆柱表面并分解成两个子涡，分别从下游圆柱两侧向尾流移动。从下游圆柱两侧通过的涡对下游圆柱产生相反的作用力，一定程度上抑制了下游圆柱的振动。