肖 姚， 秦 浩,2， 宋文武,2， 胡 帥

(1．西華大學(xué) 能源與動(dòng)力工程學(xué)院，成都 610039；2．西華大學(xué) 流體及動(dòng)力機(jī)械教育部重點(diǎn)實(shí)驗(yàn)室，成都 610039)

較大雷諾數(shù)下方柱繞流的數(shù)值模擬

肖 姚1， 秦 浩1,2， 宋文武1,2， 胡 帥1

(1．西華大學(xué) 能源與動(dòng)力工程學(xué)院，成都 610039；2．西華大學(xué) 流體及動(dòng)力機(jī)械教育部重點(diǎn)實(shí)驗(yàn)室，成都 610039)

基于k-ωSST模型和γ-Reθ轉(zhuǎn)捩模型，分析了雷諾數(shù)為1.4×104和2.2×104時(shí)二維方柱的流場(chǎng)特征，計(jì)算方柱的阻力系數(shù)和斯特勞哈爾數(shù)Sr，并對(duì)方柱平均流場(chǎng)的速度剖面進(jìn)行分析，綜合考慮不同湍流模型的計(jì)算成本及計(jì)算精度.結(jié)果表明：γ-Reθ轉(zhuǎn)捩模型比k-ωSST模型計(jì)算結(jié)果更加精確且所需計(jì)算成本較?。沪?Reθ轉(zhuǎn)捩模型能較好地預(yù)測(cè)方柱尾緣近壁回流區(qū)的位置和大小，且能較好地觀測(cè)流場(chǎng)在方柱尾緣的再附著現(xiàn)象;在求解實(shí)際大雷諾數(shù)工程問題中γ-Reθ轉(zhuǎn)捩模型為優(yōu)秀的湍流模型.

方柱繞流； 數(shù)值模擬； 阻力系數(shù)； 湍流模型； 時(shí)間步長(zhǎng)

鈍體繞流廣泛存在于流體機(jī)械、動(dòng)力工程和石油工程等領(lǐng)域中，如汽輪機(jī)葉片、換熱器管群及海洋鉆井平臺(tái)等均涉及鈍體繞流問題.鈍體繞流所引起的振動(dòng)會(huì)造成設(shè)備工作不穩(wěn)定，當(dāng)振動(dòng)頻率與結(jié)構(gòu)的固有頻率相接近時(shí)會(huì)產(chǎn)生共振，影響結(jié)構(gòu)的使用性能，振動(dòng)嚴(yán)重時(shí)甚至?xí)?duì)結(jié)構(gòu)造成破壞性損壞.方柱繞流存在復(fù)雜的渦結(jié)構(gòu)，而渦的脫落對(duì)于方柱所受的升力和阻力均有較大影響.因此，深入研究方柱繞流對(duì)解決工程實(shí)際中相關(guān)問題有著重要意義.

為了準(zhǔn)確高效地計(jì)算出方柱流場(chǎng)的分布以及了解渦的形成與脫落規(guī)律，很多學(xué)者對(duì)方柱繞流開展了研究工作[1-10].如Durao等[1]在水洞中利用多普勒激光測(cè)速儀(LDV)研究了方柱在雷諾數(shù)為1.4×104時(shí)的繞流問題，結(jié)果表明隨著湍流的波動(dòng)，在顫振速度最大區(qū)域的能量可達(dá)總能量的40%.Lyn等[2]基于LDV實(shí)驗(yàn)測(cè)試了雷諾數(shù)為2.2×104時(shí)方柱的流場(chǎng)，得出方柱尾緣的速度分布規(guī)律以及渦的分布與大小.張偉等[3]采用粒子圖像測(cè)速儀(PIV)和k-ωSST模型研究了雷諾數(shù)為2.2×104時(shí)的方柱繞流流場(chǎng)，得出k-ωSST模型計(jì)算結(jié)果與實(shí)驗(yàn)數(shù)據(jù)吻合較好.Bouris等[4]利用大渦模型數(shù)值模擬了小雷諾數(shù)(Re=250)時(shí)的二維方柱繞流，通過優(yōu)化網(wǎng)格質(zhì)量得出二維方柱的大渦模擬結(jié)果與實(shí)驗(yàn)數(shù)據(jù)吻合較好.梁思超等[6]基于分離渦(DES)模型研究了雷諾數(shù)為4×104時(shí)方柱繞流中馬蹄渦系統(tǒng)的發(fā)展過程.Franke等[7]采用k-ε渦黏模型和雷諾應(yīng)力模型研究了雷諾數(shù)為2.2×104時(shí)的二維方柱流場(chǎng),結(jié)果表明2種模型均不能較好的與實(shí)驗(yàn)結(jié)果吻合.謝志剛等[8]基于大渦模型,利用精細(xì)的邊界研究雷諾數(shù)為2.2×104時(shí)的方柱繞流，其結(jié)果可以媲美密網(wǎng)格的結(jié)果.

