等邊布置四方柱繞流渦脫落規(guī)律的數(shù)值研究

王圣潔，張 偉，劉 歡

（天津城建大學(xué) 能源與安全工程學(xué)院 天津 300384）

流體流過多個(gè)柱體的繞流現(xiàn)象，在諸多工程領(lǐng)域都有涉及，例如橋梁的支柱與吊索或者成群排列的換熱器等都涉及典型的多柱體繞流問題.其中四等邊方柱的繞流，作為較為常見的單元，有著重要的研究意義.

單方柱繞流問題因?yàn)榱鲃雍唵?、流動現(xiàn)象豐富，經(jīng)常被作為湍流模型的標(biāo)準(zhǔn)驗(yàn)證算例[1].因此，目前國內(nèi)外的相關(guān)研究大多集中于單方柱的繞流，對于多方柱的繞流研究相對較少涉及.就多柱體繞流而言，研究者們似乎對于多圓柱繞流問題給與了更多的關(guān)注，也取得了很多研究成果.張愛社等[2]對等直徑四圓柱在雷諾數(shù)Re=200時(shí)的繞流特性進(jìn)行了探究.張鵬飛等[3]針對三圓柱的繞流問題，在雷諾數(shù)Re=200時(shí)的各種流場特性進(jìn)行了實(shí)驗(yàn)和數(shù)值研究.Farrant等[4]對四圓柱和雙圓柱的繞流問題進(jìn)行研究，結(jié)果表明四圓柱與雙圓柱在繞流特性與流場阻力方面有明顯不同.李景銀等[5]針對雷諾數(shù)Re=100、間距比等于4條件下的正方形排列四圓柱繞流現(xiàn)象進(jìn)行了研究.Lam等[6]采用數(shù)值方法對二維和三維四圓柱流場流動及動力學(xué)特性進(jìn)行研究.近幾年，一些學(xué)者也陸續(xù)開展了對于多方柱繞流問題的研究，Shams等[7]采用LBM方法對低雷諾數(shù)下四等邊方柱的繞流進(jìn)行了數(shù)值模擬，但對于漩渦脫落的規(guī)律，未得出精確的結(jié)果.趙心廣[8]和張偉等[9]分別對較高雷諾數(shù)（Re=8×104）和低雷諾數(shù)（Re=100）的三方柱繞流問題進(jìn)行了數(shù)值研究.因?yàn)榉街兴膫€(gè)拐角而與圓柱的幾何特征不同，所以無論是繞流流場特性還是流體動力特性，兩者都有很大的差別.為了更加全面地了解方柱群的繞流特征，特別是對繞流柱體尾渦區(qū)漩渦生成、演化、脫落機(jī)制以及柱體間的干擾特性等尚需進(jìn)行細(xì)致的研究.

本文對雷諾數(shù)Re=100條件下呈正方形和菱形兩種幾何布置的四方柱繞流問題進(jìn)行數(shù)值研究.著重探究等邊布置四方柱在不同間距比條件下，不同的繞流流場特性、方柱間干擾效應(yīng)、方柱間流體動力特性，并嘗試準(zhǔn)確捕捉兩種幾何陣列情形下柱體繞流的臨界間距比.

1 數(shù)值方法

1.1 控制方程和參數(shù)設(shè)定

二維黏性不可壓非定常流動的控制方程為連續(xù)性方程和不可壓Navier-Stokes方程，即

式中：t為時(shí)間；vi為速度分量；p為壓力；ρ為流體密度；μ為流體動力黏度.

采用有限體積法求解二維黏性不可壓非定常Navier-Stokes方程，對流采用二階迎風(fēng)格式，壓力采用二階離散格式，壓力與速度耦合采用SIMPLEC算法求解.通過對單方柱繞流模型算例的計(jì)算，驗(yàn)證了本文算法的有效性.

升力系數(shù)、阻力系數(shù)、間距比、St數(shù)等無量綱數(shù)，由以下方程定義，即

式中：Cd為阻力系數(shù)，正方向與來流方向相反；CL為升力系數(shù)，正方向與來流方向垂直；ρ為流體密度；FD為阻力；U為來流速度；FL為升力；L為方柱中心之間的距離；D為方柱邊長；g為間距比；St為斯特勞哈爾數(shù)；f為渦脫落頻率.

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

計(jì)算域?yàn)殚L方形，長度D=1，方柱邊長為量綱1.四方柱呈正方形布置（布置1）和菱形布置（布置2）兩種情況.為了表達(dá)方便，對四個(gè)方柱分別進(jìn)行編號為A、B、C、D，具體布置形式及計(jì)算區(qū)域如圖1所示.

相關(guān)邊界條件，入口邊界為速度入口，無量綱速度U=1；出口邊界為自由出流；方柱表面為無滑移壁面；上下壁面為滑移壁面.

