基于时间解析PIV的圆柱绕流尾迹特性研究

王 勇, 郝南松, 耿子海, 王万波

(1. 中国空气动力研究与发展中心 空气动力学国家重点实验室, 四川 绵阳 621000; 2. 中国空气动力研究与发展中心 气动噪声控制重点实验室, 四川 绵阳 621000)

0 引 言

在实际工程应用中，如航空中的飞机翼间支柱、海洋开发中的深海钻井平台及其钻杆、风工程中的桩基码头和大桥桥墩等，以圆柱绕流为代表的钝物体绕流是一个广泛存在的问题。由于圆柱绕流尾迹中常常伴随着复杂的流动现象，如边界层分离、旋涡形成与脱落及耗散等[1-2]，开展圆柱体绕流尾迹特性研究具有重要的实际意义。

对圆柱绕流尾迹的研究，有理论分析、实验研究和数值模拟3种方法。由于理论简化受限、实验耗费较大，现阶段的研究以数值模拟为主、实验研究和理论分析为辅。热线风速仪[3-4](Hot-Wire Anemometer, HWA) 和粒子图像测速[5-6](Particle Image Velocimetry, PIV)是常用于圆柱绕流研究的实验设备。热线风速仪价格便宜、频率响应高(达20kHz以上)，但只能进行单点测量且与流场直接接触，会对绕流尾迹产生显著影响；另一方面，PIV能够非接触地获取较大范围内整个流场的瞬时速度信息。因此，PIV是对圆柱绕流尾迹进行定性的流场显示和定量的流场测量与分析的理想工具。

Julio[7]利用采样频率fPIV=25Hz的数字式PIV对雷诺数Re=875和769的圆柱绕流进行了测量，讨论了垂直于圆柱轴线的平面和经圆柱轴线沿流向平面上的瞬时速度场和涡量场。Krothapalli等[8]采用fPIV=10Hz的PIV对Re=(0.85～1.5)×105的圆柱近尾流区域进行了测量，实验结果给出了瞬时速度场和瞬时涡量场。Braza等[9]和Perrin等[10]采用fPIV=4Hz的2D2C PIV和2D3C PIV研究了Re=1.4×105的圆柱绕流的湍流特性，并通过相位平均分析展示了绕流尾迹中旋涡交替脱落和发展的演化过程。Sung和Yoo[11]在水洞中采用fPIV=125Hz的PIV对Re=360的圆柱绕流进行了测量，提出了一种从二维PIV数据中获取三维相位平均数据的方法。Konstantinidis等[12]采用PIV在水洞中对Re=2160的自然圆柱绕流和强迫圆柱绕流进行了测量，对比分析了2种情况下的平均速度场、雷诺应力和相位平均速度场。张玮等[13]采用fPIV=10Hz的低频PIV对雷诺数Re=500的圆柱绕流尾迹进行了测量，给出了圆柱两侧各1倍直径、圆柱后方2.2倍直径所围成区域内的一系列瞬态速度矢量场、涡量场和流线图，并着重分析了1个拟周期内旋涡的发展演化过程。张孝棣等[14]采用fPIV=8Hz的PIV系统对Re=104～105范围内的圆柱绕流尾迹进行了测量，主要分析了绕流尾迹的瞬时速度场、涡量场和流线图。涂程旭等[15]采用fPIV=15Hz的PIV对Re=630、800和950等3种雷诺数的圆柱绕流场进行了实验研究，获得了速度场、涡量场以及涡脱落的时空演化规律。