研究者們從實(shí)驗(yàn)、數(shù)值模擬和理論分析3個(gè)方面同時(shí)研究了方柱的流場(chǎng)，然而實(shí)驗(yàn)結(jié)果與數(shù)值模擬還存在一定的差異.文獻(xiàn)[7]的研究結(jié)果與實(shí)驗(yàn)結(jié)果相差較大，k-ε渦黏模型計(jì)算出的Sr和阻力因數(shù)CD均偏??；文獻(xiàn)[3]的結(jié)果雖然與文獻(xiàn)[2]的實(shí)驗(yàn)結(jié)果吻合較好，但其計(jì)算結(jié)果與實(shí)驗(yàn)結(jié)果在近方柱的尾流區(qū)存在一個(gè)明顯的跳躍，其對(duì)尾流區(qū)未能準(zhǔn)確預(yù)測(cè)；文獻(xiàn)[4]的計(jì)算結(jié)果較好，但計(jì)算的雷諾數(shù)較小且由于該湍流模型計(jì)算成本較高，將其應(yīng)用在實(shí)際工程中不現(xiàn)實(shí).實(shí)際工程多涉及大雷諾數(shù)下渦街的形成與脫落，為此需在考慮成本的情況下得出適合于工程大雷諾數(shù)的湍流模型.γ-Reθ轉(zhuǎn)捩模型[11]基于k-ωSST輸運(yùn)方程、間歇系數(shù)γ方程和動(dòng)量厚度雷諾數(shù)Reθt方程，相比k-ωSST模型，該模型可以較好地捕捉流場(chǎng)分離過程和湍流強(qiáng)度的變化.筆者基于Fluent軟件分析方柱在不同Re下的流場(chǎng)，通過對(duì)比分析方柱的Sr、阻力系數(shù)CD、升力系數(shù)CL和流場(chǎng)的速度剖面等參數(shù)來(lái)對(duì)比不同湍流模型的評(píng)估結(jié)果.

1 數(shù)值方法

1.1 控制方程

二維黏性不可壓牛頓流體運(yùn)動(dòng)的基本控制方程為連續(xù)性方程和Navier-Stokes方程：

(1)

(2)

1.2k-ωSST模型

k-ωSST模型基于Menter[12]研究成果,該湍流模型在多種流態(tài)下較standardk-ω和BSLk-ω湍流模型更加精確和可靠，其表達(dá)式如下：

Gk-Yk+Sk

(3)

Gw-Yw+Dw+Sw

(4)

式中：Gk、Gw分別代表湍動(dòng)能和w的產(chǎn)生；Γk、Γw分別為k和w的有效擴(kuò)散率；Yk、Yw分別為由湍動(dòng)產(chǎn)生的耗散率；Sk、Sw為源項(xiàng).

1.3γ-Reθ轉(zhuǎn)捩模型

Menter等[11]提出了γ-Reθ轉(zhuǎn)捩模型，該湍流模型由k-ωSST輸運(yùn)方程、間歇系數(shù)γ方程和動(dòng)量厚度雷諾數(shù)Reθt方程組成，其中間歇系數(shù)γ用來(lái)描述邊界層的狀態(tài)，并控制過渡點(diǎn)下游湍動(dòng)能的生成項(xiàng)，式(5)為間歇系數(shù)無(wú)量綱輸運(yùn)方程的守恒形式.動(dòng)量厚度雷諾數(shù)Reθt通過捕獲非局部影響的湍流強(qiáng)度來(lái)控制由邊界層內(nèi)間歇系數(shù)的增長(zhǎng)和當(dāng)?shù)販u量雷諾數(shù)構(gòu)成的轉(zhuǎn)捩觸發(fā)準(zhǔn)則，其表達(dá)式如下：

(5)

(6)

式中：μt、σθ、σγ為模型系數(shù);Pγ1、Eγ1代表轉(zhuǎn)捩源項(xiàng);Pγ2、Eγ2分別為破裂系數(shù)和再附著系數(shù);Pθ為源項(xiàng).

1.4 方柱幾何模型與網(wǎng)格劃分

圖1 方柱流場(chǎng)計(jì)算域

圖2 網(wǎng)格劃分

采用Fluent軟件進(jìn)行數(shù)值模擬，設(shè)定初始條件:u1=U0，u2=0；方柱入口定義為速度邊界條件，來(lái)流為均勻速度U0；出口定義為自由流邊界條件，方柱壁面與計(jì)算域邊界定義為無(wú)滑移壁面.

2 結(jié)果與討論

2.1 湍流模型對(duì)比

