杜曉慶， 張盛華， 王玉梁， 孫雅慧

(1. 上海大學(xué) 土木工程系，上海 200072；2. 上海大學(xué) 風(fēng)工程和氣動(dòng)控制研究中心，上海 200072)

大長(zhǎng)細(xì)比圓柱型結(jié)構(gòu)在實(shí)際工程中有大量的應(yīng)用，如橋梁并列索、超高層建筑、多分裂導(dǎo)線(xiàn)、冷卻塔、煙囪等，并且多以柱群形式出現(xiàn).多個(gè)圓柱結(jié)構(gòu)間存在強(qiáng)烈的氣動(dòng)干擾，下游圓柱在橫風(fēng)向升力作用下常發(fā)生尾流激振現(xiàn)象，導(dǎo)致結(jié)構(gòu)發(fā)生損傷或破壞[1-4]，澄清雙圓柱結(jié)構(gòu)下游圓柱的升力機(jī)理有助于理解尾流激振的發(fā)生機(jī)制.

已有研究發(fā)現(xiàn)[5-7]，對(duì)于小間距(P/D=1.2～3.5，其中P為圓心間距，D為圓柱直徑)雙圓柱，下游圓柱在小風(fēng)向角(β≈8°)下受到很大的平均升力，稱(chēng)之為內(nèi)側(cè)升力；對(duì)于中等間距雙圓柱(P/D=3.0～5.0)，下游圓柱在β=20°左右同樣受到較大的平均升力，稱(chēng)之為外側(cè)升力.

對(duì)于內(nèi)側(cè)升力，低雷諾數(shù)下的流跡顯示試驗(yàn)表明[6,8]：兩圓柱間的加速間隙流、下游圓柱的停滯點(diǎn)偏移以及上游圓柱的剪切流在下游圓柱表面發(fā)生再附是升力出現(xiàn)的三個(gè)原因.杜曉慶等[9]進(jìn)行了基于高雷諾數(shù)下的大渦模擬研究，對(duì)下游圓柱平均升力的出現(xiàn)給出了不同的解釋.

對(duì)于外側(cè)升力，Zdravkovich[7]和Ting等[10]提出，兩圓柱間隙氣流速度的增大使下游圓柱間隙側(cè)表面壓力減小，導(dǎo)致平均升力的出現(xiàn)，文獻(xiàn)[10-12]中也認(rèn)為下游圓柱風(fēng)壓停滯點(diǎn)的偏移對(duì)平均升力有貢獻(xiàn).Alam等[13-14]提出，受上游圓柱尾流旋渦的作用，下游圓柱的旋渦加速向下游流動(dòng)，這也會(huì)引起下游圓柱間隙側(cè)表面壓力的減小，進(jìn)而導(dǎo)致平均升力的出現(xiàn).

對(duì)以往文獻(xiàn)的分析可知，由于影響雙圓柱的因素眾多，干擾機(jī)理復(fù)雜，研究者對(duì)下游圓柱的升力機(jī)理尚未有統(tǒng)一的認(rèn)識(shí).值得注意的是，以往的機(jī)理解釋都是基于較低雷諾數(shù)下的流跡顯示試驗(yàn)得出的，但雙圓柱繞流有強(qiáng)烈的雷諾數(shù)效應(yīng)[3,15]，因而在高雷諾數(shù)下很可能存在不同的干擾機(jī)理.此外，以往文獻(xiàn)主要集中于對(duì)雙圓柱平均升力的研究，受到測(cè)試條件的限制，對(duì)脈動(dòng)升力及其發(fā)生機(jī)理的研究很少.

本文中采用大渦模擬(LES)方法，在高雷諾數(shù)(Re=1.4×105)下，研究了不同風(fēng)向角時(shí)中等間距(P/D=4)雙圓柱的氣動(dòng)力特性和流場(chǎng)特性，重點(diǎn)討論了下游圓柱的氣動(dòng)升力與流場(chǎng)結(jié)構(gòu)之間的內(nèi)在關(guān)系，提出了雙圓柱干擾流態(tài)隨風(fēng)向角的變化規(guī)律，進(jìn)一步澄清下游圓柱的平均升力和脈動(dòng)升力的流場(chǎng)機(jī)理.

1 數(shù)值方法和計(jì)算模型

