汪 洋,陳 立

(四川省公路規(guī)劃勘察設(shè)計研究院有限公司，四川成都 610000)

斜拉索、電纜等鈍體結(jié)構(gòu)由于其柔度大、阻尼小易導(dǎo)致各種風(fēng)致振動問題，對于這類問題一般將結(jié)構(gòu)簡化為圓柱模型進(jìn)行研究[1-5]。圓柱繞流是流體力學(xué)中的經(jīng)典問題之一，當(dāng)粘性流體流經(jīng)鈍體時，會發(fā)生邊界層分離現(xiàn)象，進(jìn)而在尾部形成周期性渦脫，同時伴隨著鈍體兩側(cè)壁面壓力的變化，從而對結(jié)構(gòu)產(chǎn)生有規(guī)律的氣動力作用。特別對于雙圓柱繞流的情況，圓柱間的耦合干擾作用會導(dǎo)致其氣動力特性有別于單圓柱繞流，因此雙圓柱繞流問題一直是專家學(xué)者研究的熱點[6-10]。

大部分研究集中在串列圓柱繞流上。Zdravkovich認(rèn)為串列圓柱之間的流體流動形態(tài)與柱間中心距有關(guān)，當(dāng)兩圓柱中心距為1.0D～1.5D時，上游圓柱分離出的自由剪切層會越過下游圓柱；當(dāng)兩圓柱中心距為1.5D～4.0D時，剪切層會附著在下游柱體上；當(dāng)中心距大于4D時，上游和下游柱體會分別形成渦脫，其中D為圓柱直徑。Kitagawa運用大渦模擬研究了Re=22000時不同中心距下的串列圓柱，結(jié)果表明：當(dāng)中心距小于3.25D時，只有下游圓柱出現(xiàn)漩渦脫落；當(dāng)中心距大于3.25D時，上下游均發(fā)生渦脫現(xiàn)象，并且發(fā)現(xiàn)中心距在3.25D附近時，上下游的平均阻力系數(shù)和升力系數(shù)會發(fā)生跳躍性變化。Tokoro與Komatsu的理論研究與風(fēng)洞試驗表明:當(dāng)串列圓柱中心距為4.3D風(fēng)向角為15 °時，下游圓柱會出現(xiàn)明顯振動。Brika和Lanevill采用風(fēng)洞試驗的方法對中心距為10D～25D的串列雙圓柱進(jìn)行了遠(yuǎn)距失穩(wěn)現(xiàn)象的研究。也有少部分學(xué)者對并列圓柱的繞流特性開展了研究，Sumner[15]分別對小中心距(1.0D～1.2D)、中等中心距(1.2D～2.2D)、大中心距(>2.2D)的并列圓柱繞流進(jìn)行了分析。在實際工程中，鈍體的布置形式相對于流體流動方向一般都是錯列的，然而錯列雙圓柱的繞流狀態(tài)卻很少受到研究者的關(guān)注。

本文采用大渦模擬的方法對Re=3900中心距3D的錯列雙圓柱繞流進(jìn)行數(shù)值模擬，得到了不同風(fēng)向角下的升阻力系數(shù)、尾流相對速度以及風(fēng)壓系數(shù)等參數(shù)，并分析了下游圓柱升力系數(shù)出現(xiàn)極大值的原因。

1 錯列雙圓柱繞流數(shù)值模型

1.1 計算工況

風(fēng)向角α的范圍為0～24 °，角度從0 °開始每2 °設(shè)置一個工況，一共13個工況，圓柱相對中心距W=3.0D(圖1)。

圖1 圓柱布置示意

1.2 計算域與邊界條件

計算域設(shè)置為：流向長度20D，橫向長度10D，展向長度3.2D，以保證湍流尾流的充分發(fā)展。計算區(qū)域劃分為三個部分，均采用結(jié)構(gòu)化網(wǎng)格：第一部分為3D×3D×3.2D矩形近模型區(qū)域，該區(qū)域為加密網(wǎng)格區(qū)，越靠近圓柱網(wǎng)格的密度越大，第一層網(wǎng)格到壁面的距離為2×10-4m，Y+值約為0.8；第二部分為直徑9D的圓形區(qū)域，方便在調(diào)整工況時改變風(fēng)向角；第三部分為外圍計算域，適當(dāng)降低網(wǎng)格密度以節(jié)省計算機(jī)資源。形成的網(wǎng)格見圖2，模型的總網(wǎng)格數(shù)量達(dá)到182×104個。

圖2 計算域網(wǎng)格劃分

