旋轉圓柱繞流流場特性分析

徐一航 陳少松

(南京理工大學能源與動力工程學院,南京 210094)

引言

旋轉導彈在飛行過程中由于攻角α的存在,此時流經彈體的橫流速度為V∞sinα,由于導彈的旋轉從而產生了馬格努斯力.因此,旋轉導彈馬格努斯力的產生從流動機理上分析是旋轉圓柱繞流問題.由于導彈飛行速度高,致使橫流速度高.所以,研究高雷諾數、低轉速比條件下的旋轉圓柱繞流問題對深入分析旋轉導彈馬格努斯力的機理有著重要意義.

為了系統地研究圓柱尾跡層流向湍流的轉變,人們進行了大量的實驗和數值研究[1-9].當圓柱體繞軸旋轉時,由于馬格努斯效應,會使圓柱后方產生非對稱尾流.雷諾數和相對轉速的變化都會對旋轉圓柱的尾流產生影響[10-13].研究發現,當相對轉速α<2 時,旋轉圓柱尾跡的脫落渦街與卡門渦街相似[14-15].當相對轉速α>2 時,旋轉圓柱后方的脫落渦被抑制[16-17].Lam[5]利用粒子圖像測速法對雷諾數Re=3600～ 5000 的旋轉圓柱進行了實驗測量.結果表明,隨著圓柱轉速的增加,尾跡變窄并向一側偏轉.Aljure 等[18]采用直接數值模擬的方法對雷諾數為5000 時不同相對轉速下的旋轉圓柱繞流進行了數值模擬,得出了渦旋停止脫落的相對轉速,分析了相對轉速變化對駐點和分離點的影響.Ray和Ezadi 等[19-20]對雷諾數為5000 到11 800 的旋轉圓柱繞流進行了數值計算和實驗分析,總結出了升力系數和阻力系數的變化規律,分析了旋轉圓柱后方的流場特征.何穎等[21-22]采用大渦模擬方法對高雷諾數下旋轉圓柱繞流問題進行了數值模擬,分析了旋轉對圓柱繞流分流流動的控制.

目前對于旋轉圓柱在高雷諾數、低轉速比條件下的實驗研究相對較少,此外關于旋轉圓柱后方流場速度變化的研究也相對較少.因此本文對雷諾數Re=20 000～ 90 000、相對轉速(ɑ=ωD/(2U∞))ɑ=0～0.72 的旋轉圓柱后方流場進行了實驗測量,分析了旋轉圓柱后方不同剖面處的速度分布規律和湍流度分布規律.采用大渦模擬方法對旋轉圓柱繞流問題進行數值模擬,分析流場變化規律.而后建立理論模型,分析后方渦位置的改變對圓柱升力和尾流區自由剪切層的影響.

1 實驗裝置與測量技術

1.1 風洞裝置

實驗在南京理工大學HG-1 風洞中進行,該風洞為開口環流式風洞.風洞實驗段截面積為700 mm×700 mm、長1.1 m、速度范圍為0.1～ 40 m/s,湍流度0.97%.

1.2 熱線風速儀

采用Dantec StreamLine Pro 型熱線風速儀對旋轉圓柱后方流場進行測量.熱線探頭由一根1.25 mm長、5 μm 直徑的鍍金鎢絲焊接于兩根不銹鋼針尖構成的支架上,并接入電路中.熱線風速儀的采樣頻率為5 kHz,采樣時間設定為2 s,每一組實驗之前對熱線風速儀進行標定.

1.3 移測坐標架控制系統

實驗采用的是奈凱移測坐標架系統,該系統由精密絲桿傳動器、步進電機、步進電機驅動器、驅動控制板、計算機和控制軟件組成.奈凱移測坐標架移測距離為600 mm×600 mm×600 mm,移測精度為0.02 mm,通過軟件控制該坐標架可以實現定步長運動.

1.4 圓柱旋轉機構

通過電機驅動圓柱旋轉,通過遙控器控制電子調速器進而控制圓柱的旋轉速度.圓柱的直徑D=50 mm、長度L=700 mm.

1.5 流場測量方法