尽管基于PIV的圆柱绕流研究已经取得了很多成果，并为理论分析和数值模拟提供了基础数据，但现有的研究存在2个方面的不足。首先，PIV的采样频率普遍偏低，只有几十赫兹，而高雷诺数下圆柱绕流的旋涡脱落频率往往达成百上千赫兹，导致PIV采样的奈奎斯特(Nyquist)频率难以超过旋涡脱落频率。这就使得现有研究通常只能给出流场的速度场、涡量场和流线图，而难以对流场进行深入的统计分析(如频谱分析等)，对圆柱绕流旋涡脱落频率的分析需要借助于数值模拟或者升力系数、阻力系数等直接与压力有关的测量。其次，圆柱绕流的相位平均分析中，需要引入1个额外的参考信号以区分不同时刻流场的相位，如Sung和Yoo[11]、Wlezien和Way[16]使用选定点或热线测量的速度信号，Braza等[9]、Perrin等[10]、Cantwell和Coles[17]使用圆柱表面的压强信号。但是，参考信号和尾流速度信号之间随机产生的相位抖动(Phase jitter)会使得相位平均结果高估湍流运动并使卡曼涡的结构更光滑[18]。

本文中使用采样率达fPIV=1000Hz的时间解析PIV系统对圆柱后方7.5倍直径、圆柱两侧各3.3倍直径所围成范围内的绕流尾迹在雷诺数Re=2.74×104下的非定常流场进行了测量。针对PIV获得的速度场数据，进行流场和频谱特性分析，探讨了圆柱绕流尾迹中的涡量场、旋转强度和旋涡脱落的频率特性等。同时，提出了一种简单有效的基于速度场之间相关性的相位平均分析方法，系统分析了圆柱上下两侧旋涡交替生成、脱落、发展并耗散的完整演化过程。

尽管关于圆柱尾流研究已经有比较丰富的研究成果，但本文在实验测量、数据处理与分析方面有一些新的创新点：(1) 借助时间解析PIV高时间分辨率和高空间分辨率的优点，对圆柱绕流尾迹进行了精细测量。高时间分辨率测量使得可以直接对PIV速度场数据进行频谱分析而得到旋涡脱落的频率特性。另一方面，全流场、高空间分辨率测量使得我们可以得到一个较大平面上的非定常流场的频率特性(如斯特劳哈尔数云图)，而不只是某几个点的频率特性。(2) 本文提出了基于速度相关性的相位平均方法来处理非定常结果，该方法简单有效，能有效地捕捉圆柱尾流中旋涡交替生成、脱落、发展并耗散的时空演化过程，在非定常流场测量(如柱翼组合模型和空腔模型的流场)方面具有普遍推广意义。

1 实验装置及方法

实验在中国空气动力研究与发展中心的0.55m×0.4m声学风洞开口试验段中进行。该风洞是一座单回流式低速风洞，具有开/闭口2个试验段，其中开口试验段长1.5m、宽0.55m、高0.4m，最大风速达100m/s。

速度场测量采用TSI公司的时间解析PIV系统，照明光源为Photonics公司的DM30单腔激光器(能量2×15mJ@1000Hz、双脉冲工作方式)，通过片光组件以水平方向照射圆柱轴向中央截面流场，片光厚度约1mm。使用Photran公司的SA4相机(配备Nikon 50/1.2镜头，分辨率1024pixel×1024pixel，最大满帧率3600帧/s)以垂直片光方向对片光照明区域进行拍摄并采集PIV图像。实验中设置粒子图像对的采样频率为1000Hz，2束激光的时间间隔为40μs。激光与相机之间的同步控制采用BNC公司的555时序控制器。PIV试验过程中采用DF-1500型烟雾发生器产生DEHS油雾粒子对流场进行示踪。

对PIV采集的粒子图像采用德国PIVTEC公司开发的PIVview2C软件对任意相邻的粒子图像对做互相关分析以获取所拍摄平面区域的速度场。处理过程中采用标准FFT相关算法、Nyquist频率滤波和9点最小二乘高斯拟合峰值算法等确定相关性，采用多次相关(Multiple correlation)和多网格查询(Multi-grid interrogation)方法提高分析精度，初始查询窗为32×32、最终查询窗为16×16，重叠因子50%。

实验使用的实心铝制圆柱模型(展长L=550mm，直径D=20mm)一端固定连接于支撑系统的顶部连接板，另一端用铁丝拉紧以消除模型振动。实验中，圆柱表面喷涂黑色哑光漆以消除强激光在模型表面反光而对PIV图像产生影响。整个实验的设置示意图如图1所示。

