圓柱繞流，你真的會算嗎？【轉(zhuǎn)發(fā)】

2017-05-01 by:CAE仿真在線來源:互聯(lián)網(wǎng)

圓柱繞流是流體力學中的經(jīng)典問題。幾乎所有的流體力學教科書在講解不可壓縮無旋流動(勢流)的時候,都以圓柱繞流為例分析平行流、偶極子以及點渦的組合。圓柱繞流的實際應用也非常廣泛,例如橋墩、煙囪、架空電纜、換熱器中的管束、飛機起落架等等。

圓柱繞流的幾何形狀看似簡單,但是其流動圖像隨著雷諾數(shù)(Re=ρUD/μ,其中ρ為流體密度,U為來流速度,D為圓柱直徑,μ為流體的動力粘性系數(shù))的變化非常豐富,至今還有許多機理并未被人們認識清楚。圖1是不同Re數(shù)下的圓柱繞流圖像。圖2是圓柱繞流阻力系數(shù)(CD=2F/ρU2A,其中F為阻力,A為迎風面積。迎風面積等于圓柱的直徑乘以圓柱的長度)隨Re數(shù)的變化情況。

越簡單=越容易?圓柱繞流,你真的會算嗎?【轉(zhuǎn)發(fā)】fluent流體分析圖片1

圖1 不同Re數(shù)下的圓柱繞流圖像。(復制自文獻[1])

越簡單=越容易?圓柱繞流,你真的會算嗎?【轉(zhuǎn)發(fā)】fluent流體分析圖片2

圖2 圓柱繞流阻力系數(shù)隨Re數(shù)的變化。實驗結(jié)果。(復制自文獻[1])

由于流動的復雜性,圓柱繞流的數(shù)值模擬是非常具有挑戰(zhàn)性的課題。在最近出版的一本流體力學教科書[2]中,作者展示了圓柱繞流模擬的三個案例。使用的程序是商業(yè)軟件FLUENT,模擬的雷諾數(shù)分別是150,104和107。使用二維求解器計算;由于圓柱的幾何對稱性,將計算域設定為一個半圓形區(qū)域(圖3);使用了三種不同密度的網(wǎng)格:粗網(wǎng)格(徑向30×圓周方向60)、中等網(wǎng)格(60×120)以及細網(wǎng)格(120×240)。圖4是中等網(wǎng)格的局部放大圖。三個雷諾數(shù)下都使用定常算法求解。

越簡單=越容易?圓柱繞流,你真的會算嗎?【轉(zhuǎn)發(fā)】fluent流體分析圖片3

圖3 計算域(復制自文獻[2])

越簡單=越容易?圓柱繞流,你真的會算嗎?【轉(zhuǎn)發(fā)】fluent分析圖片4

圖4 中等網(wǎng)格(復制自文獻[2])

當Re=150時,對于實驗測得的阻力系數(shù),不同文獻的結(jié)果有一定差異,范圍在1.1~1.4之間。實驗表明這時圓柱表面的邊界層是層流的,分離點位于82°的位置(圖5(a)),即位于圓柱最高點的前方;尾跡也是層流的。

越簡單=越容易?圓柱繞流,你真的會算嗎?【轉(zhuǎn)發(fā)】fluent分析圖片5

圖5 分離點的位置。(a) Re=150和104; (b) Re=107

對Re=150使用層流模型計算。圖6是三種不同密度網(wǎng)格的計算結(jié)果(流線圖);三種密度的網(wǎng)格的計算收斂性都很好。對于圖6(a)的粗網(wǎng)格的結(jié)果,邊界層的分離點位于圓柱最高點的后方,在120°的位置,阻力系數(shù)為1.00。對于圖6(b)的中等網(wǎng)格的計算結(jié)果,分離點稍向上游移動,位于110°的位置,似乎和實驗接近了一點,但是算出的阻力系數(shù)卻降低到0.982——與實驗值的距離更遠了。對于圖6(c)的細網(wǎng)格的結(jié)果,外觀上與中等網(wǎng)格的結(jié)果差異不大,分離點位于109°的位置,而阻力系數(shù)更小(0.977)。再加密網(wǎng)格,結(jié)果的變化就非常小了,因此可以認為已經(jīng)達到網(wǎng)格無關。

