2019, Vol. 47
2. 上海大学 风工程和气动控制研究中心，上海 200072
2. Aerodynamic Flow Control Research Center, Shanghai University, Shanghai 200072, China
大长细比圆柱型结构在实际工程中有大量的应用，如桥梁并列索、超高层建筑、多分裂导线、冷却塔、烟囱等，并且多以柱群形式出现.多个圆柱结构间存在强烈的气动干扰，下游圆柱在横风向升力作用下常发生尾流激振现象，导致结构发生损伤或破坏[1-4]，澄清双圆柱结构下游圆柱的升力机理有助于理解尾流激振的发生机制.
已有研究发现[5-7]，对于小间距(P/D＝1.2~3.5，其中P为圆心间距，D为圆柱直径)双圆柱，下游圆柱在小风向角(β≈8°)下受到很大的平均升力，称之为内侧升力；对于中等间距双圆柱(P/D=3.0~5.0)，下游圆柱在β=20°左右同样受到较大的平均升力，称之为外侧升力.
对于内侧升力，低雷诺数下的流迹显示试验表明[6, 8]：两圆柱间的加速间隙流、下游圆柱的停滞点偏移以及上游圆柱的剪切流在下游圆柱表面发生再附是升力出现的三个原因.杜晓庆等[9]进行了基于高雷诺数下的大涡模拟研究，对下游圆柱平均升力的出现给出了不同的解释.
对于外侧升力，Zdravkovich[7]和Ting等[10]提出，两圆柱间隙气流速度的增大使下游圆柱间隙侧表面压力减小，导致平均升力的出现，文献[10-12]中也认为下游圆柱风压停滞点的偏移对平均升力有贡献.Alam等[13-14]提出，受上游圆柱尾流旋涡的作用，下游圆柱的旋涡加速向下游流动，这也会引起下游圆柱间隙侧表面压力的减小，进而导致平均升力的出现.
对以往文献的分析可知，由于影响双圆柱的因素众多，干扰机理复杂，研究者对下游圆柱的升力机理尚未有统一的认识.值得注意的是，以往的机理解释都是基于较低雷诺数下的流迹显示试验得出的，但双圆柱绕流有强烈的雷诺数效应[3, 15]，因而在高雷诺数下很可能存在不同的干扰机理.此外，以往文献主要集中于对双圆柱平均升力的研究，受到测试条件的限制，对脉动升力及其发生机理的研究很少.
本文中采用大涡模拟(LES)方法，在高雷诺数(Re=1.4×105)下，研究了不同风向角时中等间距(P/D=4)双圆柱的气动力特性和流场特性，重点讨论了下游圆柱的气动升力与流场结构之间的内在关系，提出了双圆柱干扰流态随风向角的变化规律，进一步澄清下游圆柱的平均升力和脉动升力的流场机理.
1 数值方法和计算模型 1.1 控制方程和亚格子尺度模型在大涡模拟方法中，大尺度的旋涡由滤波后的Navier-Stokes方程直接求解，而小尺度的旋涡则采用亚格子尺度(SGS)模型模拟.与雷诺平均法(RANS)相比，大涡模拟方法可更好地模拟流场中的湍流旋涡，能捕捉到更丰富的流场脉动信息.因此，为了准确地模拟错列双圆柱周围的流场特性，采用大涡模拟方法进行研究.经过滤波函数的滤波，可得到大尺度旋涡的不可压缩Navier-Stokes方程，如下所示：
| $ \frac{{\partial {{\bar u}_i}}}{{\partial {x_i}}} = 0 $ | (1) |
| $ \frac{{\partial {{\bar u}_i}}}{{\partial t}} + \frac{{\partial {{\bar u}_i}{{\bar u}_j}}}{{\partial {x_j}}} = - \frac{1}{\rho }\frac{{\partial \bar p}}{{\partial {x_i}}} + \mu \frac{{{\partial ^2}\bar u}}{{\partial {x_j}\partial {x_j}}} - \frac{{\partial {\mathit{\boldsymbol{\tau }}_{ij}}}}{{\partial {x_j}}} $ | (2) |
式中：ui为滤波后速度；ρ为流体密度；p为滤波后压力；μ为流体动力黏度；τij为亚格子应力张量.亚格子应力张量τij采用Smagorinsky-Lilly的亚格子尺度模型，具有以下形式：
| $ {\mathit{\boldsymbol{\tau }}_{ij}} - \frac{1}{3}{\mathit{\boldsymbol{\tau }}_{kk}}{\mathit{\boldsymbol{\delta }}_{ij}} = - 2{\mu _t}{{\bar S}_{ij}} $ | (3) |
| $ {{\bar S}_{ij}} = \frac{1}{2}\left( {\frac{{\partial {{\bar u}_i}}}{{\partial {x_j}}} + \frac{{\partial {{\bar u}_j}}}{{\partial {x_i}}}} \right) $ | (4) |
式中:μt为亚格子尺度的湍动黏度; δij为单位张量.μt的计算式如下所示：
| $ {\mu _t} = {({C_{\rm{S}}}\mathit{\Delta })^2}\left| {\bar S} \right| $ | (5) |
| $ \left| {\bar S} \right| = \sqrt {2{{\bar S}_{ij}}{{\bar S}_{ij}}} $ | (6) |
| $ \mathit{\Delta } = {({\mathit{\Delta }_x}{\mathit{\Delta }_y}{\mathit{\Delta }_z})^{\frac{1}{3}}} $ | (7) |
式中:Δi为沿坐标轴i方向的网格尺寸；CS为Smagorinsky常数，取0.1.
1.2 计算模型和计算参数错列双圆柱的布置形式如图 1所示.图 1中:CD1、CL1分别为上游圆柱的阻力和升力；CD2、CL2分别为下游圆柱的阻力和升力；θ为圆柱表面风压测点位置的角度，以圆柱迎风点为零点，顺时针旋转为正方向.两个圆柱的直径均为D，圆柱间距P为4D，风向角β为自由来流与圆柱中心连线的夹角.为了掌握不同风向角对圆柱周围流场特性和圆柱气动力特性的影响规律，计算的风向角β分别为0°、5°、10°、15°、20°、25°、30°、45°、60°、75°、90°.
|
图 1 双圆柱计算模型 Fig.1 Computational model of two circular cylinders |
计算域网格采用结构化网格、O型域，径向长度为46D，展向长度为2D，均分为20层网格，计算模型总网格数约为266万，近壁面最小网格厚度为1×10-4D，最大阻塞率为4.3%(β=90°).量纲一时间步Δt*为0.005(Δt*=ΔtU0/D，其中Δt为实际计算时间步，U0为来流风速).计算域采用均匀速度入口边界条件，自由出口边界条件，展向两端采用周期性边界条件，圆柱表面采用无滑移壁面边界条件，如图 2所示.计算流体动力学(CFD)数值模拟基于Fluent14.0程序实现，采用压力耦合方程的半隐式(SIMPLEC)格式求解压力速度耦合方程组，空间离散采用中心差分格式，时间离散采用二阶全隐格式.
|
图 2 计算域和网格 Fig.2 Computational domain and grid |
在确定双圆柱计算模型和计算参数之前，先以单圆柱为研究对象，对六种计算模型(分别命名为Case 1~6)进行数值模拟，研究周向网格数量、时间步和展向长度等参数对计算结果的影响[16].表 1给出了单圆柱各计算模型的参数、平均阻力系数、脉动升力系数和St数，也列出了文献[17]和[18]中的风洞试验结果.通过比较可以看出，模型Case 6计算得到的气动力系数与文献的吻合情况最好，因此采用该模型参数进行错列双圆柱的大涡模拟研究.
2.2 气动力系数
图 3是下游圆柱平均气动力系数曲线.由图 3可见，下游圆柱的平均升力系数和平均阻力系数随风向角的变化趋势与文献较为吻合.与单圆柱相比，在风向角0°~30°内，下游圆柱的平均阻力系数均小于单圆柱，并且随风向角的增大而增大.从β=30°开始，下游圆柱的平均阻力系数趋于平缓，其值比单圆柱稍大.
|
图 3 下游圆柱平均气动力系数随风向角的变化 Fig.3 Mean aerodynamic coefficients as a function of incidence angles for downstream cylinder |
对于平均升力系数来说，随着风向角的增大，下游圆柱的平均升力系数先减小后增大，最后趋于稳定.在β=0°~45°时，下游圆柱的平均升力系数均为负值，并且在β=20°左右出现升力的极小值，即出现圆柱上侧升力峰值.当β≥45°时，平均升力系数便与单圆柱相似并接近于零.图 4是下游圆柱脉动气动力系数随风向角的变化曲线.由图 4可以看出，在β=0°~5°时下游圆柱的脉动升力最大，其值达到单圆柱的2倍，并且在β=20°~25°时达到极小值.当β>30°后，下游圆柱的脉动升力系数开始缓慢减小，但依然明显高于单圆柱；下游圆柱的脉动阻力则在β=10°左右达到最大值，并且在β=0°~30°范围内都明显大于单圆柱.
|
图 4 下游圆柱脉动气动力系数随风向角的变化 Fig.4 Root mean square aerodynamic coefficients as a function of incidence angles for downstream cylinder |
图 5为下游圆柱在不同风向角下的平均风压系数曲线.从图 5可以看出：当β=0°时，整个下游圆柱表面都为负压，并具有良好的对称性；当β=5°时，下游圆柱表面负压远小于其他工况，并且停滞点偏移了20°左右；在β=15°~30°时，下游圆柱表面平均风压呈现出明显的不对称性，圆柱下侧(间隙侧)表面受到的负压绝对值明显大于上侧，而平均风压分布的不均匀则导致圆柱表面受到平均升力；当风向角大于45°时，下游圆柱表面平均风压逐渐接近单圆柱，又重新表现出良好的对称性.
|
图 5 下游圆柱表面平均风压系数分布 Fig.5 Mean pressure coefficient distribution of downstream cylinder |
图 6为下游圆柱在不同风向角下的脉动风压系数曲线.由图 6可知：几乎所有风向角下下游圆柱的脉动风压系数值都不同程度地大于单圆柱，这说明受上游圆柱的影响，下游圆柱表面的脉动风压发生剧烈变化；在β=5°时，下游圆柱表面的脉动风压变化最为剧烈，也具有一定的对称性；在β=15°~30°时，下游圆柱表面的脉动风压呈现出明显的不对称性，下游圆柱上侧表面的脉动风压变化较为平缓，而下侧(间隙侧)表面的脉动风压变化很大，脉动风压的不对称性也会促使下游圆柱受到平均升力的作用；在β=45°~90°时，随着风向角的增大，下游圆柱表面的脉动风压逐渐接近于单圆柱，并且重新具有良好的对称性.
|
图 6 下游圆柱表面脉动风压系数分布 Fig.6 Root mean square pressure coefficient distribution of downstream cylinder |
对于错列双圆柱来说，不同风向角下的St差别不大，都在0.17左右，但功率谱(PSD)峰值却差别很大.图 7是不同风向角下下游圆柱的升、阻力功率谱曲线.从图 7可以看出，下游圆柱升、阻力功率谱随风向角的变化与如图 4所示下游圆柱脉动气动力系数的变化趋势一致.
|
图 7 气动力系数的功率谱峰值随风向角的变化 Fig.7 PSD peak of aerodynamic coefficient as a function of incidence angles |
当β=0°和5°时，下游圆柱的升力功率谱峰值最大，说明此时涡脱强度很大，圆柱表面出现很强的脉动风压，进而产生很大的脉动升力；当β=20°时，下游圆柱的升力功率谱峰值最小，表明此时的涡脱强度明显减弱，相应的下游圆柱表面受到的脉动升力也随之降低.
2.5 干扰流态分类通过上述分析可以发现，不同的风向角范围内，下游圆柱的气动力有着很大差异.通过观察不同风向角时下游圆柱周围的流场变化，可以总结出以下四种不同的流态(见图 8)：
|
图 8 错列双圆柱在不同风向角下的干扰流态 Fig.8 Interaction flow patterns of two staggered circular cylinders at different incidence angels |
(1) 流态1，旋涡撞击流态.当β=0°时，上游圆柱上侧脱落的旋涡撞击到下游圆柱上侧，接着上游圆柱下侧脱落的旋涡撞击到下游圆柱的下侧(间隙侧).上游圆柱脱落的旋涡对下游圆柱的不断撞击是下游圆柱表面出现较大脉动风压和脉动升力的主要原因.
(2) 流态2，旋涡撞击和剪切层干扰流态.当β=5°~10°时，上游圆柱上侧脱落的旋涡撞击到下游圆柱的迎风面，上游圆柱下侧脱落的旋涡并没有撞击到下游圆柱，而是与下游圆柱下侧的剪切层相互作用.同时存在的旋涡撞击和剪切层干扰，导致下游圆柱表面出现很大的脉动升力.
(3) 流态3，剪切层干扰流态.当β=15°~30°时，上游圆柱上侧脱落的旋涡与下游圆柱下侧(间隙侧)的剪切层和旋涡耦合产生相互作用，而上游圆柱下侧脱落的旋涡在向下游运动的过程中由于距离下游圆柱较远，故没有发生明显的相互作用.
(4) 流态4，尾流干扰流态.当β=45°~90°时，上游圆柱与下游圆柱间距较大，上游圆柱上侧脱落的旋涡在向下游运动的过程中只有一部分尾流会与下游圆柱发生相互作用，而下侧脱落的旋涡则基本不会对下游圆柱有影响.
2.6 平均升力的机理分析为了进一步澄清风向角20°附近下游圆柱出现平均升力峰值(即文献[15]中的外侧升力)的流场机理，本节从平均流场和瞬态流场两方面进行分析.
2.6.1 平均流场分析图 9是β=20°时下游圆柱的平均风压场和平均流速.从图 9可以看出，下游圆柱下侧(间隙侧)的负压强度明显大于上侧，并且从下游圆柱的表面平均风压分布可知，下游圆柱的风压停滞点相对于来流向上侧偏移了6°左右，这两个因素对下游圆柱的平均升力均有一定贡献.此外，以往的学者还认为[8]，加速间隙流也是导致下游圆柱出现平均升力的原因之一，但从平均风速图(见图 9b)和下游圆柱附近的局部风速图(见图 9c)可看出，两个圆柱间隙中的风速并未有明显的加速现象，并且下游圆柱上、下侧表面的最大平均风速也未见明显的差别.上述结果说明，在该工况下间隙流并未有明显的加速现象出现，因而对下游圆柱平均升力的贡献不大.
|
图 9 平均风压场、平均流速和下游圆柱的局部流速 Fig.9 Mean pressure coefficient field, mean velocity ratio and local velocity ratio of downstream cylinder |
从第2.6.1节分析可知，下游圆柱上、下侧表面的平均风压存在较大差异，这是下游圆柱受到较大平均升力的主要原因之一.为了进一步澄清该现象的流场机理，图 11给出了风向角为20°时九个瞬态时刻(T1~T9)的涡量图及下游圆柱表面的风压系数分布.这九个瞬态时刻在升力系数时程曲线上的具体位置如图 10所示，下游圆柱的升力系数在T1和T9时刻接近极大值，在T6时刻则接近极小值.
|
图 10 β=20°时升力系数的时程曲线及瞬态时刻 Fig.10 Lift coefficient field and instantaneous moment at β=20° |
|
图 11 β=20°时一个周期内的瞬时涡量及表面风压分布 Fig.11 Instantaneous vorticity and surface pressure distribution in a cycle at β=20° |
由图 10还可知，上、下游圆柱的升力时程曲线均有明显的周期性，但两者有一定的相位差；对于下游圆柱而言，升力系数从波峰下降至波谷的时间较长，而从波谷上升至波峰的时间较短.
由图 11可见，上游圆柱的上、下表面交替形成分离的剪切层，剪切层受另一侧旋涡的影响而发生卷起，卷起的剪切层最终脱离圆柱并在尾流中形成方向交替变化的旋涡.随着剪切层的变化，上游圆柱表面的风压也发生相应的交替变化.上游圆柱的绕流场和表面风压的变化与单圆柱相似.
下游圆柱受到上游圆柱尾流旋涡的影响，导致其气动力特性与单圆柱有很大差别.由图 11可见，从上游圆柱脱落的旋涡与下游圆柱的下侧剪切层(或旋涡)发生相互作用.相对来说，上游圆柱上侧脱落的旋涡与下游圆柱的相互作用更为剧烈(T1、T2、T8和T9时刻)，而上游圆柱下侧脱落的旋涡对下游圆柱的影响则相对较小(T4、T5和T6时刻).
特别在T8和T9时刻，当上游圆柱上侧脱落的旋涡移动至下游圆柱的下侧剪切层附近时，由于旋涡与剪切层的旋转方向不同，上游圆柱的上侧旋涡吸引下游圆柱的下侧剪切层从圆柱表面脱落而形成旋涡；随后，这一对方向相反的旋涡加速向下游移动(T9时刻).由于下游圆柱的下侧旋涡迅速离开圆柱，使得下侧旋涡对下游圆柱上侧剪切层的影响减弱，从而导致下游圆柱上侧剪切层的卷起现象不明显.在T1和T2时刻能观察到类似现象.因而，下游圆柱的上侧剪切层在所有时刻中均保持较为平顺的形态，上侧剪切层对圆柱的影响相对较小，所以下游圆柱上侧表面的负压在所有时刻中均保持在相对较低的数值，最终导致下游圆柱上、下侧表面的平均风压出现较大差异，这是下游圆柱受到较大平均升力作用的原因之一.
此外，由图 11还可见，T4~T7时刻的下游圆柱风压停滞点明显向圆柱上侧偏移，而T1和T9时刻则向下侧偏移，T2、T3和T8时刻则偏移不明显，因而总体上下游圆柱的停滞点是向上偏移，这是导致下游圆柱受到较大平均升力的另一个原因.
3 结论(1) 随着风向角的增大，两个圆柱的干扰流态依次为旋涡撞击流态、旋涡撞击和剪切层干扰流态、剪切层干扰流态、尾流干扰流态.在不同流态下，下游圆柱的气动力特性有很大差异.
(2) 在风向角为0°~10°的旋涡撞击流态以及旋涡撞击和剪切层干扰流态时，从上游圆柱脱落的旋涡与下游圆柱发生强烈的撞击，导致下游圆柱的脉动升力和迎风面的脉动风压远远大于其他流态，也远大于单圆柱.
(3) 在风向角为15°~30°的剪切层干扰流态下，下游圆柱受到显著的平均升力作用.下游圆柱风压停滞点的偏移、上游圆柱的尾流旋涡与下游圆柱间隙侧剪切层和旋涡的相互作用，是下游圆柱出现平均升力的两个主要原因.
| [1] |
沈国辉, 王宁博, 楼文娟, 等. 渡桥电厂冷却塔倒塌的塔型因素分析[J]. 工程力学, 2012, 29(8): 123 SHEN Guohui, WANG Ningbo, LOU Wenjuan, et al. Analysis of tower shape factor in the collapse of the ferrybridge cooling towers[J]. Engineering Mechanics, 2012, 29(8): 123 |
| [2] |
杜晓庆, 蒋本建, 代钦, 等. 大跨度缆索承重桥并列索尾流激振研究[J]. 振动工程学报, 2016, 29(5): 842 DU Xiaoqing, JIANG Benjian, DAICHIN, et al. Wake-induced vibration research on parallel cables of long-span cable-supported bridges[J]. Journal of Vibration Engineering, 2016, 29(5): 842 |
| [3] |
SUMNER D. Two circular cylinders in cross-flow: a review[J]. Journal of Fluids and Structures, 2010, 26(6): 849 DOI:10.1016/j.jfluidstructs.2010.07.001 |
| [4] |
CHEN S S, PAIDOUSSIS M P, AU-YANG M K, et al. Flow-induced vibration of circular cylindrical structures[J]. Journal of the Acoustical Society of America, 1998, 75(1): 299 |
| [5] |
ZDRAVKOVICH M M. The effects of interference between circular cylinders in cross flow[J]. Journal of Fluids and Structures, 1987, 1(2): 239 DOI:10.1016/S0889-9746(87)90355-0 |
| [6] |
ALAM M M, MEYER J P. Global aerodynamic instability of twin cylinders in cross flow[J]. Journal of Fluids and Structures, 2013, 41: 135 DOI:10.1016/j.jfluidstructs.2013.03.007 |
| [7] |
ZDRAVKOVICH M M. Review of flow interference between two circular cylinders in various arrangements[J]. Journal of Fluids Engineering, 1977, 99(4): 618 DOI:10.1115/1.3448871 |
| [8] |
GU Z, SUN T. On interference between two circular cylinders in staggered arrangement at high subcritical Reynolds numbers[J]. Journal of Wind Engineering and Industrial Aerodynamics, 1999, 80(3): 287 DOI:10.1016/S0167-6105(98)00205-0 |
| [9] |
杜晓庆, 张利平, 刘庆宽, 等. 缆索承重桥并列索平均气动性能的雷诺数效应[J]. 工程力学, 2017, 34(3): 189 DU Xiaoqing, ZHANG Liping, LIU Qingkuang, et al. Reynolds number effects on mean aerodynamic forces acting on parallel cables of cable-supported bridges[J]. Engineering Mechanics, 2017, 34(3): 189 |
| [10] |
TING D S K, WANG D J, PRICE S J, et al. An experimental study on the fluid elastic forces for two staggered circular cylinders in cross-flow[J]. Journal of Fluids and Structures, 1998, 12(3): 259 DOI:10.1006/jfls.1997.0140 |
| [11] |
MAIR W A, MAULL D J. Aerodynamic behaviour of bodies in the wakes of other bodies[J]. Philosophical Transactions of the Royal Society of London A: Mathematical, Physical and Engineering Sciences, 1971, 269: 425 DOI:10.1098/rsta.1971.0042 |
| [12] |
PRICE S J. The origin and nature of the lift force on the leeward of two bluff bodies[J]. The Aeronautical Quarterly, 1976, 27(2): 154 DOI:10.1017/S0001925900007642 |
| [13] |
ALAM M M, SAKAMOTO H. Investigation of Strouhal frequencies of two staggered bluff bodies and detection of multistable flow by wavelets[J]. Journal of Fluids and Structures, 2005, 20(3): 425 DOI:10.1016/j.jfluidstructs.2004.11.003 |
| [14] |
ALAM M M, MEYER J P. Two interacting cylinders in cross flow[J]. Physical Review E, 2011, 84(5): 056304 DOI:10.1103/PhysRevE.84.056304 |
| [15] |
ZHOU Y, ALAM M M. Wake of two interacting circular cylinders: a review[J]. International Journal of Heat and Fluid Flow, 2016, 62: 510 DOI:10.1016/j.ijheatfluidflow.2016.08.008 |
| [16] |
杜晓庆, 施春林, 孙雅慧, 等. 高雷诺数下串列圆柱尾流致涡激振动的机理研究[J]. 振动工程学报, 2018, 31(4): 688 DU Xiaoqing, SHI Chunlin, SUN Yahui, et al. Mechanism of wake-induced vortex vibration of two tandem circular cylinders at a high Reynolds number[J]. Journal of Vibration Engineering, 2018, 31(4): 688 |
| [17] |
NISHIMURA H, TANIIKE Y. Aerodynamic characteristics of fluctuating forces on a circular cylinder[J]. Journal of Wind Engineering and Industrial Aerodynamics, 2001, 89(7/8): 713 |
| [18] |
CANTWELL B J, COLES D. An experimental study of entrainment and transport in the turbulent near wake of a circular cylinder[J]. Journal of Fluid Mechanics, 1983, 136: 321 DOI:10.1017/S0022112083002189 |
| [19] |
SUMNER D, RICHARDS M D, AKOSILE O O. Two staggered circular cylinders of equal diameter in cross-flow[J]. Journal of Fluids and Structures, 2005, 20(2): 255 DOI:10.1016/j.jfluidstructs.2004.10.006 |