不同工況下有限長圓柱繞流實驗測量及分析

梁圣召，周軍偉，梅 蕾

(哈爾濱工業大學（威海）船舶與海洋工程學院，山東 威海 264209)

0 引 言

不同長徑比、倒立入水的有限長圓柱流體流動特性對很多工程應用，包括浮式海上結構和橋梁渦振都有重要的研究價值。因為有限長圓柱的尾流和自由端的流動相當復雜，僅用實驗方法不能滿足相似參數和相似定律的要求，而且各項數據的測量以及尾流和自由端流動狀態的觀察比較困難。實驗的這些不足可以由數值模擬來補充，但是截至目前，研究不同長徑比倒立入水的有限長圓柱繞流流動特性的文獻還不夠充分。

早期的研究以實驗為主，Okamoto 和Yagita[1]通過測量長徑比AR=1-12.5 圓柱表面的壓力分布得到了圓柱平均阻力系數，并分析了自由端的分離現象。H.Akilli和D.Rockwell[2]利用粒子圖像測速法（PIV）對倒立入水的長徑比為AR=1 和AR=2 的有限長圓柱進行了研究，研究結果形象解釋了大尺度漩渦的形成。J.R.Chaplin等[3]在雷諾數Re與傅汝德數Fr的比值恒等于2.79×105的條件下對有自由液面的倒立入水圓柱在拖曳水池中進行了阻力實驗研究。張榮譽[4]實驗研究了具有自由液面的圓柱繞流現象，測量了不同Fr數工況下尾流場的變化規律。

隨著CFD 技術的發展，Lee 等[5]采用CFD 方法（大渦模擬）對2 種長徑比AR=2.5 和10 的壁面固定的有限長圓柱進行了數值模擬研究。Afgan 等[6]也采用大渦模擬在雷諾數Re=2×104和長徑比AR=6 和10 的情況下對壁面固定的有限長圓柱進行了研究。Guilherme F.Rosetti 等[7]在雷諾數Re=4.3×104和傅汝德數Fr=0.31 的條件下，對長徑比ui=2倒立入水的有限長圓柱從實驗和數值模擬兩方面進行了研究。發現自由表面的存在會影響流動，但是占主導作用的還是自由端的影響。M.A.Benitz 等[8]在雷諾數Re=2 900 和傅汝德數Fr=0.65 的條件下，對長徑比從1～19 倒立入水的有限長圓柱進行了實驗研究，并用大渦模擬進行了對比驗證。

目前學者們已經探討了部分雷諾數，傅汝德數以及長徑比下圓柱繞流的流動特性，同時開展了以大渦模擬為主的數值方法驗證。但是Re和Fr的覆蓋范圍還不夠全面，探討還不夠充分。本文首先采用實驗測量和數值模擬相結合的方法對不同長徑比圓柱的阻力系數進行對比分析，而后對有限長圓柱的尾流場特征進行討論，為海工結構設計提供一定的數據參考。

1 數值模擬方法

1.1 控制方程

對于三維非定常不可壓縮粘性流動，其連續性方程和動量方程為：

式中：ui，uj為速度分量；p為壓力；v為總的流體運動黏性系數，表示為層流黏性系數與湍流黏性系數之和，其中湍流黏性系數由湍流模型獲得；fi為外力分量。

本文采用Wilcox 的k-ω湍流模型，具體可參照文獻[9]。基于CFX 求解器實現方程的求解。

1.2 計算域模型

幾何模型如圖1 所示，長×寬×高為4 m×2 m×1.2 m，相當于80D×40D×24D，其中特征長度D=0.05 m。坐標原點在圓柱底面的中心處。計算域入口距離圓柱中心為20D，出口距離圓柱中心60D，兩側面距離圓柱中心20D。

圖1 有限長圓柱繞流幾何模型Fig.1 Model of computational domain for finite cylinder

流體域進口設定為速度進口（inlet）；出口設定為壓力出口邊界（outlet）；流體域上表面以及2 個側面均設定為自由滑移壁面（free slip wall）,下表面設置為無滑移壁面（no slip wall）；圓柱表面包括自由端面設定為無滑移壁面（no slip wall）。

1.3 數值方法驗證