越簡單=越容易?圓柱繞流,你真的會算嗎?【轉(zhuǎn)發(fā)】fluent分析圖片6

圖6 Re=150的流線。(a) 粗網(wǎng)格;(b) 中等網(wǎng)格;(c) 細網(wǎng)格

(復制自文獻[2])

顯然,數(shù)值模擬結(jié)果并不符合實際。實驗表明分離點位于圓柱最高點前方,在82°的位置,而數(shù)值模擬的結(jié)果中,分離點位于圓柱最高點后方,在約109°的位置。此外,數(shù)值模擬算出的阻力系數(shù)也比實驗值小。

對于Re=104,實驗結(jié)果表明圓柱前部的附著邊界層是層流的,分離點也位于82°的位置,但是,圓柱的尾跡是湍流的;實驗測得的阻力系數(shù)為1.15。這個雷諾數(shù)下使用k-ε湍流模型計算;入口的湍流度設為10%,湍流長度尺度設為圓柱直徑的一半。圖7是三種不同密度網(wǎng)格的計算結(jié)果(流線圖);三種密度的網(wǎng)格的計算收斂性都很好。對于圖7(a)的粗網(wǎng)格的結(jié)果,邊界層分離點在圓柱頂部后方很遠的地方(140°處),算出的阻力系數(shù)只有0.647,幾乎比實驗值小了一半。對于圖7(b)的中等網(wǎng)格的結(jié)果,分離點前移到104°的位置,阻力系數(shù)增加到0.742——雖然離實驗值近了一些,但是仍然顯著小于實驗值。對于圖7(c)的細網(wǎng)格的結(jié)果,與中等網(wǎng)格的結(jié)果非常相似,阻力系數(shù)為0.753,分離點位于102°的位置?？梢钥闯?這次數(shù)值模擬結(jié)果仍然不符合實際;分離點的位置與實驗不符,阻力系數(shù)也顯著小于實驗值。

越簡單=越容易?圓柱繞流,你真的會算嗎?【轉(zhuǎn)發(fā)】fluent分析圖片7

圖7 Re=104的流線。 (a) 粗網(wǎng)格;(b) 中等網(wǎng)格;(c) 細網(wǎng)格

(復制自文獻[2])

對于Re=107,實驗結(jié)果表明圓柱表面的附著邊界層是湍流的(頭部駐點附近除外),尾跡當然也是湍流的;分離點位于圓柱最高點后方,在120°的位置(圖5(b));阻力系數(shù)為0.7。這個雷諾數(shù)的流動同樣使用k-ε湍流模型計算。網(wǎng)格的密度與前述的細網(wǎng)格相似,但是在壁面附近按照湍流模擬中近壁面網(wǎng)格的要求進行加密(因為這時雷諾數(shù)很高,邊界層很薄);與Re=104時一樣,入口的湍流度設為10%,湍流長度尺度設為圓柱直徑的一半。不幸的是,算出的阻力系數(shù)只有0.262——比實驗值的一半還要小。按計算結(jié)果畫出的流線如圖8所示;算出的分離點位于129°的位置,雖然與實驗值同處于圓柱最高點后方,但是與實驗值相比要往下游一些。

越簡單=越容易?圓柱繞流,你真的會算嗎?【轉(zhuǎn)發(fā)】fluent圖片8

圖8 Re=107的流線。數(shù)值模擬的結(jié)果。(復制自文獻[2])

Re=150,104和107這三種情形分別代表了圓柱分離流動的三種最典型的流態(tài),即“層流邊界層—層流尾跡”、“層流邊界層—湍流尾跡”以及“湍流邊界層—湍流尾跡”。不幸的是,三種流態(tài)的數(shù)值模擬都以徹底失敗告終。