圖1 四方柱布置形式及計(jì)算域示意

1.3 網(wǎng)格劃分

采用GAMBIT軟件對計(jì)算域進(jìn)行非結(jié)構(gòu)三角形網(wǎng)格劃分，并對各個(gè)四方柱表面進(jìn)行網(wǎng)格加密，網(wǎng)格劃分如圖2所示.

圖2 柱體表面網(wǎng)格示意

2 結(jié)果與討論

2.1 算例網(wǎng)格無關(guān)性及結(jié)果驗(yàn)證

首先針對雷諾數(shù)為100條件下的單方柱繞流問題進(jìn)行網(wǎng)格無關(guān)性驗(yàn)證.采用網(wǎng)格步長為1進(jìn)行初步試算并評估試算結(jié)果，在流場趨勢基本正確的情況下逐步加密網(wǎng)格進(jìn)行計(jì)算，并將多次計(jì)算結(jié)果進(jìn)行對比.結(jié)果表明，采用網(wǎng)格步長為0.3時(shí)阻力系數(shù)Cd隨網(wǎng)格加密的變化差異在誤差范圍0.3%以內(nèi).

對計(jì)算結(jié)果進(jìn)行分析并與經(jīng)典算例對比，見表1，以驗(yàn)證數(shù)學(xué)模型與數(shù)值模擬方法的準(zhǔn)確性.

表1 單方柱繞流模擬結(jié)果比較

2.2 流動特性

本文對雷諾數(shù)為 100，間距比為1.5，2，2.5，3，3.5，4，4.5的兩種方柱排列的繞流流場進(jìn)行了模擬.圖3給出了布置1在不同間距比下的渦量.

由圖3可見，對于低雷諾數(shù)下正方形排列的四方形柱體，隨著間距比的逐漸增大，其流場結(jié)構(gòu)呈現(xiàn)出四種特征：①當(dāng)間距比較小時(shí)，四個(gè)柱體非常接近，此時(shí)流場形態(tài)類似于一個(gè)大的孤立方柱，在尾流區(qū)形成一組渦街，如圖3a所示；②當(dāng)間距比稍微增大時(shí)，上流柱體脫落的剪切層附著于下流方柱，在兩柱體間出現(xiàn)扭曲的渦場，如圖3c所示；③當(dāng)間距比再繼續(xù)增大，大于臨界間距比時(shí)，在上游柱體后發(fā)生漩渦的脫落，但漩渦在尚未完全形成時(shí)即撞擊到下游方柱，如圖3d所示；④上游柱體后的漩渦完全脫落，擊打在下游柱體上，如圖3e所示.

Farrant[10]在模擬四圓柱繞流時(shí)認(rèn)為，漩渦的脫落都是同步的，但從相位上來說，同相脫落和反相脫落都可能出現(xiàn)并會保持一段時(shí)間穩(wěn)定.在本文的四方柱繞流中，漩渦的脫落也符合這個(gè)規(guī)律.所有穩(wěn)定的漩渦脫落過程中，漩渦的脫落都是同步的，且同相位和反相位脫落都會發(fā)生.

圖4給出了布置1在不同間距比下的渦量.

由圖4可見，對于低雷諾數(shù)下菱形排列的等邊四方柱體，隨著間距比的逐漸增大，其流場結(jié)構(gòu)呈現(xiàn)出三種特征：①當(dāng)間距比較小時(shí)，四個(gè)柱體類似于一個(gè)大的孤立方柱，流場特性與單方柱相似，如圖4a所示；②當(dāng)間距比稍微增大時(shí)，流體介質(zhì)從四個(gè)柱體中間穿過，由于通道較小而造成狹窄射流，菱形排布的方柱之間產(chǎn)生強(qiáng)烈的相互作用，如圖4d所示；③當(dāng)間距比再繼續(xù)增大，大于臨界間距比時(shí)，四個(gè)方柱后都產(chǎn)生漩渦脫落，柱B下游形成了成熟的渦脫落連續(xù)撞擊到柱C上.此間距比下，流場內(nèi)不同方柱后的渦結(jié)構(gòu)相互干擾，流動特性以及擾動最復(fù)雜，如圖3e所示.

2.3 阻力和升力系數(shù)特性

布置1的各柱阻力系數(shù)和升力系數(shù)與間距比之間的關(guān)系折線如圖5所示.

由圖5可得以下結(jié)論.

（1）上游柱體的阻力系數(shù)變化不大，而下游柱體的阻力系數(shù)變化較大.當(dāng)間距比超過臨界間距比后，上游柱體后漩渦脫落，并擊打到下游柱體上.此時(shí)，下游柱體的阻力系數(shù)將增大.但隨著間距比的進(jìn)一步增加，阻力系數(shù)將會減小，并逐漸趨向于單個(gè)柱體繞流時(shí)的阻力系數(shù).

