呂啟兵，楊 斌，楊忠超，李鵬浩

(1.重慶交通大學河海學院，重慶400074;2.溫州港集團有限公司，浙江溫州325100)

鈍體繞流問題廣泛存在于橋梁、海洋港口、航空航天、交通運輸?shù)裙こ炭茖W領域，如置于水中的橋墩、石油鉆井平臺的墩柱、港口碼頭的樁柱、飛機機翼繞流、汽車的風阻等。而方柱作為一種典型的鈍體，一直是流體力學領域繞流問題的熱點之一。關于方柱繞流問題的實驗研究也已經(jīng)進行得相當廣泛，許多研究者通過各種方法對方柱的繞流進行了實驗模擬研究，D.A.Lyn，等［1-2］提供了大量的實驗模擬資料，為實際的工程問題提供了一定的參考。

隨著計算機計算能力的提升、數(shù)值算法的改進與計算流體力學(CFD)的突飛猛進，數(shù)值模擬技術已成為鈍體繞流問題研究的一種有效手段，樁柱繞流問題的數(shù)值模擬也得到了很大的發(fā)展。早期的Hiromasa Kawai，等［3］對雷諾數(shù)為 200 的串列方柱的繞流用差分方法做了數(shù)值模擬;王遠成，等［4］人采用湍流RNG k-ε模型對鈍體繞流流場的不穩(wěn)定、非定常流動進行了數(shù)值模擬;S.Muralami，等［5］使用標準k-ε模型對鈍體繞流進行了研究，并且利用了大渦模型(LES)進行了數(shù)值計算;陳素琴，等［6］數(shù)值模擬了雷諾數(shù)為10 000情況下兩串列方柱繞流的干擾問題，并分析了不同間距下各方柱的升阻力特征及流場特性。

由于內河特殊的水文環(huán)境(大水位差變化與較大的流速)，近年來內河碼頭運用全直柱框架結構作為一種新的碼頭結構型式，在考慮施工過程、工程條件與經(jīng)濟因素等方面時，全直樁方柱框架碼頭已十分普遍。筆者利用流體軟件、采用RNG k-ε模型對串聯(lián)雙方柱繞流進行了數(shù)值模擬，分析了每根柱升、阻力特性與不同樁距下的流場演變情況，為實際工程問題提供一定參考。

1 數(shù)值模型與模型驗證

1.1 控制方程

對于不可壓縮黏性流體，在直角坐標系下，流體的運動可用Navier-Stokes方程［7］與動量方程描述：

模擬計算采用RNG k-ε模型，該模型由Yakhot和Orzag提出，模型是從瞬態(tài)N-S方程中推出的，通過大尺度運動和修正后的黏度項體現(xiàn)小尺度的影響，而使這些小尺度運動有系統(tǒng)地從控制方程中去除，所得到的k方程和ε方程與標準k-ε模型非常相似：

式中：Gk是由層流速度梯度產(chǎn)生的湍流動能k的產(chǎn)生項C1ε=1.44，η。=4.377，β =0.012，C2ε=1.68;ε 為湍動耗散率。

RNG k-ε模型來源于嚴格的統(tǒng)計技術，是對標準k-ε模型的重要改進。RNG模型在ε方程中加入了一個Rε條件，從而反映了主流的時均應變率Eij，這樣該模型中產(chǎn)生項不僅與流動情況有關，而且在同一問題中也還是空間坐標的函數(shù)，從而可以更好的處理高應變率及流線彎曲程度較大的流動。

1.2 數(shù)值算法

計算區(qū)域網(wǎng)格采用分塊結構化網(wǎng)格(圖1)，網(wǎng)格總量在4.9萬～18.9萬之間，利用有限體積法(FVM)對方程進行離散，壓力基求解器與壓力2階離散格式，動量、湍動能和湍動耗散率選項全部為2階迎風格式，計算模型選用RNG k-ε模型，PISO算法。

1.3 邊界條件

進口邊界條件：u=1，v=0。

出口邊界條件：自由出流。

上下邊界及方柱表面為壁面邊界條件，即u=0，v=0的無滑移邊界條件。

雷諾數(shù)：Re=UD/υ=22 000。式中：U為來流流速;D為方柱邊長;υ為流體的動力黏性系數(shù)。

