Fluent流固耦合基礎教程（上）

2016-09-04 by:CAE仿真在線來源:互聯網

最近用fluent UDF做了個彈性梁的流固耦合問題。這里把經驗跟大家分享一下。這個帖子主要是介紹基本操作和程序編寫,算是掃盲貼;對于CFD,有限元、梁理論、流固耦合原理不做過多探討,請大家見諒。另外作者水平有限,時間精力也有限,寫得不好也請見諒了。

作者假設讀者熟悉Fluent,對CFD的基本概念有所了解,具備一定的C/Fortran語言編程技能,并且了解有限元的基本理論和程序編寫。不過即使背景有所欠缺也無妨,點是了解基本流程和操作方法。

1. 流固耦合問題的種類

流固耦合問題可以分為很多種。按照耦合程度分類,可以分為強耦合問題和弱耦合問題。但是這種分類方法有兩個流行的版本。

第一版本根據流體和固體之間相互影響的程度來劃分,在工業應用中談得比較多。如果固體運動很小,對流場的影響不大,則認為是弱耦合。比如金屬管道內的水流引起的管壁運動,機翼的振顫等,可以屬于這類。在計算機技術不夠發達的時代,這類問題一般采用簡化的流場模型,如簡單幾何形狀的理想流體,解析出由于固體運動而引起的流場作用力,然后再將這些作用力施加到固體控制方程上,從而得到附加質量,附加粘度,附加剛度,以及其他非線性項。隨著計算機技術的發展,對這類問題的直接數值模擬逐漸成為可能。如果固體變形或運動比較大,其對流場的影響不可忽略,則屬于強耦合。強耦合問題必須對流場和固體進行細致的計算,充分考慮固體變形對流場的影響。閥門開閉,血液通過血管瓣膜,旗幟在風中飄舞,等都屬于這類問題。解決這類問題需要依靠數值模擬,傳統的簡化模型很難應用在這類問題上。

第二個版本根據流體和固體求解的模式來劃分,在學術領域用的比較多。強耦合問題被認為是流場和固體變形必須同時求解的問題。帶有滲流的多孔介質變形問題可以算這類問題。弱耦合問題則是流場和固體變形可以分開求解,但是二者的信息交換通過循環迭代完成。多數工程問題都可以算做這類問題。

這里我要講的是第一類劃分標準中的弱耦合問題,但所用的數值方法也可以應用到強耦合問題中,只是結構體的有限元方法需要采用大變形理論。按照第二種劃分標準來說,這里解決的是流場和固體變形/運動分別求解,兩場之間的相互作用通過迭代實現。

2. 流固耦合問題的數值模擬

既然是數值模擬,則需將系統方程在時域上離散。在每個時間步上分別解算流體域和固體域。固體的變形或位移導致流固邊界的運動。這個運動以動網格的形勢傳遞給流場,然后求解流場。流場的解則包括了這個運動所導致的反力。將這個反力反饋給固體,求解其下一時間步的響應。如圖1所示。

圖1

以上的概念恐怕大家都清楚,只是談到實現的時候,需要做很多細致的工作。當然,如果大家有ANSYS14 WorkBench,則流固耦合問題很容易解決,因為大部分的工作都自動化處理了。但是如果只用Fluent,則需要做一些工作。即便使用先進的數值模擬工具,如果能弄清楚整個流程以及其中的一些技巧,也會對提高效率以及保證精度有所幫助。下面就講講具體怎么做。

3. 流固耦合例題

既然是具體講,就不如找個例題來說。圖2是一根柔軟的彈性梁,浸泡在在直徑為Dt的圓截面管道里。梁的截面也是圓形,直徑為D,長度為L。梁的兩端為鉸接。整個結構是三維結構。梁可以在三維空間里彎曲。流體為水,梁中點所在截面上的平均流速為U。管道入口處的流速為Uin。上游長度為Lup,下游長度為Ldw。上下游的長度分別為10D。

圖2

4. 流體模型

流體模型的建立是關鍵的一環,因為流體模型的好壞直接影響到所得解的合理性。建立好的流體模型需要利用流體力學知識對所處理的問題進行綜合分析。根據流體的雷諾數和其它特征建立合適的網格,選擇合理的湍流模型和算法。我們希望流體模型能夠盡量準確地捕捉到梁表面的力,因此準確的near-wall處理方法很重要。采用wall function還是精確地模擬邊界層,需要給出適當的理由。由于我們處理的是動力問題,希望模型對梁表面力的變化給出準確的值。采用動態k-epsilon 或者k-omega湍流模型計算所得的激勵幅值明顯偏小,這是由于雷諾平均造成的。這里我們采用對流體動力特性捕捉比較好的大渦湍流模型(Large Eddy Simulation, LES)。LES對near-wall處理方法要求比較高,需要邊界上的第一層網格厚度滿足y+ ~ 1。在Gambit里可以方便地建立邊界層網格,達到以上要求。如何建立良好的網格是個復雜的話題,這里且不涉及。

