為何我這個流動總是算不收斂？我要砸電腦！

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

當(dāng)流動具有強(qiáng)烈的分離特性的時候,如果用定常算法計算不收斂(殘差不能降下去,monitor監(jiān)視的物理量總是在振蕩),可以改用非定常算法來計算,通常就能在每一個時間步實現(xiàn)收斂。

很多FLUENT的用戶都會遇到一個令人困惑的問題,那就是用定常流動算法計算的時候,總是不能收斂得很好:殘差總是不能降低到默認(rèn)的收斂標(biāo)準(zhǔn)(1e-03)以下(當(dāng)然,有時候殘差并不能代表收斂);更重要的是,用監(jiān)視器(monitor)監(jiān)視某些量隨著迭代次數(shù)增加的變化情況(例如某個點的速度隨迭代次數(shù)增加的變化)時,曲線總是在振蕩,不能收斂到一個確定的數(shù)值。

我們來看一個例子。如圖1所示,一個邊長為2m的立方體連接著兩根方形截面的管子。入口流速是1m/s,假定流體的密度是1kg/m3,粘性系數(shù)是10-5Pa?s。以入口流速為特征速度,入口截面的邊長為特征長度,則流動雷諾數(shù)為105,所以流態(tài)是湍流流動。坐標(biāo)原點位于立方體的中心,坐標(biāo)軸的方向如圖所示。

圖1 問題描述

網(wǎng)格使用ICEM CFD生成,為分塊結(jié)構(gòu)化網(wǎng)格(圖2),網(wǎng)格數(shù)量約37萬。對靠近壁面的網(wǎng)格加密,壁面第一層網(wǎng)格的高度為0.01m。這樣的網(wǎng)格高度,使各處壁面第一層網(wǎng)格形心的y+大約在1~60的范圍內(nèi)(通過計算結(jié)果的后處理可以看出來)。因為我們使用了FLUENT中的k-ωSST湍流模型來計算,這個模型在默認(rèn)的情況下在近壁面區(qū)域近使用了增強(qiáng)壁面處理(enhanced wall treatment),所以上述網(wǎng)格高度是合適的。

圖2 網(wǎng)格

使用FLUENT中的基于壓力的求解器計算,速度和壓力的耦合(Pressure-VelocityCoupling)采用SIMPLE算法。入口使用velocity-inlet邊界條件,流速按上述條件設(shè)置。入口的湍流參數(shù)設(shè)為:湍流強(qiáng)度=5%,水力直徑=1m。出口使用pressure-outlet邊界條件,回流的湍流參數(shù)設(shè)為與入口相同。計算前的初始流場按照入口邊界的數(shù)值初始化。

用定常算法(Steady)來計算,可以發(fā)現(xiàn)收斂得很不好。從圖3的殘差曲線看,所有的方程的殘差都不能降低到默認(rèn)的收斂標(biāo)準(zhǔn)(1e-03)以下;而使用monitor監(jiān)測點(-0.5,0.5,0.5)處的流速隨迭代次數(shù)的變化情況,其曲線振蕩得非常厲害(圖4)。調(diào)節(jié)松弛因子(Under-relaxation factors)也不能奏效。出口邊界也并未出現(xiàn)回流(reversed flow)的現(xiàn)象,所以也不是導(dǎo)致不收斂的原因。

圖3殘差曲線

圖4(-0.5,0.5,0.5)處的速度隨迭代次數(shù)的變化情況

有什么解決辦法呢?把算法改成非定常的就行了。這里我們把時間步長設(shè)為0.1s(見注1)。由圖5可以看出,非定常計算在每一個時間步都可以將殘差迭代到1e-03以下。在非定常算法下,點(-0.5,0.5,0.5)處的速度是隨著時間的推移而不斷變化的(如圖5所示)。

圖5用非定常算法計算