1.4 計算區(qū)域與網(wǎng)格劃分

計算區(qū)域如圖1，單柱計算區(qū)域為22D×13D(D為方柱邊長)，距入口邊界距離6D，考慮上下邊界對方柱周圍流場的影響，上下邊界距方柱都為6D，與出口距離為15D;雙方柱間距為L，其余設置與單方柱相同。

圖1 單方柱與雙方柱計算區(qū)域與網(wǎng)格Fig.1 Calculation area and grid of single square cylinder and double square cylinder

1.5 模型驗證

方柱繞流常用于模擬研究流體的荷載，即對建筑物和構筑物的阻力、升力和柱上渦街脫落頻率。為了驗證數(shù)值模擬的準確性，首先計算了Re=22 000情況下的單方柱繞流問題。方柱上阻力系數(shù)、升力系數(shù)與斯特羅哈數(shù)的定義為：

式中：Fd為阻力，F(xiàn)l為升力，u為來流流速，D為方柱邊長，f為渦街脫落頻率，即方柱上渦街的脫落頻率，、Cl分別代表阻力系數(shù)均值與升力系數(shù)。

從表1可以看出，模擬計算結果與 D.A.Lyn，等［1-2］試驗結果非常接近，模擬計算誤差在7%以內，模擬得到的單方柱情形下的升、阻力系數(shù)如圖2。

表1 單方柱繞流模擬計算結果比較Table 1 Comparison results of flow around square cylinder simulation

圖2 Re=22 000時單方柱升、阻力系數(shù)Fig.2 Lift and drag coefficients of single square cylinder when Re=22 000

2 模擬結果及分析

2.1 升、阻力分析

筆者主要關注不同間距情況下各方柱的阻力、升力系數(shù)與變化情況，簡略分析了流場演變情況。

表2給出了不同間距下上游方柱的阻力、升力與斯特勞哈數(shù)。在L/D=1時，兩方柱的升力都非常的小，表明兩方柱此時的脈動微弱;在L/D=2時，上游方柱阻力較其他情況時有所減小，說明間距較小時，兩柱存在干擾作用。隨著間距比的不斷增大，上游方柱阻力(圖3)與斯特勞哈數(shù)St趨于單柱繞流情況下的情形;在L/D=3時上下方柱的升力都變化較大，此種間距比下游方柱增加的升力有較大增加，主要是因為：上游方柱的尾流對下游方柱的來流有很大影響，此間距下上游方柱的尾渦已開始脫落，擴散出去，改變了下游方柱的來流(圖4)，導致上下游方柱的升力都有大幅度的增加。

表2 上下游方柱阻力計算結果Table 2 Calculation results of drags of square cylinders at upstream and downstream

圖3 上下游方柱阻力與升力系數(shù)走勢Fig.3 Trend of drag and lift coefficient of upstream and downstream square cylinder

圖4 流線圖局部放大Fig.4 Local amplification of flow line

從表2可以看出，在L/D=1與2間距下時下游方柱阻力為負值，此時方柱受力為水流的相反方向，其主要原因是：當兩方柱間距較小時(如L/D=1與2)，下游方柱受上游方柱尾流漩渦影響，柱前形成為負壓，從而導致這樣的結果，此結果與一般的實驗情況基本吻合。此時兩方柱的升力系數(shù)非常小，斯特勞哈數(shù)St較單柱情況時變小，說明方柱的脈動較小，而渦的形成頻率較單柱時加快。隨后方柱脈動開始變大，上游方柱阻力有所減小，但減小得并不十分明顯，而后兩方柱的脈動放緩。兩方柱的間距比L/D值的變化對上游方柱所受阻力影響不大，主要為上游方柱對下游方柱的影響;當L/D≥5時，上游方柱的所受阻力、升力、斯特勞哈數(shù)均與單柱情況幾乎一致。隨著方柱間距的增大，下游方柱所受阻力呈變大的趨勢，升力較單柱時大了近一倍左右。其原因為下游方柱位于上游方柱分離的脈動高速區(qū)附近，上游方柱脫落的渦直接作用于下游方柱，產(chǎn)生更小尺度的渦，在其前角產(chǎn)生額外的角動量，從而產(chǎn)生這樣的結果。