首先以文獻(xiàn)[1]和文獻(xiàn)[2]的實(shí)驗(yàn)結(jié)果作為數(shù)值模擬的參照對(duì)象，計(jì)算不同Re下各湍流模型最佳的模擬時(shí)間步長(zhǎng)，得出方柱阻力系數(shù)和升力系數(shù)的時(shí)程曲線，將升力系數(shù)時(shí)程曲線進(jìn)行快速傅氏變換(FFT)得到渦街脫落頻率f.其中阻力系數(shù)CD、升力系數(shù)CL及Sr的定義分別為：

(a) 方柱上壁面Y+

(b) 方柱尾緣壁面Y+

(7)

式中：FD和FL分別為方柱所受的阻力和升力.

表1 不同Re下方柱繞流的Sr和

(a) k-ω SST模型

(b) γ-Reθ轉(zhuǎn)捩模型

2.2 流場(chǎng)分析

由第2.1節(jié)可知γ-Reθ轉(zhuǎn)捩模型可以較好地評(píng)估實(shí)際工程大Re下的流場(chǎng)變化，為了與實(shí)驗(yàn)進(jìn)行對(duì)比，下文只分析Re為2.2×104時(shí)方柱尾緣渦的脫落過程.圖5為方柱平均流場(chǎng)速度分布曲線與水平湍動(dòng)能變化曲線.從圖5(a)可以看出，2種湍流模型計(jì)算結(jié)果在方柱前緣的速度場(chǎng)與文獻(xiàn)[1]的實(shí)驗(yàn)結(jié)果吻合.在1.8

(a) 水平速度剖面

(b) 豎向速度剖面

(c) 水平剖面湍動(dòng)能

綜上，γ-Reθ轉(zhuǎn)捩模型與三維大渦模擬結(jié)果以及實(shí)驗(yàn)結(jié)果較為接近，但由于尾流近壁回流區(qū)是三維流動(dòng)，本文結(jié)果還存在一定的偏差.

2.3 尾渦脫落

由表1可知，雷諾數(shù)為2.2×104時(shí)采用γ-Reθ轉(zhuǎn)捩模型計(jì)算得出的方柱渦脫落頻率為17 Hz，周期為0.059 s.圖6給出了Re為2.2×104時(shí)方柱尾緣渦脫落的半個(gè)周期過程，其中τ為時(shí)間.從圖6可以看出，在方柱前緣-0.05 m左右流場(chǎng)開始分離，并在方柱兩側(cè)產(chǎn)生回流區(qū)，在方柱尾緣也有短暫再附著現(xiàn)象.在這半個(gè)周期內(nèi)尾渦從方柱尾緣上方逐漸形成，并在尾緣擴(kuò)大，在一定程度時(shí)渦開始脫落，其形成尾渦的大小與方柱的特征長(zhǎng)度相近.

(a) τ=0 s

(b) τ=0.000 2 s

(c) τ=0.000 4 s

(d) τ=0.000 6 s

(e) τ=0.000 8 s

(f) τ=0.001 s

(g) τ=0.001 2 s

(h) τ=0.001 4 s

(i) τ=0.001 6 s

(j) τ=0.001 8 s

圖6 方柱尾渦脫落過程

Fig.6 Shedding process of vortex from the square cylinder

3 結(jié) 論

(1) 與其他湍流模型相比，γ-Reθ轉(zhuǎn)捩模型可以精確地評(píng)估大Re下方柱的Sr、CD和CL.

(2) 與相關(guān)實(shí)驗(yàn)結(jié)果對(duì)比，γ-Reθ轉(zhuǎn)捩模型與實(shí)驗(yàn)結(jié)果的速度剖面吻合較好，可以較好地預(yù)測(cè)方柱尾緣回流區(qū)的大小和位置，計(jì)算結(jié)果與三維大渦模擬結(jié)果較為接近.

(3)γ-Reθ轉(zhuǎn)捩模型可以準(zhǔn)確地分析方柱壁面渦的形成與脫落過程以及流場(chǎng)再附著現(xiàn)象，很好地揭示了繞流的流動(dòng)機(jī)理.

[1] DURAO D, HEITOR M V, PEREIRA J C F. Measurements of turbulent and periodic flows around a square cross-section cylinder[J]. Experiments in Fluids, 1988, 6(5): 298-304.

[2] LYN D A, EINAV S, RODI W, et al. A laser-Doppler velocimetry study of ensemble-averaged characteristics of the turbulent near wake of a square cylinder[J]. Journal of Fluid Mechanics, 1995, 304: 285-319.

[3] 張偉, 葛耀君. 方柱繞流粒子圖像測(cè)速試驗(yàn)與數(shù)值模擬[J]. 同濟(jì)大學(xué)學(xué)報(bào)(自然科學(xué)版), 2009, 37(7): 857-861, 892.

ZHANG Wei, GE Yaojun. Particle image velocimetry study and numerical simulation of turbulent near wake of square cylinder[J]. Journal of Tongji University(Natural Science), 2009, 37(7): 857-861, 892.