計算域左側(cè)為速度入口，流速為U0=0.57m/s；右側(cè)為壓力出口；上下前后流域邊界為對稱邊界；圓柱表面為無滑移壁面。

坐標(biāo)系選擇：雙圓柱中心連線的中點作為坐標(biāo)原點，x軸方向與來流方向一致，y軸垂直于計算域上下邊界，z軸由右手定則確定。

1.3 參數(shù)設(shè)置

流體為空氣，密度ρ=1.225 kg/m3，動力粘度μ=1.7894×10-5pa·s。亞格子尺度模型采用WALE模型，壓力速度耦合求解算法采用對大多數(shù)問題都適用的SIMLPE算法，梯度離散采用節(jié)省計算量的最小二乘法，壓力離散格式采用標(biāo)準(zhǔn)格式，動量離散格式采用具有高階格式的有界中心差分格式。瞬態(tài)求解時時間步長Δt=0.005s，以保證庫朗數(shù)小于5且每個時間步迭代收斂。

1.4 網(wǎng)格獨立性驗證

為驗證網(wǎng)格獨立性，對雷諾數(shù)Re=3900的單圓柱在同樣的網(wǎng)格參數(shù)(計算域尺寸、網(wǎng)格密度、時間步長與錯列圓柱模擬相同)條件下進(jìn)行了數(shù)值模擬，并與Lourenco[17]與Norberg[18]的風(fēng)洞試驗進(jìn)行對比。從圖3與圖4中可見時均速度與表面壓力系數(shù)的曲線與試驗曲線吻合較好。

圖3 單柱平均流速

圖4 單柱壓力系數(shù)

2 模擬結(jié)果

2.1 尾流平均風(fēng)速

在中心平面上(z=0)本文選取了3個剖面來監(jiān)測尾流平均風(fēng)速，分別是x/D=0，x/D=2.1，x/D=3。因工況較多，本文僅列出風(fēng)向角α=0°與α=24°情況下的模擬結(jié)果。

圖5給出了α=0°時平均流速沿y軸的分布規(guī)律。從圖中可以看出速度關(guān)于x軸對稱，這與結(jié)構(gòu)的對稱性相一致，隨著湍流的發(fā)展，速度分布曲線的形狀由“U”形逐漸變?yōu)椤癡”形，說明湍流尾流在周邊流體的影響下，速度逐漸趨近于周邊流體的速度。在y=0處，速度大小出現(xiàn)極小值，風(fēng)速大小從圓柱中心向兩側(cè)迅速增長，快速增長區(qū)域大約為2D，之后風(fēng)速趨于穩(wěn)定。圖6為α=0°平均流速云圖，可以直觀地看到下游圓柱完全浸沒在上游圓柱的尾流之中，氣流在圓柱之間形成了較大范圍的回流區(qū)，在下游圓柱尾部也形成了長度約為1D的回流區(qū)。

圖5 α=0°平均流速分布

圖6 α=0°平均流速云圖

圖7給出了α=24°時平均流速沿y軸的分布規(guī)律。由于圓柱對氣流的阻擋作用，在y方向上與圓柱位置對應(yīng)處的流速會出現(xiàn)兩個極小值，下游圓柱尾流流速的極值出現(xiàn)在圓柱的正后方，而上游圓柱尾流因為間隙流的擠壓效應(yīng)，速度極值的出現(xiàn)位置會向y軸負(fù)向偏移。圖8為α=24°的平均流速云圖，下游圓柱不再浸沒在上游圓柱的尾流之中，圓柱之間出現(xiàn)間隙流并分別在尾部出現(xiàn)回流區(qū)，下游圓柱回流區(qū)長度約為1.5D，上游圓柱回流區(qū)略小于上游圓柱，長度約為1.3D。

圖7 α=24°平均流速云圖

圖8 α=24°平均流速云圖

2.2 平均風(fēng)壓系數(shù)

沿圓柱周向布置36個監(jiān)測點，監(jiān)測點間隔為10 °。圖9為α=0°時的圓柱時均壓力系數(shù)分布。從圖中可以發(fā)現(xiàn)上游圓柱表面正對來流處(θ1=0°)的壓力系數(shù)最大，大約為1.0。隨著氣流向圓柱兩側(cè)擴(kuò)展，Cp值迅速減小，并逐漸轉(zhuǎn)變?yōu)樨?fù)壓，在θ1=70°處達(dá)到最大負(fù)壓，大約為-1.0，說明氣流邊界層在該處與上游圓柱分離，分離角θsep=70°，邊界層分離后壓力系數(shù)Cp逐漸回升并增大到一個較為穩(wěn)定的值，在圓柱背后形成穩(wěn)定的負(fù)壓分布；下游圓柱受到柱間回流區(qū)域的影響，表面壓力系數(shù)均為負(fù)值，θ2=30°與θ2=330°處的負(fù)壓值最大，最大負(fù)壓值約為-0.77，在θ2=70°與θ2=290°處負(fù)壓值最小，最小值為-0.14，說明氣流從上游圓柱分離之后在該處進(jìn)行了再附，隨著氣流向圓柱后方延伸，負(fù)壓值略微增大之后趨于穩(wěn)定，并在下游圓柱后方形成穩(wěn)定的壓力分布區(qū)域。