实验中，自由来流风速U∞=20m/s，基于圆柱直径的雷诺数为Re=ρU∞D/μ=2.74×104，此时圆柱绕流处于亚临界状态，边界层为层流分离而尾流已转变为湍流涡街。PIV采样时间为2.7s，共获取2700幅粒子图像对。

图1 实验设置示意图(俯视图)

(1)

其中，cov(Xk,Xl)为Xk和Xl的协方差，var(Xk)和var(Xl)分别为Xk和Xl的方差。

2 实验结果及分析

2.1 平均流场特性分析

图2给出了根据PIV测量结果计算得到的平均速度场的流线图，从图中可以看出，在圆柱壁面正后方约1.4D的区域内存在2个与圆柱半径同尺度的反向旋转的旋涡，旋涡中心位于0.8D的位置附近，圆柱上下两侧交替脱落的旋涡进入尾流场后形成了低速回流区。

图2 平均速度场的流线图

图3给出了平均速度场的归一化涡量图以及根据λci[19]准则确定的旋转强度(Swirling strength)图。从图中可以看出：(1) 涡量模值较大的区域主要位于圆柱的上、下两侧，“带状”分布的涡量沿圆柱后方回流区的边缘向下游扩散、向旋涡区聚集并逐渐减弱。(2) 圆柱上、下两侧的涡量强度相同但符号相反，顺时针旋转的旋涡对应于负的涡量，而逆时针旋转的旋涡对应于正的涡量。(3) 旋转强度沿着圆柱后方回流区的边缘逐渐增强，在0.8D的位置附近旋转强度达最大值，随后旋转强度迅速减弱。(4) 涡量和旋转强度均表明，回流中心区域(0.8D的位置附近)是流动结构变化最剧烈的区域。

图3 平均速度场的归一化涡量图和旋转强度图

2.2 脉动流场特性分析

(2)

图4 脉动速度场的湍动能和归一化雷诺剪切应力

2.3 频谱特性分析

由于实验中采用的时间解析PIV采样频率达1000Hz，大于圆柱尾涡脱落频率的2倍，因此可以对PIV速度场数据进行频谱分析。采用pWelch方法计算每个PIV网格点上速度场的功率谱密度(Power Spectra Density, PSD)。pWelch方法中引入Hamming窗以减小信号分析的“频谱泄露”，分块大小为1024，重叠因子为50%，频率分辨率为0.98Hz。

图5给出了圆柱绕流速度场的斯特劳哈尔数St=fD/U∞云图，其中f为功率谱密度的最大值对应的频率，即旋涡脱落频率。从图中可以看出，除了视场的上下边缘和前缘外，圆柱绕流尾迹中的斯特劳哈尔数稳定在0.2左右，这与文献[20-26]中报道的亚临界雷诺数条件下圆柱绕流的St≈0.20～0.21是一致的。

图5 斯特劳哈尔数云图

2.4 相位平均分析