[4] BOURIS D, BERGELES G. 2D LES of vortex shedding from a square cylinder[J]. Journal of Wind Engineering and Industrial Aerodynamics, 1999, 80(1/2): 31-46.

[5] 李慧君, 王炯. 基于水平圓管外加裝排液板降膜流動(dòng)的數(shù)值模擬[J]. 動(dòng)力工程學(xué)報(bào), 2016, 36(4): 265-270.

LI Huijun, WANG Jiong. Numerical simulation of falling film flow on horizontal circular tubes with drainage strips[J]. Journal of Chinese Society of Power Engineering, 2016, 36(4): 265-270.

[6] 梁思超, 張曉東, 康順, 等. 方體繞流的DES模擬及其旋渦結(jié)構(gòu)分析[J]. 工程熱物理學(xué)報(bào), 2015, 36(3): 535-540.

LIANG Sichao, ZHANG Xiaodong, KANG Shun, et al. DES simulation of the flow around a cube and analysis of the vortex structure[J]. Journal of Engineering Thermophysics, 2015, 36(3): 535-540.

[7] FRANKE R, RODI W. Calculation of vortex shedding past a square cylinder with various turbulence models[M]//DURST F, FRIEDRICH R, LAUNDER BE, et al. Turbulent Shear Flows 8.Berlin Heidelberg:Springer, 1993.

[8] 謝志剛, 許春曉, 崔桂香, 等. 方柱繞流大渦模擬[J]. 計(jì)算物理, 2007, 24(2): 171-180.

XIE Zhigang, XU Chunxiao, CUI Guixiang, et al. Large eddy simulation of flows around a square cylinder[J]. Chinese Journal of Computational Physics, 2007, 24(2): 171-180.

[9] OKAJIMA A. Strouhal numbers of rectangular cylinders[J]. Journal of Fluid Mechanics, 1982, 123: 379-398.

[10] SAMANI M, BERGSTROM D J. Effect of a wall on the wake dynamics of an infinite square cylinder[J]. International Journal of Heat and Fluid Flow, 2015, 55: 158-166.

[11] MENTER F R, LANGTRY R B, LIKKI S R, et al. A correlation-based transition model using local variables—Part I: model formulation[J]. Journal of Turbomachinery, 2006, 128(3): 413-422.

[12] MENTER F R. Two-equation eddy-viscosity turbulence models for engineering applications[J]. AIAA Journal, 1994, 32(8): 1598-1605.

[13] BREUER M, POURQUIé M. First experiences with LES of flows past bluff bodies[C]//Proceedings of the 3rd International Symposium of Engineering Turbulence Modelling and Measurements. Amsterdam,Holland:Elsevier,1996: 177-186.

[14] SAHA A K, MURALIDHAR K, BISWAS G. Experimental study of flow past a square cylinder at high Reynolds numbers[J]. Experiments in Fluids, 2000, 29(6): 553-563.

Numerical Simulation of the Flow Around a Square Cylinder at Higher Reynolds Numbers

XIAOYao1,QINHao1,2,SONGWenwu1,2,HUShuai1

(1.School of Energy and Power Engineering, Xihua University, Chengdu 610039, China; 2. Key Laboratory of Fluid and Power Machinery, Ministry of Education, Xihua University, Chengdu 610039, China)

Based on thek-ωSST model andγ-Reθtransition model, a study was conducted on the flow characteristics around a square cylinder respectively at the Reynolds number of 1.4×104and 2.2×104, so as to calculate the drag coefficient and Strouhal number of the square cylinder, analyze its average velocity profiles, and finally to find out the calculation cost and precision of different turbulent models. Results show that compared withk-ωSST model, theγ-Reθtransition model is more superior in calculation precision and cost, which is more effective in predicting the position and magnitude of the recirculation region behind the square cylinder, and in analyzing the reattachment of separated flows on side walls near the rear corner, and therefore is a practical model for engineering applications at high Reynolds numbers.

flow around square cylinder; numerical simulation; drag coefficient; turbulent model; time step

1674-7607(2017)08-0679-06

O357

A

470.10

2016-08-04

2016-10-24

流體及動(dòng)力機(jī)械教育部重點(diǎn)實(shí)驗(yàn)室開放基金資助項(xiàng)目(szjj2016-003);西華大學(xué)自然科學(xué)重點(diǎn)基金資助項(xiàng)目(Z1510418);國(guó)家自然科學(xué)基金資助項(xiàng)目(51408505)

肖 姚(1990-)，男，四川資陽(yáng)人，碩士研究生，研究方向?yàn)榱黧w機(jī)械. 秦 浩(通信作者)，男，講師，博士，電話(Tel.)：13547904041;E-mail:qh.email.@foxmail.com.