圖9 α=0°平均風(fēng)壓系數(shù)

圖10為α=24°時的圓柱時均壓力分布。在正對來流處(θ=0°)兩圓柱的壓力系數(shù)最大，最大值均為0.25。隨著氣流在圓柱壁面上的流動，壓力系數(shù)迅速減小，并逐漸轉(zhuǎn)變?yōu)樨?fù)壓，上游圓柱的氣流分離角為70 °與290 °，在該角度處的壓力系數(shù)均為-0.21，壓力系數(shù)在圓柱截面上下對稱分布；下游圓柱的氣流分離角為70 °與300 °，對應(yīng)的壓力系數(shù)分別為-0.20、-0.24，由于兩圓柱之間間隙流的影響，壓力系數(shù)在下游圓柱截面的分布是非對稱的。

圖10 α=24°平均風(fēng)壓系數(shù)

2.3 升、阻力系數(shù)與Strouhal數(shù)

圖11到圖14分別給出了α=0°與α=24°的升、阻力系數(shù)時程曲線。α=0°時力系數(shù)隨時間的變化不平穩(wěn)，但仍存在一定周期性特征。上游圓柱升力系數(shù)約為0，阻力系數(shù)約為0.93；下游圓柱升力系數(shù)約為0，阻力系數(shù)約為-0.24。由于上游圓柱尾流的撞擊作用，下游圓柱力系數(shù)的波動幅值要大于上游圓柱。

圖11 α=0°上游圓柱力系數(shù)

圖12 α=0°下游圓柱力系數(shù)

圖13 α=24°上游圓柱力系數(shù)

圖14 α=24°下游圓柱力系數(shù)

α=24°時力系數(shù)隨時間的變化具有較強的周期性。上游圓柱升力系數(shù)約為-0.02，阻力系數(shù)約為1.05；下游圓柱升力系數(shù)約為-0.38，阻力系數(shù)約為0.98。

圖15展示了時均升力系數(shù)Cl與風(fēng)向角α之間的關(guān)系，可以發(fā)現(xiàn)風(fēng)向角的變化不會影響上游圓柱升力系數(shù)，其值穩(wěn)定在0附近。下游圓柱的升力系數(shù)均小于零，說明下游圓柱有向上游圓柱靠近的趨勢，升力系數(shù)曲線形狀近似于“V”字形，風(fēng)向角在8 °以下時，升力系數(shù)絕對值隨著風(fēng)向角的增加而增大，即下游圓柱向上游圓柱靠近的趨勢愈加明顯，當(dāng)風(fēng)向角大于8 °時，升力系數(shù)絕對值又隨著風(fēng)向角的增加而減小。

圖15 升力系數(shù)與風(fēng)向角的關(guān)系

圖16展示了時均阻力系數(shù)Cd與風(fēng)向角α之間的關(guān)系，上游圓柱阻力系數(shù)隨著風(fēng)向角的增大而緩慢增加。下游圓柱的阻力系數(shù)隨著風(fēng)向角的增加而明顯增大，阻力系數(shù)由負(fù)值轉(zhuǎn)變?yōu)檎?，說明下游圓柱的運動趨勢由向上游運動逐漸轉(zhuǎn)變?yōu)橄蛳掠芜\動。

圖16 阻力系數(shù)與風(fēng)向角的關(guān)系

將升力系數(shù)時程進(jìn)行傅里葉變換并進(jìn)行無量綱化處理可以得到圓柱的Strouhal數(shù)。圖17給出了Strouhal數(shù)與角度α的關(guān)系，可以發(fā)現(xiàn)上下游圓柱的Strouhal數(shù)始終保持一致，α=0°時Strouhal數(shù)最小，最小值為0.169，對應(yīng)的渦脫頻率0.963 Hz；α=6°時Strouhal數(shù)最大，最大值為0.348，對應(yīng)的渦脫頻率1.985 Hz。