1.1 控制方程和亞格子尺度模型

在大渦模擬方法中，大尺度的旋渦由濾波后的Navier-Stokes方程直接求解，而小尺度的旋渦則采用亞格子尺度(SGS)模型模擬.與雷諾平均法(RANS)相比，大渦模擬方法可更好地模擬流場(chǎng)中的湍流旋渦，能捕捉到更豐富的流場(chǎng)脈動(dòng)信息.因此，為了準(zhǔn)確地模擬錯(cuò)列雙圓柱周?chē)牧鲌?chǎng)特性，采用大渦模擬方法進(jìn)行研究.經(jīng)過(guò)濾波函數(shù)的濾波，可得到大尺度旋渦的不可壓縮Navier-Stokes方程，如下所示：

(1)

(2)

(3)

(4)

式中:μt為亞格子尺度的湍動(dòng)黏度;δij為單位張量.μt的計(jì)算式如下所示：

(5)

(6)

(7)

式中:Δi為沿坐標(biāo)軸i方向的網(wǎng)格尺寸；CS為Smagorinsky常數(shù)，取0.1.

1.2 計(jì)算模型和計(jì)算參數(shù)

錯(cuò)列雙圓柱的布置形式如圖1所示.圖1中:CD1、CL1分別為上游圓柱的阻力和升力；CD2、CL2分別為下游圓柱的阻力和升力；θ為圓柱表面風(fēng)壓測(cè)點(diǎn)位置的角度，以圓柱迎風(fēng)點(diǎn)為零點(diǎn)，順時(shí)針旋轉(zhuǎn)為正方向.兩個(gè)圓柱的直徑均為D，圓柱間距P為4D，風(fēng)向角β為自由來(lái)流與圓柱中心連線(xiàn)的夾角.為了掌握不同風(fēng)向角對(duì)圓柱周?chē)鲌?chǎng)特性和圓柱氣動(dòng)力特性的影響規(guī)律，計(jì)算的風(fēng)向角β分別為0°、5°、10°、15°、20°、25°、30°、45°、60°、75°、90°.

圖1 雙圓柱計(jì)算模型

計(jì)算域網(wǎng)格采用結(jié)構(gòu)化網(wǎng)格、O型域，徑向長(zhǎng)度為46D，展向長(zhǎng)度為2D，均分為20層網(wǎng)格，計(jì)算模型總網(wǎng)格數(shù)約為266萬(wàn)，近壁面最小網(wǎng)格厚度為1×10-4D，最大阻塞率為4.3%(β=90°).量綱一時(shí)間步Δt*為0.005(Δt*=ΔtU0/D，其中Δt為實(shí)際計(jì)算時(shí)間步，U0為來(lái)流風(fēng)速).計(jì)算域采用均勻速度入口邊界條件，自由出口邊界條件，展向兩端采用周期性邊界條件，圓柱表面采用無(wú)滑移壁面邊界條件，如圖2所示.計(jì)算流體動(dòng)力學(xué)(CFD)數(shù)值模擬基于Fluent14.0程序?qū)崿F(xiàn)，采用壓力耦合方程的半隱式(SIMPLEC)格式求解壓力速度耦合方程組，空間離散采用中心差分格式，時(shí)間離散采用二階全隱格式.

a 計(jì)算域和邊界條件

b 局部的計(jì)算網(wǎng)格

2 計(jì)算結(jié)果及分析

2.1 計(jì)算參數(shù)驗(yàn)證

在確定雙圓柱計(jì)算模型和計(jì)算參數(shù)之前，先以單圓柱為研究對(duì)象，對(duì)六種計(jì)算模型(分別命名為Case 1～6)進(jìn)行數(shù)值模擬，研究周向網(wǎng)格數(shù)量、時(shí)間步和展向長(zhǎng)度等參數(shù)對(duì)計(jì)算結(jié)果的影響[16].表1給出了單圓柱各計(jì)算模型的參數(shù)、平均阻力系數(shù)、脈動(dòng)升力系數(shù)和St數(shù)，也列出了文獻(xiàn) [17]和[18]中的風(fēng)洞試驗(yàn)結(jié)果.通過(guò)比較可以看出，模型Case 6計(jì)算得到的氣動(dòng)力系數(shù)與文獻(xiàn)的吻合情況最好，因此采用該模型參數(shù)進(jìn)行錯(cuò)列雙圓柱的大渦模擬研究.