當兩方柱靠得較近時，阻力與升力變化都比較劇烈，其變化的趨勢各不相同，當距離較遠時，兩方柱的阻力趨于一致，升力也趨于穩(wěn)定，這在工程設計建設中對水流荷載的考慮提供一定參考價值。

2.2 流場演變

圖5給出5種情形下的瞬時流線圖。間距較小時，只在方柱背后產(chǎn)生漩渦，兩方柱附近的流場改變不大，尾流區(qū)較遠處沒有渦的擺動;隨著L的增大，下游方柱前后流態(tài)變化明顯，尾流區(qū)形成周期性的渦街，兩柱都有規(guī)律性渦脫落，并且下游方柱渦街擺動越來越大，說明L的增大增加了兩方柱尾流的不穩(wěn)定性，而上游方柱周圍流場幾乎不受下游方柱影響。從圖中對比可以看出，串列雙方柱情形下流場只有在L/D=1間距比下才出現(xiàn)了對稱情況，并伴有對稱尾流，其余情況主要為偏轉尾流。在L/D=2時兩柱之間有較大回流區(qū)，隨著L的增大，回流區(qū)有一定的減弱。

圖5 流線圖Fig.5 Flow line

圖6給出了 L/D 為1，2，3，6，10情況下的渦量圖。可以發(fā)現(xiàn)，隨著間距的增大，渦量較大的地方由方柱兩側移向方柱的后方，而后方柱卻剛好相反。

圖6 渦量圖Fig.6 The vorticity

3 結論

筆者主要研究了兩串列方柱擾流的干擾情況。模擬計算過程中發(fā)現(xiàn)，兩方柱間網(wǎng)格精度尤為重要，直接關系到下游方柱繞流特性模擬所得的數(shù)值。經(jīng)研究得到以下結論：

1)當兩方柱串列排列時，上游方柱所受阻力幾乎無變化，與單方柱情形時大致相同。

2)下游方柱所受阻力在隨著間距的增大而先減小后增大，且均小于上游方柱所受阻力，屏蔽效應明顯;由于上游方柱的存在，下游方柱升力明顯增大。

3)間距較大時，上下游方柱的脈動頻率幾乎相等。

4)在串聯(lián)時兩方柱之間間距L的變化對流場的改變作用明顯。

［1］Lyn D A，Einav S，Rodi W.Laser doppler velocimetry study of ensemble-averaged characteristic of the turbulent flow near of a square cylinder［J］.Fluid Mechanics，1995，304：285-319.

［2］Lyn D A，Rodi W.The flapping shear layer formed by flow separation from the forward corner of a square cylinder［J］.Fluid Mesh，1994，267：353-376.

［3］Hiromasa Kawai，Kiyoshi Fujinamic.Numerical Simulation of Flow around Square Prisms in Tandem Arrangements［C］//9thICME.Houston Texas：Gulf Publishing，1995：185-186.

［4］王遠成，吳文權.方柱繞流流場的RNG方法模擬研究［J］.水動力學研究與進展：A，2004，19(增刊1)：916-920.Wang Yuancheng，Wu Wenquan.Numerical simulation of flow around square cylinder using RNG k-ε turbulence model［J］.Journal Hydrodynamics：A，2004，19(S1)：916-920.

［5］Muralami S，Mochida A.3-D numerical simulation of airflow around a cubic model by means of the k-ε model［J］.Journal of Wind Engineering and Industrial Aerodynamics，1988，31(1)：327-338.

［6］陳素琴，黃自萍，沈劍華，等.兩串列繞流的干擾數(shù)值模擬研究［J］.同濟大學學報：自然科學版，2001，29(3)：320-325.Chen Suqin，Huang Ziping，Shen Jianhua，et al.Numerical computation of the flow around two square cylinders in tandem arrangement［J］.Journal of Tongji University：Natural Science，2001，29(3)：320-325.

［7］張遠君.流體力學大全［M］.北京：北京航空航天大學出版社，1991.Zhang Yuanjun.Fluid Mechanics［M］.Beijing：Beijing University of Aeronautics and Astronautics Press，1991.

［8］Rodi W.Simulation of flow past buildings with statistical turbulence models［J］.Wind Climate in Cities，1995，277：649-668.