图7和8分别给出了根据速度场之间的相关性确定的相位平均流线图和涡量图，从图中可以看到：(1) 圆柱上、下两侧不断地有旋涡交替地生成，从圆柱表面脱落向下游发展并随着能量的衰竭而耗散，形成较规则的Karman涡街脱落。(2) 在t=0T(T为圆柱尾涡脱落的周期)时刻，圆柱正后方约1.3D处存在1个顺时针旋转的大尺度旋涡，在旋涡附近存在负的涡量聚集区。在紧邻圆柱表面的右下方，1个逆时针旋转的小旋涡开始形成。(3) 在t=1/5T时刻，顺时针旋转的大尺度旋涡沿流向明显拉伸并向四周扩散，其涡量降低并开始衰竭。圆柱下方的小旋涡进一步长大、涡量逐渐聚集并增大到最大值。(4) 在t=2/5T时刻，圆柱右下方逆时针旋转的旋涡已经发展得非常饱满，成为大尺度旋涡，其中心位于圆柱正后方约1.0D的位置，正的涡量聚集区进一步发展。圆柱右上方t=0T时刻生成的大尺度旋涡已经扩散并向下游移动，涡量进一步降低、聚集区逐渐消失。此时，与圆柱表面相邻的右上方，1个新的旋涡结构开始酝酿形成。(5) 在t=3/5T时刻，逆时针旋转的旋涡随着向下游流动，旋涡中心后移、涡量降低，开始进入被主流拉伸的初始阶段。圆柱右上方形成的小尺度旋涡继续长大。(6) 在t=4/5T时刻，圆柱上方顺时针旋转的旋涡已经发展得比较饱满、其中心位于圆柱正后方约0.9D的位置，负的涡量聚集区进一步扩大。t=2/5T时刻发展得非常饱满的大尺度旋涡已经扩散并向下游移动，一个新的逆时针旋转的旋涡结构开始在紧邻圆柱的右下方酝酿形成。(7)t=T时刻的圆柱尾流特征，不论是流线还是涡量分布，都与t=0T时刻的尾流特征完全一致，圆柱绕流的旋涡演化过程完成1个完整的周期。

图6 速度场之间的相关性分析

图7 相位平均流线图

图8 相位平均涡量场

3 结 论

本文使用采样频率为fPIV=1000Hz的时间解析PIV系统对圆柱后方7.5倍直径、圆柱两侧各3.3倍直径所围成范围内的绕流尾迹在雷诺数Re=2.74×104下的非定常流场进行了细致地测量。PIV采样频率高于旋涡的变化频率，一方面使得可以直接利用速度场数据分析圆柱绕流尾迹中的旋涡脱落频率，另一方面由于每个涡脱落周期能捕获到若干幅PIV图像从而使得基于速度场之间相关性的相位平均分析成为可能。实验研究结果表明：

(1) 在圆柱后方1.4D的范围内存在1个低速回流区，其中心0.8D位置附近是流动结构变化最剧烈的区域；

(2) 圆柱后方1.9D位置附近是上/下两侧脱落旋涡交汇、耦合的区域，湍流强度和雷诺剪切应力较大，湍流脉动最强；

(3) 每个PIV网格点上速度场的频谱分析表明，圆柱绕流尾迹中的旋涡脱落频率对应的斯特劳哈尔数稳定在0.2左右；

(4) PIV速度场之间相关性系数的频谱分析表明，功率谱密度的峰值频率与圆柱后方的旋涡脱落频率是完全一致的。基于速度场之间相关性的相位平均分析方法简单有效，可以准确地识别绕流尾迹中旋涡交替生成、脱落、发展并耗散的时空演化过程，并且该方法可以很容易地推广应用到其它非定常流场的测量和分析中。

[1]Von Karman T H, Rubach H. On the mechanism of resistance in fluid[J]. Physikalische Zeitschrift, 1912, 13(5): 351-358.

[2]Williamson C H K. Vortex shedding in the cylinder wake[J]. Annu Rev Fluid Mech, 1996, 28: 477-539.

[3]Raffel M, Willert C, Wereley S, et al. Particle image velocimetry: a practical guide[M]. 2nd ed. Berlin Heidelberg: Springer-Verlag, 2007.

[4]Adrian R J. Twenty years of particle image velocimetry[J]. Exp Fluids, 2005, 39(2): 159-169.

[5]Comte-Bellot G. Hot-wire anemometry[J]. Annu Rev Fluid Mech, 1976, 8: 209-231.

[6]Perry A E. Hot-wire anemometry[M]. Oxford: Clarendon Press, 1982.

[7]Julio S. An investigation of the near wake of a circular cylinder using a video-based digital cross-correlation particle image velocimetry technique[J]. Experimental Thermal and Fluid Science, 1996, 12(2): 221-233.

[8]Krothapalli A, Shih C, Lourenco L. The near wake of a circular cylinder at 0.3