實際上,對于Re=150的情形,圓柱尾跡中存在非定常的卡門渦街,附著渦交替地從圓柱的兩側(cè)脫落。而在計算模型中,根據(jù)圓柱的幾何對稱性天真地將計算域設為實際流動區(qū)域的一半,并在x軸上應用了對稱邊界條件;再者,計算過程中使用了定常的求解器;這些都是違背實際發(fā)生的物理機制的。計算結(jié)果表明,如果將計算域擴展到整個流動區(qū)域,并且使用非定常求解算法,就可以算出非定常的卡門渦街來,對算得的非定常阻力系數(shù)取時間平均之后,約為1.14,這時就落在實驗值的范圍(1.1~1.4)里面了。計算值接近實驗值范圍下限的原因可能是計算時采用了二維模型,使得任何傾斜的渦脫落或者其它三維不均勻性都被抑制了。

對于Re=104的情形,我們?yōu)榱四M湍流尾跡,使用了湍流模型進行計算。但是,實驗結(jié)果表明,圓柱表面的邊界層是層流的,這是計算結(jié)果不符合實際的重要原因。在我們的計算結(jié)果中,邊界層的分離點位于>90°的位置;實際上只有當圓柱表面為湍流邊界層時,分離點才會移動到這樣的位置?？傊?層流/湍流混合流動對于數(shù)值模擬來說,仍然是非常困難的問題。另外,實驗表明這個雷諾數(shù)下尾跡是非定常的,而我們使用定常求解算法并在x軸使用了對稱邊界條件,這也是不符合實際的。

對于Re=107的情形,實際上圓柱下游存在湍流的卡門渦街,而我們使用定常求解算法并在x軸使用了對稱邊界條件,這是不符合實際的。

對Re=104這一類“層流邊界層—湍流尾跡”的流態(tài)以及Re=107這一類“湍流邊界層—湍流尾跡”的流態(tài),要采用更細致的模擬方法才能得到合理的結(jié)果。例如,文獻[3]采用了二維非定常RANS和三維DES對這兩種流態(tài)進行了模擬(見注1),其結(jié)果表明,對于后一種流態(tài),兩種方法都能得到比較合理的結(jié)果;而對于前一種流態(tài),只有三維DES才能得到合理的結(jié)果。

順便說一句,從圖2可以看出,Re=2.5×105時阻力系數(shù)突然下降。層流邊界層分離點位置與湍流邊界層分離點位置的截然不同,正是造成這時阻力系數(shù)突然下降的原因。當雷諾數(shù)超過此值時,邊界層從層流轉(zhuǎn)捩為湍流,由于湍流邊界層對抗分離的能力比層流的強,所以分離點突然從<90°的位置移動到>90°的位置,使得尾跡顯著變窄,壓差阻力大大降低。這就是著名的圓柱繞流“阻力危機”(drag crisis)現(xiàn)象。所謂的危機并不是流動現(xiàn)象有什么很壞的影響,而是說這個現(xiàn)象使傳統(tǒng)的流體力學理論受到了挑戰(zhàn),學術界出現(xiàn)了危機[4]。

從這幾個例子還可以看出,要正確地模擬流動問題,必須對流動的物理機制有一定的認識;不懂物理機制而盲目地進行計算,常常會得到錯誤的結(jié)果。

未經(jīng)許可,不得轉(zhuǎn)載

參考文獻

[1]莊禮賢, 尹協(xié)遠, 馬暉揚. 流體力學. 第二版. 中國科學技術大學出版社,2009

[2] Yunus A. ?engel, John M. Cimbala. Fluid Mechanics: fundamental and applications. McGraw-Hill, 2014

[3] ANDREI TRAVIN, MICHAEL SHUR, MICHAEL STRELETS and PHILIPPE SPALART. Detached-eddy simulations Past a Circular Cylinder. Flow, Turbulence and Combustion, 63, 1999

[4] 王洪偉. 我說理解的流體力學. 北京:國防工業(yè)出版社,2014

注1:RANS是“雷諾平均方法”(Reynolds averaged Navier-Stokes),即我們在商業(yè)軟件中常用的湍流模型;DES是“脫體渦數(shù)值模擬”,即在流動附著于壁面的區(qū)域采用RANS,在其它區(qū)域自動轉(zhuǎn)變?yōu)榇鬁u模擬。

轉(zhuǎn)自公眾號:流體那些事兒

開放分享：優(yōu)質(zhì)有限元技術文章,助你自學成才