展向剪切流作用下斜置圆柱气动特性研究

杜晓庆，费陈杰，况中华，赵 燕

(上海大学 土木工程系，上海 200072)

圆柱绕流是一个经典的流体力学问题，具有极为重要的研究意义和应用价值。在过去一个世纪内，虽然对这个问题已取得了大量研究成果［1］，但某些特定条件下的圆柱绕流问题仍需要进一步的研究。文献［2－8］对均匀流作用下的斜圆柱绕流问题进行了研究，而文献［9－14］研究了展向剪切流作用下直立圆柱的绕流特性。目前尚未见到研究展向剪切流作用下的斜置圆柱绕流的成果。

在很多实际工程中，如悬索桥主缆、斜拉桥拉索、桅杆拉索、海洋工程中的管道和锚缆等圆柱形结构，流体并非垂直作用于圆柱。并且来流也很少以均匀流的形式出现，在大气边界层和海底边界层内均存在相当大的剪切作用，局部地貌也会增加来流的不均匀性。因此，对展向剪切流作用下的斜置圆柱绕流问题进行研究，具有研究和工程意义。

研究人员通过试验方法［2－5］和计算流体力学方法(CFD)［6－8］对均匀流作用下斜置圆柱的绕流问题进行了研究。研究发现:沿斜置圆柱轴线方向的流体速度与垂直于圆柱截面方向的流体速度具有相同的量级［5，7－8］。Mastumoto 等［5］认为沿圆柱轴向方向的流体运动(或称为轴向流)会干扰并减弱卡门涡强度。

对于剪切流垂直作用于直立圆柱的绕流问题，研究者一般也是通过试验方法［9－11］和 CFD 方法［12－14］进行研究。研究发现:在展向剪切流作用下，圆柱展向会产生倾斜的胞状旋涡脱落［9，14］。

本文针对展向剪切流作用下的斜置圆柱绕流问题，采用大涡模拟(LES)方法，在雷诺数为3 900时(雷诺数根据平均来流速度和圆柱直径计算得到)，研究了来流与圆柱轴向夹角为45°的斜置圆柱在均匀流和两种不同强度展向剪切流作用下，斜置圆柱的表面风压分布、Strouhal数、气动力功率谱和相关性等气动性能，并探讨了剪切流的作用机理。

1 计算模型

1.1 坐标系及测点布置

坐标系和圆柱倾角定义如图1所示，图中α为圆柱倾角，D为圆柱直径，L为圆柱长度，U为来流速度，Uo为来流在垂直于圆柱轴线方向的速度分量，Ua为来流在沿圆柱轴线方向的速度分量。图中XYZ为整体坐标，xyz为圆柱轴截面局部坐标，其中Y轴和y轴都垂直于平面向内。

图1 坐标系及测压截面Fig.1 Coordinates and pressure measuring sections

在本文中，圆柱长度L为20D，斜置圆柱的倾角α为45°。为了与直立圆柱进行比较，本文还考虑了α为0°的直立圆柱在均匀流和剪切流作用下情况。为了降低计算工作量，计算雷诺数取为3 900(雷诺数根据平均来流速度U和圆柱直径D计算得到)。无量纲的时间步长为 Δt*=0.02(Δt*=UΔt/D)。

沿圆柱展向(即圆柱轴线方向)布置了4个测压截面(见图1)，每个截面沿圆柱周向均匀布置32个压力测点，在计算结果稳定可靠后开始采集数据。在本文中，除特别注明外，文中结果均来自B截面计算数据。

1.2 剪切流定义

展向线性剪切流指来流沿圆柱高度是线性均匀变化的如图2所示，风速沿着整体坐标系的Z轴方向线性增大。剪切流的强度可用剪切系数k来表示，k=AD/(其中 A= U/Z，D 为圆柱直径，为来流平均速度，见图2所示)。在本文中，剪切系数k取为0(即均匀流)、0.02 和 0.05。

图2 展向剪切流定义Fig.2 Definition of span-wise shear flow

1.3 气动力定义