如圖1 所示,將熱線風速儀固定在移測坐標架的支壁上,通過軟件控制移測坐標架的移動,從而對旋轉圓柱的后方流場進行準確的測量.熱線風速儀分別對旋轉圓柱后方四處截面的速度進行測量,每個截面到圓柱軸心的距離分別為1D,1.5D,2D和3D.

如圖2 所示,以圓柱軸線為中心向上、下共選取87 個測量點,每兩個測量點之間的距離為2.5 mm.為了保證測量的準確性,在每一組實驗之前對熱線風速儀進行標定,并對表1 中的工況進行了實驗測量.

圖2 旋轉圓柱測量示意圖Fig.2 Schematic diagram of measuring a rotating cylinder

表1 實驗工況Table 1 Experiment cases

2 實驗結果分析

2.1時均速度

圖3 和圖4 給出了雷諾數為20 000 和60 000 時旋轉圓柱后方4 個剖面處無量綱時均速度分布.當雷諾數一定時,旋轉對圓柱的尾流區產生了一定的影響,致使尾流向一側偏離.旋轉圓柱尾跡區域下方速度突變處的位置隨著相對轉速的增加而上移,而旋轉圓柱尾跡區域上方速度突變處的位置隨著相對轉速的增加并未發生太大變化.例如:當Re=20 000,x/D=1 時,隨著相對轉速由0 增大到0.4、再增大到0.72,其時均速度突變位置分別上移了0.25D和0.11D.時均速度的測量是圓柱后方多個周期渦脫落的時均結果,速度突變位置的變化說明了后方渦時均位置的改變.對于不旋轉的圓柱,其后方渦在Y方向上的時均位置是對稱的.旋轉圓柱尾跡區域下方速度突變處位置的上移說明圓柱的旋轉導致了尾跡區域下方渦的上移.而上方速度突變處位置幾乎不變說明上方渦的位置并未發生太大的變化.圖5 和圖6 給出了ɑ=0.2 和ɑ=0.4 時旋轉圓柱后方4 個剖面無量綱時均速度分布,隨著雷諾數的增加,旋轉圓柱尾跡區域下方速度突變處位置有小幅度的下移(在ɑ=0.2,x/D=1.5 處,隨著雷諾數的增加,其時均速度突變位置分別下移0.05D和0.05D),旋轉圓柱尾跡區域上方速度突變處位置變化不明顯.雷諾數增加對尾跡區域的影響隨著x/D的增加表現得越來越明顯.

2.2 湍流強度

圖7 和圖8 給出了不同雷諾數、不同相對轉速下圓柱后方4 個剖面上湍流度的分布.與速度分布相似,圓柱尾流區下方湍流強度突變點隨著相對轉速的增加而上移,尾流區上方湍流強度突變點隨著相對轉速的增加變化很小.隨著雷諾數的增加,旋轉圓柱尾跡區域下方湍流度突變處位置有小幅度的下移,旋轉圓柱尾跡區域上方湍流度突變處位置變化較小.圖9 給出了不同工況下圓柱后方上、下輪廓線湍流強度分布(上、下輪廓線位置如圖2 所示),相對轉速的增加致使下輪廓線的湍流強度減小,而對上輪廓線的影響不大.從圖9(b)中也可以看出,由于圓柱的旋轉,下輪廓線湍流強度隨x/D增加的減小速度要大于上輪廓線.這說明圓柱的旋轉對圓柱后方流場下側的影響較大,而對上側影響相對較小.

圖3 Re=20 000 時旋轉圓柱后方4 個剖面無量綱時均速度分布Fig.3 Dimensionless time averaged velocity distribution of four sections behind the rotating cylinder when Re=20 000

圖4 Re=60 000 時旋轉圓柱后方4 個剖面無量綱時均速度分布Fig.4 Dimensionless time averaged velocity distribution of four sections behind the rotating cylinder when Re=60 000

圖5 ɑ=0.2 時旋轉圓柱后方4 個剖面無量綱時均速度分布Fig.5 Dimensionless time averaged velocity distribution of four sections behind the rotating cylinder when ɑ=0.2

圖6 ɑ=0.4 時旋轉圓柱后方4 個剖面無量綱時均速度分布Fig.6 Dimensionless time averaged velocity distribution of four sections behind the rotating cylinder when ɑ=0.4