（2）方柱A與方柱C的平均升力系數(shù)折線關(guān)于x=0這條線呈對稱分布.方柱B與方柱D的平均升力系數(shù)折線關(guān)于x=0這條線呈對稱分布.

（3）間距比越大，平均升力系數(shù)越趨向于零，此時(shí)流動狀態(tài)越趨向于單個(gè)柱體的繞流.而間距比較小時(shí)，上下兩側(cè)柱體的平均升力系數(shù)之和接近于零，表明此時(shí)上方和下方的兩個(gè)柱體更像作為一個(gè)整體受到升力的作用.

（4）結(jié)合繞流流場和升阻力系數(shù)變化規(guī)律可知，臨界間距比在2.5到3之間.

布置2的各柱阻力系數(shù)和升力系數(shù)與間距比之間的關(guān)系折線如圖6所示.

由圖6可得以下結(jié)論.

（1）間距比小于2.5時(shí)，四方柱的阻力系數(shù)變化都不大，柱C的阻力系數(shù)比較小，柱A柱D的阻力系數(shù)比較大且升力系數(shù)對稱分布變化較大.這是因?yàn)橹鵄柱D對稱分布阻力特性相同，柱C幾乎完全位于其他三根方柱的尾渦之中，類似被屏蔽，阻力較小.

（2）隨著間距比的逐漸增大，柱A的阻力系數(shù)大幅度減小，柱C的阻力系數(shù)大幅度增加.這是因?yàn)橹鵄的自由剪切層逐漸附著到柱C上，一部分阻力也隨之作用于柱C上.

（3）當(dāng)間距比介于3.0到5.0之間時(shí)，柱A的阻力系數(shù)又大幅度回升，柱C的阻力系數(shù)大幅度減小.這是因?yàn)榇碎g距比下，流場處于渦沖擊模式，每個(gè)方柱都受到了來流沖擊.

（4）柱A和柱D的平均阻力系數(shù)值要比柱B和柱C的大很多.這是因?yàn)樵诹庑闻帕邢?，柱A和柱D的上游為均勻來流，其受到的干擾相對較小.

（5）隨著間距比的增大，柱A和柱D的平均升力系數(shù)一直比較穩(wěn)定，接近于零.這是因?yàn)榱庑闻挪枷轮鵄和柱D上下對稱分布，所受升力相互抵消.

（6）結(jié)合繞流流場和升阻力系數(shù)變化規(guī)律可知，臨界間距比在3左右.

2.4 渦脫落頻率特性

本文計(jì)算結(jié)果表明，布置1時(shí)，上下游柱體的漩渦脫落同步，因此St數(shù)相同.這與相關(guān)文獻(xiàn)的計(jì)算結(jié)果一致[4，6].而當(dāng)間距比較小時(shí)，四個(gè)方柱的渦脫落頻率存在微小差異.St數(shù)隨間距比變化的規(guī)律折線如圖7所示.

圖7 布置1 St隨g變化曲線

由圖7可以看出：①隨著間距比的增加，St先增大后減小，最后接近單方柱繞流的St數(shù)；②隨著間距比的增大，增加同樣間距比時(shí)，St數(shù)的變化幅度越來越??；③在臨界間距比附近，St數(shù)的變化幅度較大.

圖8為布置2時(shí)，St數(shù)隨間距比變化的規(guī)律折線.由圖8可以看出：①隨著間距比的增加，St先減小后逐漸增大，最后接近單方柱繞流的St數(shù)；②在臨界間距比之前，四個(gè)方柱的漩渦脫落同步，St相同，在臨界間距比之后，流場渦脫落頻率發(fā)生較大變化，四根方柱St出現(xiàn)較大差異，柱A和柱D的St近似，都增加幅度較大后再趨于一穩(wěn)定值，柱B和柱C的St近似，并平緩趨于一穩(wěn)定值；③隨著間距比的增大，增加同樣間距比時(shí)，St數(shù)的改變幅度越來越??；④在臨界間距比附近，St數(shù)的改變幅度較大，臨界間距比之后漩渦脫落不同步.

圖8 布置2 St隨g變化曲線

2.5 臨界間距比的驗(yàn)證

綜上所述，臨界間距比是四等邊方柱的繞流模擬中極為重要的一個(gè)參數(shù)，漩渦的脫落，阻力系數(shù)，升力系數(shù)等參數(shù)都在跨過臨界間距比時(shí)發(fā)生改變.現(xiàn)對布置 1，取 g=2.6，2.7，2.8，2.9；布置 2 取 g=2.8，2.9，3.1，3.2進(jìn)行模擬分析，旨在能夠較為準(zhǔn)確地捕捉臨界間距比的大小.

對于布置1，通過進(jìn)一步模擬得到在臨界間距比附近，平均阻力系數(shù)隨間距比的變化規(guī)律，如圖9所示.