2.2 氣動(dòng)力系數(shù)

圖3是下游圓柱平均氣動(dòng)力系數(shù)曲線(xiàn).由圖3可見(jiàn)，下游圓柱的平均升力系數(shù)和平均阻力系數(shù)隨風(fēng)向角的變化趨勢(shì)與文獻(xiàn)較為吻合.與單圓柱相比，在風(fēng)向角0°～30°內(nèi)，下游圓柱的平均阻力系數(shù)均小于單圓柱，并且隨風(fēng)向角的增大而增大.從β=30°開(kāi)始，下游圓柱的平均阻力系數(shù)趨于平緩，其值比單圓柱稍大.

表1 單圓柱計(jì)算結(jié)果對(duì)比[16]

a 阻力系數(shù)

b 升力系數(shù)

對(duì)于平均升力系數(shù)來(lái)說(shuō)，隨著風(fēng)向角的增大，下游圓柱的平均升力系數(shù)先減小后增大，最后趨于穩(wěn)定.在β=0°～45°時(shí)，下游圓柱的平均升力系數(shù)均為負(fù)值，并且在β=20°左右出現(xiàn)升力的極小值，即出現(xiàn)圓柱上側(cè)升力峰值.當(dāng)β≥45°時(shí)，平均升力系數(shù)便與單圓柱相似并接近于零.圖4是下游圓柱脈動(dòng)氣動(dòng)力系數(shù)隨風(fēng)向角的變化曲線(xiàn).由圖4可以看出，在β=0°～5°時(shí)下游圓柱的脈動(dòng)升力最大，其值達(dá)到單圓柱的2倍，并且在β=20°～25°時(shí)達(dá)到極小值.當(dāng)β>30°后，下游圓柱的脈動(dòng)升力系數(shù)開(kāi)始緩慢減小，但依然明顯高于單圓柱；下游圓柱的脈動(dòng)阻力則在β=10°左右達(dá)到最大值，并且在β=0°～30°范圍內(nèi)都明顯大于單圓柱.

圖4 下游圓柱脈動(dòng)氣動(dòng)力系數(shù)隨風(fēng)向角的變化

Fig.4 Root mean square aerodynamic coefficients as a function of incidence angles for downstream cylinder

2.3 下游圓柱表面風(fēng)壓

圖5為下游圓柱在不同風(fēng)向角下的平均風(fēng)壓系數(shù)曲線(xiàn).從圖5可以看出：當(dāng)β=0°時(shí)，整個(gè)下游圓柱表面都為負(fù)壓，并具有良好的對(duì)稱(chēng)性；當(dāng)β=5°時(shí)，下游圓柱表面負(fù)壓遠(yuǎn)小于其他工況，并且停滯點(diǎn)偏移了20°左右；在β=15°～30°時(shí)，下游圓柱表面平均風(fēng)壓呈現(xiàn)出明顯的不對(duì)稱(chēng)性，圓柱下側(cè)(間隙側(cè))表面受到的負(fù)壓絕對(duì)值明顯大于上側(cè)，而平均風(fēng)壓分布的不均勻則導(dǎo)致圓柱表面受到平均升力；當(dāng)風(fēng)向角大于45°時(shí)，下游圓柱表面平均風(fēng)壓逐漸接近單圓柱，又重新表現(xiàn)出良好的對(duì)稱(chēng)性.

圖6為下游圓柱在不同風(fēng)向角下的脈動(dòng)風(fēng)壓系數(shù)曲線(xiàn).由圖6可知：幾乎所有風(fēng)向角下下游圓柱的脈動(dòng)風(fēng)壓系數(shù)值都不同程度地大于單圓柱，這說(shuō)明受上游圓柱的影響，下游圓柱表面的脈動(dòng)風(fēng)壓發(fā)生劇烈變化；在β=5°時(shí)，下游圓柱表面的脈動(dòng)風(fēng)壓變化最為劇烈，也具有一定的對(duì)稱(chēng)性；在β=15°～30°時(shí)，下游圓柱表面的脈動(dòng)風(fēng)壓呈現(xiàn)出明顯的不對(duì)稱(chēng)性，下游圓柱上側(cè)表面的脈動(dòng)風(fēng)壓變化較為平緩，而下側(cè)(間隙側(cè))表面的脈動(dòng)風(fēng)壓變化很大，脈動(dòng)風(fēng)壓的不對(duì)稱(chēng)性也會(huì)促使下游圓柱受到平均升力的作用；在β=45°～90°時(shí)，隨著風(fēng)向角的增大，下游圓柱表面的脈動(dòng)風(fēng)壓逐漸接近于單圓柱，并且重新具有良好的對(duì)稱(chēng)性.