圖7 Re=20 000 時旋轉圓柱后方4 個剖面湍流強度Fig.7 Turbulence intensity distribution of four sections behind the rotating cylinder when Re=20 000

圖8 ɑ=0.4 時旋轉圓柱后方4 個剖面湍流強度Fig.8 Turbulence intensity distribution of four sections behind the rotating cylinder when ɑ=0.4

圖9 (a) Re=20 000 和(b) ɑ=0.4 時旋轉圓柱后方上、下輪廓線湍流強度分布Fig.9 Turbulence intensity distribution of upper and lower contour behind the rotating cylinder when (a) Re=20 000 and (b) ɑ=0.4

3 數值模擬

3.1 控制方程

采用大渦模擬方法求解物理空間盒式濾波形式的三維納維?斯托克斯方程

3.2 網格劃分與邊界條件

選取來流方向沿x軸正向,圓柱直徑為D,如圖10 所示.計算區域的長為50D、寬為40D,圓柱中心與上、下邊界的距離為20D,與左、右邊界的距離分別為20D和30D.圓柱長度為πD.圓柱的旋轉是通過滑移網格區域的網格運動實現的,如圖11所示.滑移網格技術是將整個計算域分為兩個或多個單元區域.在計算過程中,單元區域沿著網格分界面相互滑動,而各自區域的網格不會發生變化,區域邊界之間進行數據交換.圓柱的旋轉區域采用旋轉坐標系方法,壁面第一層網格高度由y+～ 1 確定,近壁面網格增長率為1.02.進口邊界條件為速度入口,出口邊界條件為壓力出口,圓柱壁面為無滑移條件,展向設為周期邊界條件.通過對網格進行獨立性驗證,最終選擇的網格數量為1340 萬.計算時間步長由確定,最終選擇時間 步長為5×10?5s,每一次步長內迭代60 次.空間離散采用通量差分分裂格式,時間推進方法采用雙時間步隱式迭代推進方法.

圖10 流動示意圖Fig.10 Schematic diagram

圖11 網格示意圖Fig.11 Grid diagram

3.3 數值計算驗證

圖12 是雷諾數Re=60 000,x/D=1 時不同亞格子模型計算結果對比圖,通過將不同亞格子模型得到的計算結果與實驗值進行對比得出,采用VWLES S-Omega 模型和Smagorinsky-Lilly 模型的計算結果與實驗值最為接近,雖然在最小速度值方面數值模擬得到的結果要高于實驗值,但是其變化趨勢完全吻合.圖13 為雷諾數Re=60 000,x/D=3 時VWLES S-Omega 和Smagorinsky-Lilly 模型的計算結果,兩者與實驗值均吻合較好.相比之下,VWLES S-Omega模型的計算結果更為精確,因此在后續的數值模擬中均采用VWLES S-Omega 亞格子模型進行計算.

圖12 Re=60 000,x/D=1 時不同亞格子模型計算結果對比圖Fig.12 Calculation results of different subgrid model when Re=60 000 and x/D=1

3.4 數值計算分析

3.4.1 旋轉對圓柱流場特性的影響

圖14 和圖15 為Re=20 000,ɑ=0 和ɑ=0.72時圓柱后方三維渦量圖,給出了不同轉速下渦的發展和自由剪切層的差異.當圓柱不旋轉時,圓柱后方的渦脫落位置是上、下對稱的,上下的渦交替生成.但是當圓柱旋轉之后,圓柱后方下側渦的位置明顯上移,且幅度較大.在下方渦的擠壓下,上方渦的位置有小幅度的上移.不旋轉的圓柱其自由剪切層上下對稱,但當圓柱旋轉后,下方的自由剪切層有明顯的上移,上方的自由剪切層位置變化較小.圖16 和圖17 為Re=20 000,ɑ=0 和ɑ=0.78 時圓柱后方煙流圖.與數值模擬結果相同,旋轉后的圓柱下方自由剪切層的位置明顯上移.圖14～ 圖17 的結果說明,圓柱旋轉后下方自由剪切層的上移導致了圓柱后方流場下側速度突變位置的上移.圓柱后方流場上側速度突變位置基本不變也是因為圓柱后側上方自由剪切層位置變化不大.

圖13 Re=60 000,x/D=3 時兩種亞格子模型計算結果對比圖Fig.13 Calculation results of different subgrid model when Re=60 000 and x/D=3