由于固體-液體邊界要移動,流體網格必須改變,需要運用動網格技術。此時建模的時候需要將動網格區域和不動網格區域分開,如圖3所示。

圖3

固體和液體的交界面附近的區域是動網格區域(fluid_dyn),遠離固體的區域可以定義為不動網格區域。動網格區域的大小應當視具體問題而定,但起碼要能夠囊括固體的最大位置變化范圍。梁表面的網格節點位置隨著梁的運動而變化,因此需要將梁的表面單獨定義為一個區域 wall_mov_beam。梁的兩個端點也需要進行控制,因此也可以分別定義獨立區域:wall_mov_end1 和 wall_mov_end2。在劃分流體網格的時候,接近固體表面的地方需要建立邊界層。Fluent幫助文檔上說邊界層應該單獨作為一個動網格區域,這是為了便于控制邊界層的變形,防止出現邊界層的異化。實際上如果將邊界層的厚度增加到物理邊界層厚度的兩倍以上,并且只采用smoothing網格控制,則可以不必將邊界層單獨定義為一個區域。但是如果變形較大,流體網格需要重新劃分,則最好將邊界層單獨分離為一個區域,并且在這個區域內使用smoothing控制。對于管壁附近的邊界層,我們并不關心管壁上的壓力和剪力,而且考慮到出了邊界層以后這個區域的局部解對整體解影響不是很大,因此可以適當降低網格要求,甚至采用wall function近似逼近,這樣可以節省不少計算量。下面兩張圖分別是圖4幾何模型(部分)和圖5網格劃分(部分)。

圖4

圖5

流體為不可壓縮流體,求解器采用一次隱式瞬態算法。LES的動網格計算上Fluent 不支持二次格式,這就需要時間步設定要小一些。時間步長和求解器參數的選取也是個復雜的問題,這里也不展開討論了。在這個算例中,我們取較大的步長,以加快計算速度。具體參數如下:

Time step size: 0.0005 s

Solver settings: unsteady / pressure-based / 1st order implicit

Viscous model: LES / Smagorinsky-Lilly / no dynamic stress

Solution settings:

velocity-pressure coupling: SIMPLE

Relaxation factor: pressure 0.3 / momentum 0.7 / density 1.0 / body force 1.0

Discretization: pressure standard / momentumn bounded central difference

在不加入流固耦合的情況下,計算結果收斂很好,穩定以后每個時間步上的循環(iteration)為兩步。每個時間步循環終止時的殘差為:

continuity x-velocity y-velocity z-velocity

4.2946e-04 8.3873e-06 8.1793e-06 1.3922e-04

流體模型解算成功之后就需要考慮動網格和固體變形的問題了。梁表面流體網格節點的運動需要利用UDF(User Defined Function)來控制。這將在下面的章節里說明。這里首先設置動網格的參數。

動網格區域fluid_dyn應該被設置為Deforming,采用的網格控制方法為Smoothing,但是對于復雜結構,可能需要Remeshing。具體選項的含義請參見Fluent 幫助文檔中的用戶手冊第11章(Modeling Flows Using Sliding and Deforming Meshes)。這里只說明大體思路。這個區域的網格將會隨著這個區域的邊界(梁表面wall_mov_beam)變化。我們只要控制wall_mov_beam,則Fluent可以自動計算出這個區域內的流體網格的新位置。參數設置如圖6所示。

圖6

面區域wall_mov_end 1和wall_mov_end2也設置為動網格區域,跟fluid_dyn類似,具體的網格變形由Fluent自動處理。這是由于我們的梁是兩端鉸接的,因此端面上的網格變化不大。對于其他的情況,比如懸臂梁,這兩個面上的點應該跟據梁的運動用UDF來控制,不能交給Fluent自動處理。具體設置參數如圖7所示。

圖7

梁表面網格節點的位移必須通過UDF來控制。因此在選擇的時候應當選User-Defined。在Motion Attributes里面選擇相應的UDF。這一點在下面的章節中說明。這里先給出Meshing Options的值,參看圖8。這里提醒讀者,wall_mov_beam是驅動網格變化的“源頭”。我們通過UDF控制wall_mov_beam的運動,而周邊的網格將隨著wall_mov_beam的運動而改變。