1）分別進行網格數量和湍流模型對計算結果和計算時間的對比分析。計算條件為：流場無窮遠處的均勻來流速度為U=0.6 m/s，流體為水，密度ρ=1 000 kg/m3，運動粘性系數ν=10-6m2/s，計算得雷諾數Re=3×104，時間步長0.001 s。網格在圓柱周向均布100 個節點，法向擠出20 層，線性增長率1.3。近壁面處保證y+＜2。3 種網格密度采用同樣的加密區域7D×3D，采用四面體、棱柱型與六面體網格相結合的方式，以改善網格質量，計算域網格劃分如圖2 所示。3 套網格分別命名為Coarse，Medium 和Fine，網格數分別為49 萬、119 萬和200 萬。表1 為k-ω湍流模型不同網格密度條件下，以及在Medium 網格下分別采用DES，k-ω和k-ε三種湍流模型所得到的升阻力系數以及計算所耗時間。可以看出當采用k-ω模型時，使用Medium 和Fine 網格，Cd計算結果相差2.7%，Cl計算結果相差9.93%，但是計算時間相差1 倍，所以Medium 網格數計算精度足夠。而在Medium 網格數下，k-ε的計算結果與DES 和k-ω的結果相差較大，DES 和k-ω的計算精度最大誤差在15%以內，但是計算時間相差近1 倍。綜合考慮計算精度與計算時間成本，本文選取Medium 網格數119 萬，湍流模型選擇k-ω。

圖2 不同密度的網格劃分Fig.2 Grid with different densities

表1 不同網格密度和湍流模型的比較Tab.1 Numerical results of different mesh refinement and different turbulence model

2）鑒于文獻在研究相關問題時多采用大渦模擬（LES）湍流模型，為了進一步驗證本文計算模型的可用性，在同等初始條件下，對不同工況及長徑比的有限長圓柱進行計算，并與前人所做相近工況的實驗和數值模擬結果進行對比。首先，本文計算和實驗結果與前人所做相近工況的實驗和數值模擬結果進行對比，如表2 所示。可以看出：在長徑比AR=6，雷諾數Re=3×104工況下，本文計算結果比H.Sakamoto[10]的實驗值大6.59%，實驗結果比H.Sakamoto[10]的實驗值小12.3%；在長徑比AR=10，雷諾數Re=3×104工況下本文的計算結果比LEE[5]的計算值大5.18%，實驗值比其小2.56%。與相似工況下的文獻結果的誤差均在可以接受范圍內。同時，在長徑比AR=10，雷諾數Re=3×104工況下，本文的實驗結果與計算結果也進行了對比，如圖3 所示。可知實驗值與計算值的升阻力系數時程圖比較接近，而橫流向的振動頻率約是順流向振動頻率的2 倍，這與王曉聰[11]的研究結論相符，也側面佐證了本文研究的準確性。最終選取k-ω湍流模型進行數值分析。

圖3 實驗與計算升阻力系數時程圖Fig.3 Distribution of the time-dependent lift and drag coefficient for EXP and CFD

表2 不同作者相似工況計算結果Tab.2 Calculation results of similar conditions for different authors

2 實驗測試

實驗在哈爾濱工業大學（威海）海洋工程學院的循環水槽實驗室進行，三維模型如圖4 所示。實驗所用圓柱材質為均質6061 鋁合金，直徑50 mm，圓柱表面經過高速車削處理，以達到較小的粗糙度。實驗中，測力單元處于水面之上，而底板略低于水面，以盡量減小自由液面對測量結果的影響。

圖4 三維模型Fig.4 3D model

實驗來流速度U選擇為0.2，0.4 和0.6 m/s，對應的雷諾數分別是1×104，2×104和3×104。測試3 個長徑比的圓柱，其高度分別為H=150，300 和500 mm，對應長徑比為AR=3，6 和10。實驗中對測力計進行細致的標定，以保證測量結果的準確度。3 個不同的來流速度以及3 個長徑比組合共9 組實驗，每組實驗進行4 次測量，取平均值。

3 結果與分析

3.1 阻力系數分析

將不同長徑比下，數值模擬與實驗研究所得阻力系數的均方根值進行對比，其與雷諾數的關系如圖5所示。研究發現，當長徑比AR=3，6 和10 時，阻力系數均隨流速的增加有1 個逐漸減小的趨勢，這與R.T.Gon?alves[13]的實驗結果相一致，原因是流速越大，自由端對流動的影響越大。但阻力系數總體變化不大，表明在長徑比大于3 的工況下，有限長圓柱繞流的阻力系數趨于平穩，其值在0.61 與0.86 之間。但同時也發現，本文的實驗結果與數值模擬結果存在5%～25%的誤差，其產生的原因有2 個，一是在實驗過程中循環水槽存在槽體共振現象。而實驗中阻力的最大值只有4N 左右，所以輕微的槽體振動也會對實驗結果產生一定的影響。第2 個原因是本文數值模擬采用的是壁面固定式圓柱，而實驗則是采用的“倒立入水”式。為了盡量減小自由液面對實驗結果影響，實驗中在與自由液面相切的位置增加了一塊底板，防止液面波動對圓柱受力產生影響，但仍不能完全避免自由液面的影響。

圖5 雷諾數與阻力系數的關系Fig.5 Relationship between Re and drag coefficient