圖18 為Re=20 000,ɑ=0.72 時圓柱周圍速度矢量圖,上、下表面附近的流場流動有所不同.由于圓柱的旋轉,圓柱下表面的運動方向與主流方向相同,甚至下表面的運動速度要大于主流的運動速度.因此,圓柱下表面會帶動周圍的流體向右上方移動,這導致了圓柱后方下側渦位置和自由剪切層的上移.

圖14 Re=20 000,ɑ=0 時圓柱后方渦量圖Fig.14 Vorticity diagram behind cylinder at Re=20 000,ɑ=0

圖15 Re=20 000,ɑ=0.72 時圓柱后方渦量圖Fig.15 Vorticity diagram behind cylinder at Re=20 000,ɑ=0.72

圖16 Re=20 000,ɑ=0 時圓柱后方煙流圖Fig.16 Smoke plume behind cylinder at Re=20 000,ɑ=0

圓柱上表面的運動方向與主流方向相反,這就使得上表面附近的流體運動方向與主流相反.但是上表面并沒有像下表面一樣對尾跡區產生嚴重的影響.這是因為當雷諾數較高而相對轉速較低時,上表面帶動附近流體向著與主流相反方向運動的現象被較大的主流速度所覆蓋.但是下表面則不同,它帶動著流體順著主流方向運動,最終由于逆壓梯度的作用導致流體無法繼續隨壁面運動而產生分離渦,這就使得下表面的分離渦上移,也導致了自由剪切層的上移.

3.4.2 雷諾數對圓柱流場特性的影響

從圖19 中可以看出,由于圓柱的旋轉,致使圓柱上、下表面時均壓力系數分布不對稱,進而產生了馬格努斯力.雷諾數的變化對圓柱表面時均壓力分布影響較小,只是對上、下表面最小壓力系數的極值產生了一定的影響.從表面時均摩擦力系數分布可知,隨著雷諾數的增加,上表面的流動分離點后移,下表面的流動分離點也后移.但是,上表面的流動分離點位置的變化要比下表面明顯.

圖17 Re=20 000,ɑ=0.78 時圓柱后方煙流圖Fig.17 Smoke plume behind cylinder at Re=20 000,ɑ=0.78

圖18 Re=20 000,ɑ=0.72 時速度矢量圖Fig.18 Velocity vector diagram at Re=20 000,ɑ=0.72

圖19 a=0.4 時不同雷諾數下圓柱表面(a)時均壓力系數和(b)時均摩擦力系數分布Fig.19 Distribution of (a) the average pressure coefficient and (b) the average wall shear stress coefficient on the cylinder surface under different Reynolds numbers when ɑ=0.4

4 理論模型

對于繞圓柱的有環量流動,傳統上采用均勻直線流動、偶極流、點渦三者疊加后的復合流動來表示,其速度勢函數為

其中r2=x2+y2,θ=arctan(y/x).

考慮一個有限的點渦列(圖20 所示),每一個渦的強度都是K,沿一條直線分布在點0,a,2a,···,na處,2n+1 個渦的復勢為

圖20 單渦列Fig.20 Single vortex column

對數項的有限和等于一個正弦項和一個對數項的積

因此其流函數為

在流場中任一點處速度為

假設有兩個平行的點渦列,兩個點渦列之間的距離是b,同一列中連續的兩個渦之間的距離為a,兩個點渦列中渦的強度分別為+K和?K.假設參考點渦分別位于z1和z2,其復勢為

因此,將繞圓柱的有環量流動與兩個平行的點渦列進行疊加,就得到了如圖21 所示的流動分布.它可以最大程度上還原圓柱背風區域的流動特征.圓柱后方流場中任一點處速度為

圖21 帶有兩個平行點渦列的有環量圓柱繞流Fig.21 Circular cylinder with two parallel vortex column

其中,x1,x2,y1,y2的位置由Foppl 等[23-25]提出的圖像法獲得,,經過與實驗值對比,取K1=K2=0.0001Re?0.5.