圆柱表面的气动力系数CD、CL定义如下，方向如图3所示:

式中:ρ为空气密度;Fx和Fy为作用在单位长圆柱模型上的x和y方向的气动力，Cpi为圆柱表面第i个测点的风压系数;Δθi为测点i和i+1之间交角;θcyli为第i个测点的位置角。

图3 气动力系数定义Fig.3 Definition of aerodynamic forces

2 计算方法

2.1 控制方程和亚格子模型

大涡模拟(LES)方法中，大尺度涡通过滤波后的Navier-Stokes方程直接求解，小尺度的涡则采用亚格子尺度模型(SGS)模拟。与雷诺平均法(RANS)相比，大涡模拟方法可更好地模拟流场中的湍流漩涡，能捕捉到丰富的流场脉动信息。而与直接数值模拟(DNS)相比，大涡模拟又能节省很大的计算量。本文采用Fluent软件的大涡模拟湍流模型进行数值分析。

经过滤波函数的滤波，可得到大尺度涡的不可压缩Navier-Stokes方程:

式中:为滤波后速度分量;ρ为流体密度;p—滤波后压力;ν为流体动力粘度;τij为亚格子应力张量。

亚格子应力τij采用Smagorinsky-Lilly的亚格子尺度模型，具有以下形式:

式中，νt是亚格子尺度的湍动粘度，为:

式中，Δi代表沿坐标轴i方向的网格尺寸，Cs是Smagorinsky常数，本文取0.1。

2.2 计算域和边界条件

计算域为:长度×宽度×展向=30 D×20 D×20 D，详见图4。网格类型为六面体网格，圆柱周向均匀布置128个点，沿圆柱展向间距为0.1D，模型最小网格4 ×10－4D，最大网格为 1.05 D，总网格数量约为200万。

计算域的边界条件如下:进口采用速度边界条件;出口采用出流边界条件;圆柱上下面采用自由滑移壁面条件;圆柱展向端部采用对称边界条件;圆柱表面采用无滑移壁面条件。

图4 平面网格及计算域Fig.4 Computational domain and grids

3 计算结果分析

3.1 计算结果验证

为了验证本文计算结果和计算方法的正确性和适用性，将本文直立圆柱的平均阻力系数和Strouhal数(St)与文献［15］的大涡模拟结果和文献［16］的风洞试验结果进行了比较，结果列于表1内。从表1可见，本文的圆柱展向长度与文献［16］的风洞试验模型相似;而三者的雷诺数相近，均在亚临界区内。从数值上看，本文的计算结果与其他研究者的结果相近，但平均阻力系数稍偏大而St数偏小。这种偏差可能是由于本文的计算条件和文献［15］的计算参数以及文献［16］的风洞试验条件的差异所造成的。Zdravkovich［1］指出圆柱绕流对圆柱的表面粗糙度、来流湍流度、雷诺数、模型长度和阻塞率等条件非常敏感。因此，从总体上看，本文的计算结果是正确的，本文所采用的计算方法是合适的。

表1 直立圆柱的计算结果比较Tab.1 Comparison of results of non-inclined cylinders

3.2 风压系数分布

图5为倾角α为45°的斜置圆柱分别在均匀流(k为0)和两种不同剪切流(k为0.02和0.05)作用下的平均风压系数和脉动风压系数(RMS)，图中横坐标为从停滞点开始的测压点角度(角度定义见图3)。为了进行比较，图中也列出了直立圆柱(α为0°)在均匀流作用下的计算结果。风压系数的计算均采用测压截面处垂直于圆柱轴线方向的来流风速分量Uo。

由图5(a)可见，对于均匀流作用下的圆柱，当倾角从0°增大为45°时，圆柱表面最大平均负压以及平均基底风压的绝对值均减小。有趣的是:对于倾角为45°的斜置圆柱，最大平均负压和平均基底风压的绝对值又随着剪切系数的增大而逐渐增大;在k为0.05的剪切流作用下，斜置圆柱表面的平均风压系数与均匀流作用下直立圆柱的表面风压系数非常接近。