3.2 流場分析

流速U=0.6 m/s 工況下，長徑比AR=3，6 和10 的三維流線圖如圖6 所示。長徑比AR=3 的圓柱在展向幾乎沒有渦流，也沒有漩渦脫落。而長徑比AR=6 和10 的圓柱在展向則有明顯的渦流以及漩渦脫落現象。同時發現，當長徑比AR=3 時，圓柱后方有極其混亂的湍流，其影響范圍從圓柱自由端一直到固定壁面。分析原因應該是當流體流經圓柱自由端時產生“下洗”現象，當長徑比較小時，“下洗”運動直接沖擊到圓柱的根部，同時也破壞了圓柱后方有規律的漩渦脫落，產生了極其混亂的湍流；當長徑比AR為6 和10 時，圓柱后方均有一個較小的回流區，并且也會產生“下洗”現象，但因為長徑比較大，“下洗”運動只能影響圓柱的自由端，破壞自由端附近的漩渦脫落，而對圓柱中間以及根部沒有影響。而且相比固定壁面處，圓柱中間部分的漩渦脫落更為明顯，說明固定壁面的存在一定程度上抑制了漩渦脫落的產生。即圓柱的長徑比越長，“下洗”運動影響的長度占展向長度的比例就越小，所以隨著長徑比的增加，阻力系數會緩慢增加。而有限長圓柱繞流的阻力系數比無限長圓柱低的主要原因同樣可能是有限長圓柱自由端的存在使得流體可以繞過自由端流動，而在無限長的情況下流體只能從兩側繞著圓柱流動。

圖6 不同長徑比三維流線圖Fig.6 3D streamline with different AR

為了能清晰地得到圓柱后方漩渦脫落的流動，采用Q準則（Q-criterion）對尾渦進行可視化。Q的定義為：

不同長徑比下三維渦對比如圖7 所示。可以看出，當長徑比AR=3 時，漩渦結構沒有震蕩，因為其長徑比較小，所以自由端產生的高渦量控制了整個圓柱體的流動行為；當長徑比AR=6 和10 時，出現了不同的漩渦脫落區域；漩渦在圓柱的中下部周期性震蕩，并且幾乎與圓柱平行；在圓柱的中上部，漩渦的周期性震蕩仍然存在，但是不再與圓柱平行，而是呈現出一定的角度。這種現象可能是由圓柱中上部分不同展向位置處的振幅不同引起的。長徑比AR=10 時不同流速下三維渦對比如圖8 所示。隨著流速的增加，渦黏度最大值也在增加。當從0.2 m/s 增加到0.4 m/s時，渦黏度最大值增加了92.4%，從0.4 m/s 增加到0.6 m/s時，渦黏度最大值增加了12.4%。說明隨著流速的增加，渦黏度最大值的增速在逐漸降低。從圖中也可以看出，隨著流速的增加，圓柱后方不同位置處渦的相互作用也增強了。

圖7 不同長徑比下三維渦粘度圖（Q=7）Fig.7 3D Eddy viscosity with different AR（Q=7）

圖8 長徑比AR=10 時不同流速下三維渦粘度圖（Q=7）Fig.8 3D Eddy viscosity at different velocity when AR=10（Q=7）

4 結 語

本文對雷諾數Re=1×104-3×104和長徑比AR=3，6 和10 的三維有限長圓柱繞流進行了數值模擬和實驗研究，得到如下結論：

1）在所涉及的工況范圍內，采用k-ω模型對三維有限長圓柱繞流數值分析是可行的，計算穩定，耗時少，計算誤差在可以接受范圍內。

2）當長徑比AR=3，6 和10 時，數值模擬和實驗的結果均表明，阻力系數隨流速的增加有一個逐漸減小的趨勢，數值模擬結果穩定在0.8 左右，實驗分析得到的阻力系數落在0.61～0.93 區間內。表明在長徑比大于3 的工況下，有限長圓柱繞流的阻力系數趨于平穩。

3）圓柱自由端的存在使流動具有更明顯的三維特性，自由端的“下洗”作用會繞過自由端沖擊圓柱后方的漩渦脫落，導致阻力系數小于無限長圓柱。在長徑比AR=6 和10 時，出現了不同的漩渦脫落區域。漩渦在圓柱的中下部周期性震蕩，并且幾乎與圓柱平行。在圓柱的中上部，漩渦的周期性震蕩仍然存在，但是不再與圓柱平行，而是呈現出一定的角度。

4）渦黏度最大值隨著流速的增加也在增加。當從0.2 m/s 增加到0.4 m/s 時，渦黏度最大值增加了92.4%，從0.4 m/s 增加到0.6 m/s 時，渦黏度最大值增加了12.4%。