圖5 下游圓柱表面平均風(fēng)壓系數(shù)分布

圖6 下游圓柱表面脈動(dòng)風(fēng)壓系數(shù)分布

2.4 氣動(dòng)力功率譜

對(duì)于錯(cuò)列雙圓柱來(lái)說(shuō)，不同風(fēng)向角下的St差別不大，都在0.17左右，但功率譜(PSD)峰值卻差別很大.圖7是不同風(fēng)向角下下游圓柱的升、阻力功率譜曲線(xiàn).從圖7可以看出，下游圓柱升、阻力功率譜隨風(fēng)向角的變化與如圖4所示下游圓柱脈動(dòng)氣動(dòng)力系數(shù)的變化趨勢(shì)一致.

圖7 氣動(dòng)力系數(shù)的功率譜峰值隨風(fēng)向角的變化

Fig.7 PSD peak of aerodynamic coefficient as a function of incidence angles

當(dāng)β=0°和5°時(shí)，下游圓柱的升力功率譜峰值最大，說(shuō)明此時(shí)渦脫強(qiáng)度很大，圓柱表面出現(xiàn)很強(qiáng)的脈動(dòng)風(fēng)壓，進(jìn)而產(chǎn)生很大的脈動(dòng)升力；當(dāng)β=20°時(shí)，下游圓柱的升力功率譜峰值最小，表明此時(shí)的渦脫強(qiáng)度明顯減弱，相應(yīng)的下游圓柱表面受到的脈動(dòng)升力也隨之降低.

2.5 干擾流態(tài)分類(lèi)

通過(guò)上述分析可以發(fā)現(xiàn)，不同的風(fēng)向角范圍內(nèi)，下游圓柱的氣動(dòng)力有著很大差異.通過(guò)觀察不同風(fēng)向角時(shí)下游圓柱周?chē)牧鲌?chǎng)變化，可以總結(jié)出以下四種不同的流態(tài)(見(jiàn)圖8)：

(1) 流態(tài)1，旋渦撞擊流態(tài).當(dāng)β=0°時(shí)，上游圓柱上側(cè)脫落的旋渦撞擊到下游圓柱上側(cè)，接著上游圓柱下側(cè)脫落的旋渦撞擊到下游圓柱的下側(cè)(間隙側(cè)).上游圓柱脫落的旋渦對(duì)下游圓柱的不斷撞擊是下游圓柱表面出現(xiàn)較大脈動(dòng)風(fēng)壓和脈動(dòng)升力的主要原因.

(2) 流態(tài)2，旋渦撞擊和剪切層干擾流態(tài).當(dāng)β=5°～10°時(shí)，上游圓柱上側(cè)脫落的旋渦撞擊到下游圓柱的迎風(fēng)面，上游圓柱下側(cè)脫落的旋渦并沒(méi)有撞擊到下游圓柱，而是與下游圓柱下側(cè)的剪切層相互作用.同時(shí)存在的旋渦撞擊和剪切層干擾，導(dǎo)致下游圓柱表面出現(xiàn)很大的脈動(dòng)升力.

(3) 流態(tài)3，剪切層干擾流態(tài).當(dāng)β=15°～30°時(shí)，上游圓柱上側(cè)脫落的旋渦與下游圓柱下側(cè)(間隙側(cè))的剪切層和旋渦耦合產(chǎn)生相互作用，而上游圓柱下側(cè)脫落的旋渦在向下游運(yùn)動(dòng)的過(guò)程中由于距離下游圓柱較遠(yuǎn)，故沒(méi)有發(fā)生明顯的相互作用.

a 流態(tài)1(旋渦撞擊流態(tài)，β=0°)

b 流態(tài)2(旋渦撞擊和剪切層干擾流態(tài)，β=5°～10°)