对于圆柱表面的脉动风压系数(见图5(b))，在均匀流作用下，斜置圆柱表面的脉动风压系数要比直立圆柱小。这可能是由于斜置圆柱受到尾流区轴向流的影响，从而导致卡门涡脱强度降低所造成的。斜置圆柱脉动风压系数随剪切系数的变化规律与平均风压相似:随着剪切系数的增大，斜圆柱的脉动风压系数逐渐增大并接近直立圆柱的情况。

3.3 Strouhal数沿展向变化

图6为直立和斜置圆柱在剪切流作用下的Strouhal数(St)沿圆柱展向的变化规律。图中纵坐标(z/D)为展向高度与圆柱直径的比值。Strouhal数定义为:St=fD/，为来流平均风速，f为卡门涡脱频率。

图5 倾角和剪切系数对风压系数的影响Fig.5 The effects of inclined angle and shear parameters on wind pressure coefficients

图6 Strouhal数沿展向的变化Fig.6 The change of Strouhal numbers along span-wise

由图6(a)可得，对于直立圆柱体，在剪切流作用下，St沿着圆柱展向发生明显的变化，并且可以明显观察到在一定长度范围内存在St锁定的现象，即沿着展向形成所谓的胞状涡，这与文献［9，13］中所得结论基本一致。然后，对于剪切流作用下的斜置圆柱体，除了在圆柱上端外，St沿着展向没有观察到明显的变化，见图6(b)。这说明在斜圆柱尾流区形成的轴向流会破坏胞状涡的形成。

3.4 气动力系数

图7为直立和斜置圆柱在均匀流和剪切流作用下的气动力系数时程曲线。气动力系数根据式(1)和(2)计算得到。

图7 倾角和剪切系数对气动力时程的影响Fig.7 Effects of inclined angle and shear parameter on time history of aerodynamic forces

由图7可得，均匀流作用下，直立圆柱的平均阻力系数为 1.15，与 Breuer［15］采用大涡模拟方法在相同雷诺数下得到的结果基本一致。斜置圆柱在均匀流作用下的平均阻力系数为0.98，与直立圆柱相比有所减少。在剪切流作用下，斜置圆柱的平均阻力系数随着剪切系数的增大而增大，并逐渐接近均匀流作用下直立圆柱的平均阻力系数。平均阻力系数的这种变化趋势与平均风压系数的变化规律相同。

在圆柱倾角和剪切流影响下，升力系数时程较阻力系数时程的变化更为剧烈。在均匀流作用下，斜置圆柱的升力时程曲线与直立圆柱相比，升力的波动幅度大大降小，并且升力波动的规律性也降低了。这可能是因为斜置圆柱的卡门涡脱的强度和规律性受到轴向流的影响导致的。对于斜置圆柱，随着剪切流剪切系数的增大，升力时程曲线的波动幅度逐渐增大，规律性也逐步增强。这说明剪切流可增强斜置圆柱尾流区的卡门涡脱强度。

图8为上述四种情况下升力系数的频率谱(PSD)。由图可得，在均匀流作用下，斜圆柱的卡门涡脱强度比直立圆柱明显降低。这与Matsumoto等［5］所得结论基本一致，在斜置圆柱的尾流区会形成强烈的轴向流，而轴向流会干扰尾流区的卡门涡脱，从而导致卡门涡脱的下降。而在剪切流作用下，斜置圆柱体的卡门涡脱再次得到增强，并且卡门涡的强度随剪切系数的增大而增强。

图8 升力系数功率谱Fig.8 Power spectra density of lift coefficients

对于St，均匀流作用下的直立圆柱为0.19，而斜置圆柱为0.15。对于斜置圆柱，剪切流的剪切系数基本不影响St数的大小。

3.5 气动力展向相关性

图9是圆柱展向不同截面的气动力系数之间的相关系数。各截面之间的无量纲间距的计算方法见表2，截面位置见图1。

图9 气动力系数的相关性Fig.9 Correlation coefficients of lift coefficients