圖17 St數(shù)與風(fēng)向角的關(guān)系

2.4 渦脫形式

Q判據(jù)[19-20]為識別旋渦的一種方法，形式簡單，并且運用廣泛，所以本文以Q判據(jù)作為漩渦出現(xiàn)的判斷準(zhǔn)則，并給出了α=0°與α=24°時圓柱三維渦結(jié)構(gòu)。當(dāng)α=0°時，下游圓柱會阻礙上游尾渦的發(fā)展，導(dǎo)致兩圓柱之間形成許多小型漩渦，并且多個方向的漩渦交錯在一起使得流態(tài)十分復(fù)雜見圖18；當(dāng)α=24°時，上游圓柱的尾渦會繞過下游圓柱，并與下游圓柱尾渦相結(jié)合，結(jié)合后的尾流渦結(jié)構(gòu)比α=0°的更寬且更加有規(guī)律性見圖19，這也解釋了24 °風(fēng)向角比0 °風(fēng)向角的的升阻力系數(shù)更具有規(guī)律性的原因。

圖18 α=0°三維渦結(jié)構(gòu)

圖19 α=24°三維渦結(jié)構(gòu)

2.5 大升力原因分析

觀察圖15可以發(fā)現(xiàn)，下游圓柱在α=8°時受到了很大的升力，為了進(jìn)一步探討下游圓柱受到很大升力的機(jī)理，選取α=8°的某個時刻來分析局部流場特性，見圖20中右側(cè)圓柱。圖21為α=8°時的流線圖/風(fēng)速云圖，圖22為α=8°時的風(fēng)壓云圖。

圖20 α=8°渦量云圖

圖21 α=8°流線圖/風(fēng)速云圖

圖22 α=8°風(fēng)壓云圖

由圖21可看出，靠近下游圓柱下表面的地方，出現(xiàn)了高速間隙流，最大速度系數(shù)達(dá)到了1.75，對應(yīng)圖22中相同的地方出現(xiàn)了最大負(fù)壓區(qū)域，該負(fù)壓區(qū)域使圓柱受到較大的升力作用，這是出現(xiàn)較大升力的第一個原因；其次，觀察圖21中圓柱下表面的流線，形成了一個分離泡，說明在間隙流的影響下，圓柱下側(cè)邊界層在此分離，分離后會形成多個小型漩渦，小型漩渦與分離泡會在圓柱表面造成局部強負(fù)壓，這就是導(dǎo)致圓柱出現(xiàn)較大升力的第二個原因；此外，觀察圖21圓柱上表面的流線，也存在邊界層分離現(xiàn)象，在圖22中對應(yīng)位置出現(xiàn)了負(fù)壓極值，分離角為79 °，根據(jù)之前的分析得知α=0 °時分離角為70 °，比較可知α=8°時負(fù)壓極值點上移，導(dǎo)致下游圓柱會受到更大的升力作用。

以上為下游圓柱在α=8°時受到很大的升力的三個原因，并總結(jié)如下：

(1)高速間隙流引起的強負(fù)壓。

(2)小型漩渦與分離泡產(chǎn)生的局部強負(fù)壓。

(3)邊界層分離點上移。

3 結(jié)論

在Re=3900的情況下采用大渦模擬的方法對中心距為3D、風(fēng)向角0～24 °的錯列雙圓柱進(jìn)行了數(shù)值模擬，研究了空氣動力特性隨風(fēng)向角的變化，分析了大升力出現(xiàn)的原因。主要結(jié)論如下：

(1)隨著風(fēng)向角的改變，錯列圓柱會出現(xiàn)兩種繞流形態(tài)：0～8 °時，圓柱之間出現(xiàn)大范圍的回流區(qū)域，8～24 °時，柱間回流區(qū)域轉(zhuǎn)變?yōu)楦咚匍g隙流，這種轉(zhuǎn)變對圓柱氣動力有顯著影響。

(2)隨著風(fēng)向角的增大，升、阻力曲線會變得更有規(guī)律性，且下游圓柱升、阻力曲線波動幅值始終大于上游圓柱，這是因為在大風(fēng)向角時，上游圓柱尾渦會繞過下游圓柱并與下游圓柱尾渦相結(jié)合，形成更寬更規(guī)則的渦結(jié)構(gòu)。

(3)在8 °風(fēng)向時，下游圓柱升力系數(shù)會出現(xiàn)一個峰值，分析其原因主要有三點：一是高速間隙流引起的強負(fù)壓；二是小型漩渦與分離泡產(chǎn)生局部強負(fù)壓；三是邊界層分離點由70 °上升到79 °。