在這個例子中,從入口管道到立方體是一個截面突然擴(kuò)大的流動,存在強(qiáng)烈的分離流動。經(jīng)驗表明,當(dāng)流動具有強(qiáng)烈的分離特性的時候,使用定常算法往往就無法計算收斂。來自波音公司的湍流模型大師Spalart(著名的SA模型的提出者)在一篇關(guān)于湍流模擬策略的文章[1]中也提到了非定常模擬的問題。他在文章中寫道:

“More frequently, a user that is after massive separation will find that the code simply cannot find a steady state, and that the only course is to operate in a time-accurate mode and analyse an unsteady solution.”(更為常見的情形是,當(dāng)計算具有強(qiáng)烈分離特性的流動時,計算程序無法算出一個穩(wěn)定的解,所以唯一的出路就是使用非定常算法來計算。)

另一個典型的例子是鈍體繞流;高雷諾數(shù)鈍體繞流伴隨著強(qiáng)烈的流動分離,常常存在大尺度的非定常現(xiàn)象。圖6顯示了三維方柱繞流的非定常計算結(jié)果(垂直于方柱的某個截面的速度云圖)。該算例的雷諾數(shù)為105(特征長度取柱的截面的邊長)。所使用的湍流模型也是k-ωSST。可以看出方柱兩側(cè)的交替的、周期性的渦脫落現(xiàn)象(卡門渦街)。對于這個算例,使用定常算法的時候,與前面的例子一樣,是不能算收斂的。這個例子中,使用非定常算法的物理基礎(chǔ)非常明顯:流動中存在大尺度的非定常現(xiàn)象——周期性渦脫落現(xiàn)象。

圖6 方柱繞流

實際上,對于定常算法不能收斂的題,改用非定常算法,這也是文獻(xiàn)中常見的做法。例如,文獻(xiàn)[2]在模擬一級方程式賽車的車輪周圍的氣流流動時(見圖7),當(dāng)遇到定常算法不收斂的時候,就使用非定常算法進(jìn)行了計算,如圖8所示。可以看出,在計算車輪CI的時候,對于幾種湍流模型,使用定常算法都得到了收斂的結(jié)果。但是,計算車輪CII的時候,大部分湍流模型在使用定常算法的時候都不收斂,因此使用非定常算法來計算。該文獻(xiàn)得出的結(jié)論是,與其它模型相比較,使用非定常k-ωSST模型的結(jié)果,以及使用大渦模擬(Large eddy simulation,簡稱LES)的結(jié)果,與實驗數(shù)據(jù)吻合得最好。(當(dāng)然,對于非定常算法,與實驗測得的時間平均量比較的時候,應(yīng)當(dāng)對足夠長的物理時間內(nèi)的計算結(jié)果取平均值。)

圖7模擬一級方程式賽車的車輪周圍的氣流流動(復(fù)制自文獻(xiàn)[2])

圖8文獻(xiàn)[2]中使用定常算法和非定常算法的情況(復(fù)制自文獻(xiàn)[2])

湍流的模擬仍然是一個非常開放的問題。如果讀者有什么想法,歡迎聯(lián)系公眾號。

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

參考文獻(xiàn)

[1] P. R. Spalart. Strategies forturbulence modelling and simulations. International Journal of Heat and FluidFlow, 21, pp. 252-263, 2000

[2] John Axerio-Cilies, Emin Issakhanian, Juan Jimenez, Gianluca Iaccarino. AnAerodynamic Investigation of an Isolated Stationary Formula 1 Wheel Assembly.Journal of fluids engineering. 134, 2012

注1:根據(jù)FLUENT的User’s Guide,時間步長應(yīng)該至少比系統(tǒng)中的最小的時間常數(shù)小一個量級;一種判斷方法是看FLUENT在每個時間步之內(nèi)需要迭代多少次才收斂,理想的狀態(tài)是5至10次。從圖5可以看出,FLUENT在每個時間步之內(nèi)迭代5次左右就能收斂,說明所設(shè)定的時間步長是合適的。另外,FLUENT的User’s Guide建議在后處理中繪制Cell Courant Number的云圖,如果非定常運動的區(qū)域的Cell Courant Number不超過20-40,就說明所取的時間步長是合適的。

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

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