从总体上看，阻力系数的展向相关性较升力系数差，倾角和剪切流对阻力系数的影响较升力系数小。从图9(a)可见，在均匀流作用下，直立圆柱的阻力系数相关性比斜置圆柱略好;而在剪切流作用下，斜置圆柱阻力系数的相关性略有增强。从图9(b)可见，在均匀流作用下，直立圆柱与斜置圆柱的升力相关性比较接近，这说明升力系数对圆柱的倾角不是敏感。但在剪切流作用下，斜置圆柱的升力相关性随着剪切系数的增大而增强，与直立圆柱的升力相关性相比有很大的提高。

表2 各截面无量纲间距Tab.2 Dimensionless distance between different sections

3.6 机理分析

图10为斜置圆柱在均匀流和剪切流作用下的流线图。由图10(a)可见，在均匀流作用时，流线在接近圆柱迎风侧时，会沿着圆柱轴线方向移动;然后再向圆柱横截面方向弯曲并绕经圆柱;流线绕过圆柱并沿圆柱轴线方向移动一段距离后分为两部分，一部分继续沿着轴向流动，另一部分则沿着与来流方向相近的方向向圆柱后移动。这说明:在均匀流作用下，倾角为45°的斜置圆柱尾流区会产生强烈的轴向流。这一现象与Zhao等［8］在雷诺数为1 000时采用DNS计算得到的结果相似。斜置圆柱尾流区的轴向流会干扰卡门涡脱，降低卡门涡脱强度(见图8(b))，从而导致升力时程波动幅度的大幅减小(见图7(b))。这些结论与Matsumoto等［5］的研究结论相一致。由图10(b)和图10(c)可见，在剪切流作用下，随着剪切流强度的增大，斜置圆柱尾流区的轴向流强度有逐渐减弱的趋势，并且当剪切流系数达到0.05时，圆柱尾流区的轴向流现象基本消失。

图10 斜置圆柱流线图Fig.10 The oblique cylindrical streamline

为了进一步分析剪切流对斜置圆柱的作用机理。图11列出了直立圆柱在均匀流和剪切流作用下的流线图。由图11(a)可见，在均匀流作用下，直立圆柱尾流区的流线几乎与来流方向相同，流线绕过圆柱后沿着与来流方向相近的方向向圆柱后移动。而在剪切流作用下，见图11(b)，流线绕过圆柱后分为两部分，一部分沿着与来流方向相近的方向向圆柱后移动;另一部分则沿着圆柱展向向上流动，流动一定距离后再沿着与来流相近的方向向圆柱后侧移动。这说明在剪切流作用下，在圆柱后侧也会形成轴向流，但剪切流引起的轴向流的流动方向与图10(a)中斜置圆柱的轴向流方向相反。

从以上的分析可以推断:对于展向剪切流作用下的斜置圆柱，由于剪切流引起的轴向流与斜置圆柱尾流区的剪切流方向相反，随着剪切流强度的增大，剪切流形成的轴向流可减弱或甚至抑制斜置圆柱的轴向流。这也解释了斜置圆柱的卡门涡脱强度会随着剪切流剪切系数的增大而增强的原因(见图8):在均匀流作用下，斜置圆柱尾流区会产生沿圆柱轴向流动的轴向流，轴向流会干扰并减弱卡门涡脱的强度;但展向剪切流则会破坏轴向流在尾流区的形成，消除了导致卡门涡脱减弱的因素，使得斜置圆柱的卡门涡脱的强度得到增强。

图11 直立圆柱流线图Fig.11 The vertical cylindrical streamline

4 结论

本文采用大涡模拟方法，对展向剪切流作用下斜置圆柱的气动性能和流场特征进行了研究，得到以下结论:

(1)在均匀流作用下，斜置圆柱尾流区会产生强烈的轴向流，而轴向流的出现会减弱卡门涡脱强度，并导致圆柱表面的脉动风压系数降低、圆柱平均阻力系数减小、升力时程的波动幅度减弱;