[9]Braza M, Perrin R, Hoarau Y. Turbulence properties in the cylinder wake at high Reynolds numbers[J]. Journal of Fluids and Structures, 2006, (22): 757-771.

[10]Perrin R, Cid E, Cazin S, et al. Phase-averaged measurements of the turbulence properties in the near wake of a circular cylinder at high Reynolds number by 2C-PIV and 3C-PIV[J]. Exp Fluids, 2007, (42): 93-109.

[11]Sung J, Yoo J Y. Three-dimensional phase averaging of time-resolved PIV measurement data[J]. Meas Sci Technol, 2001, (12): 655-662.

[12]Konstantinidis E, Balabani S, Yianneskis M. Conditional averaging of PIV plane wake data using a cross-correlation approach[J]. Exp Fluids, 2005, (39): 38-47.

[13]张玮, 王元, 徐忠, 等. 圆柱绕流涡系演变的DPIV测试[J]. 空气动力学学报, 2002, 20(4): 379-387.

Zhang W, Wang Y, Xu Z, et al. Experimental investigation of vortex evlovement around a circular cylinder by digital particle image velocimetry (DPIV)[J]. Acta Aerodynamica Sinca, 2002, 20(4): 379-387.

[14]张孝棣, 蒋甲利, 贾元胜, 等. 圆柱体绕流尾迹的PIV测量[J]. 实验流体力学, 2005, 19(2): 74-78.

Zhang X D, Jiang J L, Jia Y S, et al. Measurement of cylinder’s wake by PIV[J]. Journal of Experiments in Fluid Mechanics, 2005, 19(2): 74-78.

[15]涂程旭, 王昊利, 林建忠. 圆柱绕流的流场特性及涡脱落规律研究[J]. 中国计量学院学报, 2008, 19(2): 98-102.

Tu C X, Wang H L, Lin J Z. Experimental research on the flow characteristics and vortex shedding in the flow around a circular cylinder[J]. Journal of China Jiliang University, 2008, 19(2): 98-102.

[16]Wlezien R W, Way J L. Techniques for the experimental investigation of the near wake of a circular cylinder[J]. AIAA J, 1979, 17(6): 563-570.

[17]Cantwell B, Coles D. An experimental study of entrainment and transport in the turbulent near wake of a circular cylinder[J]. J Fluid Mech, 1983, 136: 321-374.

[18] Perrin R, Braza M, Cid E, et al. Coherent and turbulent process analysis in the flow past a circular cylinder at high Reynolds number[J]. Journal of Fluids and Structures, 2008, 24(8): 1313-1325.

[19]Zhou J, Adrian R J, Balachandar S, et al. Mechanisms for generating coherent packets of hairpin vortices in channel flow[J]. J Fluid Mech, 1999, 387: 353-396.

[20]Kravchenko A G, Moin P. Numerical studies of flow over a circular cylinder atRe=3900[J]. Physics of Fluids, 2000, 12(2): 403-417.

[21] Parnaudeau P, Carlier J, Heitz D, et al. Experimental and numerical studies of the flow over a circular cylinder at Reynolds number 3900[J]. Physics of Fluids, 2008, 20(8): 261-271.

[22] Schlichting H. Boundary layer theory[M]. New York: McFraw-hill Book Company, 1968.

[23] Franke J, Frank W. Large eddy simulation of the flow past a circular cylinder atRe=3900[J]. Journal of Wind Engineering and Industrial Aerodynamics, 2002, 90(10): 1191-1206.

[24]Breuer M. A challenging test case for large eddy simulation: high Reynolds number circular cylinder flow[J]. Int J Heat Fluid Flow, 2000, 21: 648-654.

[25] Dong S, Karniadakis E, Ekmekci A, et al. A combined direct numberical simulation-partical image velocimetry study of the turbulent near wake[J]. Journal of Fluid Mechanics, 2006, 569(12): 185-207.

[26]Fey U, Konig M, Eckelmann H. A new Strouhal-Reynolds number relationship for the circular cylinder in the range 47