c 流態(tài)3(剪切層干擾流態(tài)，β=15°～30°)

d 流態(tài)4(尾流干擾流態(tài)，β=45°～90°)

圖8 錯(cuò)列雙圓柱在不同風(fēng)向角下的干擾流態(tài)

Fig.8 Interaction flow patterns of two staggered circular cylinders at different incidence angels

(4) 流態(tài)4，尾流干擾流態(tài).當(dāng)β=45°～90°時(shí)，上游圓柱與下游圓柱間距較大，上游圓柱上側(cè)脫落的旋渦在向下游運(yùn)動(dòng)的過(guò)程中只有一部分尾流會(huì)與下游圓柱發(fā)生相互作用，而下側(cè)脫落的旋渦則基本不會(huì)對(duì)下游圓柱有影響.

2.6 平均升力的機(jī)理分析

為了進(jìn)一步澄清風(fēng)向角20°附近下游圓柱出現(xiàn)平均升力峰值(即文獻(xiàn)[15]中的外側(cè)升力)的流場(chǎng)機(jī)理，本節(jié)從平均流場(chǎng)和瞬態(tài)流場(chǎng)兩方面進(jìn)行分析.

2.6.1平均流場(chǎng)分析

圖9是β=20°時(shí)下游圓柱的平均風(fēng)壓場(chǎng)和平均流速.從圖9可以看出，下游圓柱下側(cè)(間隙側(cè))的負(fù)壓強(qiáng)度明顯大于上側(cè)，并且從下游圓柱的表面平均風(fēng)壓分布可知，下游圓柱的風(fēng)壓停滯點(diǎn)相對(duì)于來(lái)流向上側(cè)偏移了6°左右，這兩個(gè)因素對(duì)下游圓柱的平均升力均有一定貢獻(xiàn).此外，以往的學(xué)者還認(rèn)為[8]，加速間隙流也是導(dǎo)致下游圓柱出現(xiàn)平均升力的原因之一，但從平均風(fēng)速圖(見(jiàn)圖9b)和下游圓柱附近的局部風(fēng)速圖(見(jiàn)圖9c)可看出，兩個(gè)圓柱間隙中的風(fēng)速并未有明顯的加速現(xiàn)象，并且下游圓柱上、下側(cè)表面的最大平均風(fēng)速也未見(jiàn)明顯的差別.上述結(jié)果說(shuō)明，在該工況下間隙流并未有明顯的加速現(xiàn)象出現(xiàn)，因而對(duì)下游圓柱平均升力的貢獻(xiàn)不大.

a 平均風(fēng)壓場(chǎng)

b 平均流速

c 下游圓柱附近的平均流速

圖9 平均風(fēng)壓場(chǎng)、平均流速和下游圓柱的局部流速

Fig.9 Mean pressure coefficient field, mean velocity ratio and local velocity ratio of downstream cylinder

2.6.2瞬態(tài)流場(chǎng)分析

從第2.6.1節(jié)分析可知，下游圓柱上、下側(cè)表面的平均風(fēng)壓存在較大差異，這是下游圓柱受到較大平均升力的主要原因之一.為了進(jìn)一步澄清該現(xiàn)象的流場(chǎng)機(jī)理，圖11給出了風(fēng)向角為20°時(shí)九個(gè)瞬態(tài)時(shí)刻(T1～T9)的渦量圖及下游圓柱表面的風(fēng)壓系數(shù)分布.這九個(gè)瞬態(tài)時(shí)刻在升力系數(shù)時(shí)程曲線(xiàn)上的具體位置如圖10所示，下游圓柱的升力系數(shù)在T1和T9時(shí)刻接近極大值，在T6時(shí)刻則接近極小值.

由圖10還可知，上、下游圓柱的升力時(shí)程曲線(xiàn)均有明顯的周期性，但兩者有一定的相位差；對(duì)于下游圓柱而言，升力系數(shù)從波峰下降至波谷的時(shí)間較長(zhǎng)，而從波谷上升至波峰的時(shí)間較短.

圖10 β=20°時(shí)升力系數(shù)的時(shí)程曲線(xiàn)及瞬態(tài)時(shí)刻