(2)随着展向剪切流强度的增大，斜置圆柱的平均阻力系数逐渐增大，升力时程的波动幅度也会增强，但剪切流对Strouhal数的影响很小;

(3)在展向剪切流作用下，斜置圆柱升力的展向相关性随着剪切系数的增大而增强，剪切流对阻力的展向相关性则影响较小;

(4)对于斜置圆柱，展向剪切流会减弱或抑制轴向流在圆柱尾流区的形成，并且当剪切流的剪切系数足够大时，可完全抑制轴向流的出现。从而减弱或消除了抑制卡门涡脱的因素，使斜置圆柱的卡门涡脱强度得到恢复。

［1］Zdravkovich M M.Flow around circular cylinders［M］.Oxford:Oxford University Press，1997.

［2］Ramberg SE.The effects of yaw and finite length upon the vortex wakes of stationary and vibrating circular cylinders［J］.Journal of Fluid Mechanics，1983，128:81 －107.

［3］杜晓庆，顾明.斜置拉索表面脉动风力特性研究［J］.振动与冲击，2012，31(4):139－144.DU Xiao-qing，GU Ming.Characteristics of fluctuating forces on inclined stay cable［J］.Journal of Vibration and Shock，2012，31(4):139 －144.

［4］杜晓庆，顾明.临界雷诺数下斜拉桥拉索的平均风压和风力特性［J］.空气动力学学报，2010，28(6):639－644.DU Xiao-qing，GU Ming.Wind pressure distribution and aerodynamic characteristic of stay cable in the critical Reynolds number regime［J］.Acta Aerodynamica Sinica，2010，28(6):639 －644.

［5］ Matsumoto M，Yagi T，Hatsuda H，et al.Dry galloping characteristics and its mechanism of inclined/yawed cables［J］. Journal of Wind Engineering and Industrial Aerodynamics，2010，98(6－7):317－327.

［6］Marshall J S.Wake dynamics of a yawed cylinder［J］.ASME Journal of Fluids Engineering，2003，125(1):97－103.

［7］Yeo D H，Jones N P.Investigation on 3－D characteristics of flow around a yawed and inclined circular cylinder［J］.Journal of Wind Engineering and Industrial Aerodynamics，2008，96(10－11):1947－1960.

［8］Zhao M，Cheng L，Zhou T M.Direct numerical simulation of three-dimensional flow past a yawed circular cylinder of infinite length ［J］.Journal of Fluids and Structures，2009，25(5):831－847.

［9］Griffin O M.Vortex shedding from bluff bodies in a shear flow:A review ［J］.ASME Journal of Fluids Engineering，1985，107(3):298－306.

［10］Woo H G C，Cermak JE，Peterka JA.Secondary flows and vortex formation around a circular cylinder in constant-shear flow ［J］.Journal of Fluid Mechanics，1989，204:523－542.

［11］ Kappler M，Rodi W，Szepessy S，et al.Experiments on the flow past long circular cylinders in a shear flow ［J］.Experiment in Fluids，2005，38(3):269－284.

［12］ Xu Y，Dalton C.Computation of force on a cylinder in a shear flow［J］.Journal of Fluids and Structures，2001，15(7):941－954.

［13］Mukhopadhyay A，Venugopal P，Vanka SP.Oblique vortex shedding from a circular cylinder in linear shear flow ［J］.Computer＆ Fluids，2002，31(1):1－24.

［14］ Parnaudeau P， Heitz D， Lamballais E，et al. Direct numerical simulations of vortex shedding behind cylinders with spanwise linear nonuniformity ［J］. Journal of Turbulence，2007，8(13):1－13.

［15］Breuer M.Numerical and modeling influences on large eddy simulations for the flow past a circular cylinder ［J］.International Journal of Heat and Fluid Flow，1998，19(5):512－521.

［16］ Cardell G S.Flow past a circular cylinder with a permeable splitter plate［D］，Ph.D.Thesis，Graduate Aeronautical Laboratory，California Institute of Technology，1993.