圖9 布置1 Cd隨g變化曲線

由圖9可見，正方形排列四等邊方柱繞流的臨界間距比為2.8.當(dāng)間距比超過2.8時(shí)，開始在上游柱體后發(fā)生漩渦的脫落，上下游柱體的阻力系數(shù)都發(fā)生躍升.其中，對于下游方柱，超過臨界間距比后，阻力系數(shù)大約增加了4.1倍.

對于布置2，通過進(jìn)一步模擬得到在臨界間距比附近，平均阻力系數(shù)隨間距比的變化規(guī)律，如圖10所示.

由圖10可見，菱形排列四等邊方柱繞流的臨界間距比為3.0.當(dāng)間距比等于3.0時(shí)，流場特性發(fā)生較大變化，柱D的阻力系數(shù)達(dá)到最小值，而柱C的阻力系數(shù)達(dá)到最大值.其中，對于下游方柱，超過臨界間距比后，阻力系數(shù)大約增加了1.8倍.

圖10 布置2 Cd隨間距比g變化曲線

3 結(jié)論

本文采用有限體方法，對正方形排列和菱形排列下的四方柱繞流問題進(jìn)行了數(shù)值模擬研究.重點(diǎn)分析了漩渦的脫落規(guī)律、升阻力系數(shù)隨著間距比的變化規(guī)律.得出的主要結(jié)論如下.

（1）四方柱繞流中也存在著臨界間距比現(xiàn)象.并準(zhǔn)確地得到了臨界間距比的數(shù)值.當(dāng)Re=100時(shí)，布置1情形的臨界間距比為2.8，布置2情形的臨界間距比為3.0.

（2）對于布置1，在臨界間距比之前，下游柱體的阻力系數(shù)隨間距比增加而增加，但幅度較??；當(dāng)達(dá)到臨界間距比時(shí)，下游柱體的阻力系數(shù)出現(xiàn)了陡然的增加.在達(dá)到臨界間距比后，下游柱體的阻力系數(shù)開始逐步減小，最終趨近于單柱體繞流的結(jié)果.對于布置2，四個(gè)方柱的阻力系數(shù)總體變化不大，只是在臨界間距比附近出現(xiàn)了大幅度的變化.兩種布置下，升力系數(shù)折線均關(guān)于x=0軸線呈對稱分布，并隨間距比的變化而逐漸趨向于零.

（3）數(shù)值模擬中也觀測到了漩渦的同相位和反相位脫落現(xiàn)象，這與實(shí)驗(yàn)是相符合的.對于布置1，四方柱漩渦脫落頻率相同.對于布置2，四方柱在臨界間距比之前漩渦脫落頻率相同，在臨界間距比之后，漩渦脫落頻率出現(xiàn)較大差異.此外，當(dāng)間距比較小時(shí)，流動趨向于不穩(wěn)定，易于發(fā)生相位之間的轉(zhuǎn)換；而當(dāng)間距比增加到3以上時(shí)，相位轉(zhuǎn)換則較少發(fā)生.

（4）St數(shù)隨間距比的增加而增加，且在臨界間距比附近時(shí)，其增幅最大.

參考文獻(xiàn)：

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

[2]張愛社，張 陵.等邊布置三圓柱繞流的數(shù)值分析[J].應(yīng)用力學(xué)學(xué)報(bào)，2003，20（1）：31-36.

[3]張鵬飛，蘇中地，廣林端，等.直徑變化對三圓柱繞流的影響[J].水動力學(xué)研究與進(jìn)展，2012，27（5）：554-560.

[4]FARRANT T，TAN M，PRICE W G.A cell boundary element method applied to laminar vortex-shedding from arrays of cylinders invarious arrangements[J].Journal of Fluids and Structures，2001，30（2）：211-236.

[5]李景銀，CHANK T.繞正方形排列的順排的四個(gè)圓柱的流動研究[J].工程熱物理學(xué)報(bào)，2004，25（1）：59-63.

[6]LAM K，ZOU L.Three-dimensional numerical simulations of cross-flow around four cylinders in an in-line square configuration[J].Journal of Fluids and Structures，2010，26（3）：482-502.

[7]SHAMS U.Numerical simulations of cross-flow around[J].International Journal of Mathematical&Physical and Engineering Sciences，2009，3（2）：76-85.

[8]趙心廣.三方柱繞流的大渦模擬及頻譜分析[J].水科學(xué)與工程技術(shù)，2011，28（4）：28-30.

[9]張 偉，戴玉滿.低雷諾數(shù)下等邊布置三方柱繞流的數(shù)值研究[D].天津：天津城建大學(xué)，2015.

[10]FARRANT T.A cell boundary element method applied to laminar vortex-shedding from arrays of cylinders in various arrangements[J].Journal of Fluids and Structures，2000，14：375-402.