由圖11可見(jiàn)，上游圓柱的上、下表面交替形成分離的剪切層，剪切層受另一側(cè)旋渦的影響而發(fā)生卷起，卷起的剪切層最終脫離圓柱并在尾流中形成方向交替變化的旋渦.隨著剪切層的變化，上游圓柱表面的風(fēng)壓也發(fā)生相應(yīng)的交替變化.上游圓柱的繞流場(chǎng)和表面風(fēng)壓的變化與單圓柱相似.

下游圓柱受到上游圓柱尾流旋渦的影響，導(dǎo)致其氣動(dòng)力特性與單圓柱有很大差別.由圖11可見(jiàn)，從上游圓柱脫落的旋渦與下游圓柱的下側(cè)剪切層(或旋渦)發(fā)生相互作用.相對(duì)來(lái)說(shuō)，上游圓柱上側(cè)脫落的旋渦與下游圓柱的相互作用更為劇烈(T1、T2、T8和T9時(shí)刻)，而上游圓柱下側(cè)脫落的旋渦對(duì)下游圓柱的影響則相對(duì)較小(T4、T5和T6時(shí)刻).

特別在T8和T9時(shí)刻，當(dāng)上游圓柱上側(cè)脫落的旋渦移動(dòng)至下游圓柱的下側(cè)剪切層附近時(shí)，由于旋渦與剪切層的旋轉(zhuǎn)方向不同，上游圓柱的上側(cè)旋渦吸引下游圓柱的下側(cè)剪切層從圓柱表面脫落而形成旋渦；隨后，這一對(duì)方向相反的旋渦加速向下游移動(dòng)(T9時(shí)刻).由于下游圓柱的下側(cè)旋渦迅速離開(kāi)圓柱，使得下側(cè)旋渦對(duì)下游圓柱上側(cè)剪切層的影響減弱，從而導(dǎo)致下游圓柱上側(cè)剪切層的卷起現(xiàn)象不明顯.在T1和T2時(shí)刻能觀察到類(lèi)似現(xiàn)象.因而，下游圓柱的上側(cè)剪切層在所有時(shí)刻中均保持較為平順的形態(tài)，上側(cè)剪切層對(duì)圓柱的影響相對(duì)較小，所以下游圓柱上側(cè)表面的負(fù)壓在所有時(shí)刻中均保持在相對(duì)較低的數(shù)值，最終導(dǎo)致下游圓柱上、下側(cè)表面的平均風(fēng)壓出現(xiàn)較大差異，這是下游圓柱受到較大平均升力作用的原因之一.

圖11 β=20°時(shí)一個(gè)周期內(nèi)的瞬時(shí)渦量及表面風(fēng)壓分布

此外，由圖11還可見(jiàn)，T4～T7時(shí)刻的下游圓柱風(fēng)壓停滯點(diǎn)明顯向圓柱上側(cè)偏移，而T1和T9時(shí)刻則向下側(cè)偏移，T2、T3和T8時(shí)刻則偏移不明顯，因而總體上下游圓柱的停滯點(diǎn)是向上偏移，這是導(dǎo)致下游圓柱受到較大平均升力的另一個(gè)原因.

3 結(jié)論

(1)隨著風(fēng)向角的增大，兩個(gè)圓柱的干擾流態(tài)依次為旋渦撞擊流態(tài)、旋渦撞擊和剪切層干擾流態(tài)、剪切層干擾流態(tài)、尾流干擾流態(tài).在不同流態(tài)下，下游圓柱的氣動(dòng)力特性有很大差異.

(2)在風(fēng)向角為0°～10°的旋渦撞擊流態(tài)以及旋渦撞擊和剪切層干擾流態(tài)時(shí)，從上游圓柱脫落的旋渦與下游圓柱發(fā)生強(qiáng)烈的撞擊，導(dǎo)致下游圓柱的脈動(dòng)升力和迎風(fēng)面的脈動(dòng)風(fēng)壓遠(yuǎn)遠(yuǎn)大于其他流態(tài)，也遠(yuǎn)大于單圓柱.

(3)在風(fēng)向角為15°～30°的剪切層干擾流態(tài)下，下游圓柱受到顯著的平均升力作用.下游圓柱風(fēng)壓停滯點(diǎn)的偏移、上游圓柱的尾流旋渦與下游圓柱間隙側(cè)剪切層和旋渦的相互作用，是下游圓柱出現(xiàn)平均升力的兩個(gè)主要原因.