圖22 和圖23 分別為雷諾數Re=20 000 時旋轉圓柱后方x/D=1 處理論模型與實驗結果對比結果和相對轉速ɑ=0.2 時旋轉圓柱后方x/D=1 處理論模型與實驗結果對比結果.從圖中可以看出,理論模型的計算結果與實驗結果吻合良好,能夠較為準確地預測相對轉速變化時旋轉圓柱后方速度突變的位置,證明了理論模型的正確性.

由伯努利方程求得圓柱表面的壓強分布

圓柱在Y方向上的受力為

圖22 Re=20 000 時旋轉圓柱后方x/D=1 處理論模型與實驗結果對比Fig.22 Comparison of theoretical and experimental results at x/D=1 behind a rotating cylinder when Re=20 000

圖23 ɑ=0.2 時旋轉圓柱后方x/D=1 處理論模型與實驗結果對比Fig.23 Comparison of theoretical and experimental results at x/D=1 behind a rotating cylinder when ɑ=0.2

式中ds=Rdθ,為單位長度圓柱側表面的微面積.圖24是相對轉速ɑ=0.4 時理論模型和LES 升力計算結果對比.隨著雷諾數的增加旋轉圓柱的升力逐漸增加.理論模型的計算結果與LES 的計算結果相對吻合,其變化趨勢相同,但是兩者之間存在一定的差別.LES 計算得到的升力隨雷諾數變化曲線越來越陡,而理論模型計算得到的升力隨雷諾數變化曲線要比其平緩.通過分析渦量圖可以得出如下兩個原因是造成理論模型與LES 存在一定差距的原因.首先,隨著相對轉速和雷諾數的增加,圓柱后方的分離渦是向上方偏移的,而在理論模型中圓柱后方的分離渦是對稱分布的.其次,圓柱后方上、下兩列分離渦的強度是不同的,而在理論模型中圓柱后方的分離渦強度是相等的.

圖24 ɑ=0.4 時理論模型和LES 升力計算結果對比Fig.24 Comparison between theoretical and LES lift calculation results when ɑ=0.4

為了判斷旋轉圓柱后方渦的位置對旋轉圓柱升力的影響,采用理論模型進行了如圖25 所示的靈敏度分析,其中Δy為正表示渦向上移動.可以看出當圓柱逆時針旋轉時,圓柱后側下方點渦列位置的改變對升力的影響顯著,升力的變化與下方點渦列位置的改變呈負相關,而上方點渦列位置的改變對升力的影響相對較小.

圖25 渦位置靈敏度分析Fig.25 Sensitivity analysis of vortex position

通過理論分析發現,圓柱后側下方渦位置的上移對圓柱升力影響十分顯著,圖24 所示數值模擬值要小于理論值,而理論值恰好是隨著下方渦的上移而減小.根據數值模擬的結果可知,圓柱后側下方渦的位置明顯向上移動.因此,足以說明在高雷諾數、低相對轉速的條件下,旋轉圓柱后側下方渦位置的改變對旋轉圓柱的升力、尾流區自由剪切層的變化起到了重要的影響.

5 結論

本文采用熱線風速儀對雷諾數為20 000～ 90 000的旋轉圓柱后方流場進行了測量,通過大渦模擬和理論建模對旋轉圓柱繞流問題進行了分析,得出如下結論:

(1)當圓柱逆時針旋轉時,同一雷諾數下隨著相對轉速的增加,旋轉圓柱尾跡區域下方速度突變處的位置隨著相對轉速的增加而上移,而圓柱尾跡區域上方速度突變處的位置隨著相對轉速的增加并未發生太大變化.同一相對轉速下雷諾數的增加使旋轉圓柱尾跡區域下方速度突變處位置有小幅度的下移,旋轉圓柱尾跡區域上方速度突變處位置變化不明顯.

(2)通過數值模擬發現,圓柱旋轉之后,圓柱后方下側渦的位置明顯上移,且幅度較大.在下方渦的擠壓下,上方渦的位置有小幅度的上移.不旋轉的圓柱其自由剪切層上下對稱,但當圓柱旋轉后,下方的自由剪切層有明顯的上移,上方的自由剪切層位置變化較小.

(3)通過理論分析發現,圓柱后側下方渦位置的上移對圓柱升力影響十分顯著,在高雷諾數、低相對轉速的條件下,旋轉圓柱后側下方渦位置的改變對旋轉圓柱的升力、尾流區自由剪切層的變化起到了重要的影響.