FLUENT算例 —— Turbulent Pipe Flow (LES) 圆管湍流流动(大涡模拟)
來自https://confluence.cornell.edu/pages/viewpage.action?pageId=191794314
將原本的非結(jié)構(gòu)網(wǎng)格轉(zhuǎn)化成了結(jié)構(gòu)網(wǎng)格,同時添加了關(guān)于LES相關(guān)設(shè)置的說明。
Turbulent Pipe Flow (LES) 圓管湍流流動(大渦模擬)
以ANSYS 17.0為例
問題描述
考慮通過圓形截面直管道的流動問題,圓管直徑D=0.0127mD=0.0127mD=0.0127m,長度L=5D=0.0635mL=5D=0.0635mL=5D=0.0635m。管道進口處的平均流速為Ubulk=6.58m/sU_{bulk}=6.58m/sUbulk?=6.58m/s,假設(shè)流體密度為定值,ρ=1.331kg/m3\rho=1.331kg/m^3ρ=1.331kg/m3,流體動力粘性系數(shù)μ=2.34×10?5kg/(m?s)\mu=2.34\times10^{-5}kg/(m\cdot s)μ=2.34×10?5kg/(m?s)。那么基于圓管直徑、平均流速、流體密度、動力粘性系數(shù)算得該問題的Reynold數(shù)(Re)為
Re=ρUbulkDμ=1.331?6.58?0.01272.34×10?5≈4750Re=\frac{\rho U_{bulk} D}{\mu}=\frac{1.331*6.58*0.0127}{2.34\times10^{-5}}\approx 4750Re=μρUbulk?D?=2.34×10?51.331?6.58?0.0127?≈4750
接下來咱們用ANSYS FLUENT中的LES方法來求解該流動問題,繪制在距離進口x/Dx/Dx/D處下游截面上隨著半徑變化的平均速度和均方根速度,并比較由LES方法和k??k-\epsilonk??方法模擬得到的平均速度。
1 預(yù)分析和準備工作
預(yù)分析
在大渦模擬中,瞬時速度U?(x?,t)\vec{U}(\vec{x},t)U(x,t)被分解為濾波后的分量U?ˉ(x?,t)\bar{\vec{U}}(\vec{x},t)Uˉ(x,t)以及剩余的殘差分量U?′(x?,t)\vec{U}'(\vec{x},t)U′(x,t),濾波后的速度分量表征了大尺度的非定常運動。在LES中,大尺度的湍流運動被直接表征,而小尺度的湍流運動則用模型近似。關(guān)于濾波速度的濾波方程可以從Navier-Stokes方程推出,由于殘差操作,動量方程中的非線性對流項引入了一個應(yīng)力張量的殘差項,該殘差應(yīng)力張量需要通過構(gòu)造模型來完成方程組的封閉,而FLUENT中提供了從易到難的多種模型。
既然咱們要求解U?ˉ(x?,t)\bar{\vec{U}}(\vec{x},t)Uˉ(x,t),那么LES就是個非定常的模擬過程,需要在時域內(nèi)向前推進。為了收集統(tǒng)計平均量,比如平均和均方根(root mean square(r.m.s.))速度,咱們需要首先達到統(tǒng)計上的穩(wěn)定狀態(tài)(然后再開展統(tǒng)計平均的處理)。作為對比,k??k-\epsilonk??模型求得的平均速度也一并給出。
關(guān)于LES的詳細理論和方程可以再很多湍流的書籍中找到。
準備工作
LES是三維非定常計算(只能適用于三維問題和非定常問題),那么計算域是全部的管道。
在打開ANSYS之前,先創(chuàng)建一個文件夾turbulent_pipe_LES,然后里面在創(chuàng)建一個ICEM文件夾和FLUENT文件夾,分別用來存放ICEM的建模和畫網(wǎng)格文件,以及FLUENT的計算文件。
2 構(gòu)建幾何模型
打開ICEM CFD 17.0軟件,在其中完成建模工作,咱們計算域是圓管內(nèi)部流道,也就是一個圓柱體,讓圓柱體的軸線沿著xxx方向,進口截面位于yOz平面yOz平面yOz平面上,圓心位于坐標原點OOO。
注意:在建模過程中所用長度單位為m。
創(chuàng)建點
Geometry -> Create Point -> Explicit Coordinates
依次創(chuàng)建坐標(x,y,z)(x,y,z)(x,y,z)為(0,0,0)(0,0,0)(0,0,0)、(0,0.00635,0)(0,0.00635,0)(0,0.00635,0)、(0,0,0.00635)(0,0,0.00635)(0,0,0.00635)、(0.0635,0,0)(0.0635,0,0)(0.0635,0,0)的點0、1、2、3。
創(chuàng)建線
Geometry -> Create/Modify Curve -> Cricle or arc from Center point and 2 points on plane. Optional Radius
依次選中前3個點,即圓心、圓環(huán)上兩個點,OK,生成進口圓環(huán)曲線crv.00
Geometry -> Create/Modify Curve -> From Points
選中第1和第4個點,即軸線上兩個點,OK,生成軸線crv.01
創(chuàng)建面
Geometry -> Create/Modify Surface -> Simple Surface
將Simple Surface Method中由默認的From 2-4 Curves改為From Curves
選中進口的圓環(huán),OK,生成進口面Suf.00
Geometry -> Create/Modify Surface -> Curve Driven
選擇Driving Cruve為crv.01,Driven Curves為crv.00,OK,生成環(huán)面即為圓管的壁面Suf.01。
Geometry -> Create/Modify Surface -> Simple Surface
將Simple Surface Method中由默認的From 2-4 Curves改為From Curves
選中出口的圓環(huán),OK,生成進口面Suf.02
創(chuàng)建體
Geometry -> Create Body -> Material Point
選中pnt.03和pnt.01兩個點,則創(chuàng)建好了圓柱體Body
創(chuàng)建邊界面
Parts 右擊 Create Part,Part中命名為inlet,選中進口面,OK;
同樣,將出口面命名為outlet;將側(cè)面命名為pipeWall。
刪除模型創(chuàng)建過程中多余的點、線
Geometry -> Delete Point,按下鍵盤中的A刪除所有的點;
Geometry -> Delete Curve,按下鍵盤中的A刪除所有的線;
Model -> Parts -> GEOM右擊,Delete,刪除沒有幾何元素的空Part。
現(xiàn)在就只剩下來了一個三維實體BODY,還有三個邊界面INLET、OUTLET和PIPEWALL。
創(chuàng)建幾何模型的拓撲結(jié)構(gòu)
Geometry -> Repair Geometry,保持默認,Apply,再次生成表征幾何體所必須的Point和Curve。
系統(tǒng)會根據(jù)體和面生成所必須的點和線。
之所以這么做,是為了保證體、面與點、線的關(guān)聯(lián)性,便于后續(xù)網(wǎng)格剖分。
3 剖分網(wǎng)格
還是在ICEM軟件中操作,接下來剖分幾何體為分區(qū)結(jié)構(gòu)化網(wǎng)格。
為確保完全,先保存下,Save Project到咱們最早新建的文件夾…\turbulent_pipe_LES\ICEM中,名字叫Pipe就好了。
創(chuàng)建整體Block
Blocking -> Create Block,
Part命名為Fluid,保持默認的Initialize Blocks,選擇默認的Type為3D Bounding Box,OK。
創(chuàng)建O-Block
Blocking -> Split Block,
選擇Ogrid Block,
Select Block(s),選擇需要劃分的Block;
Select Face(s),選擇進口與出口的兩個面;
OK,即可看到切分好的O-型分區(qū)。
建立映射關(guān)系
Blocking -> Associate -> Snap Project Vertices,選擇All Visible,OK。
則Block中的點自動和幾何體中最近的點映射到一起了。
Blocking -> Associate -> Associate Edge to Curves -> Associate Edge to Curves,
將對應(yīng)的Block中的Edge映射到幾何體的Curves上。
定義網(wǎng)格節(jié)點數(shù)目
Blocking -> Pre-Mesh Params -> Edge Params
將所有軸向Edge設(shè)置為均分160個單元,這樣單元尺度為h=0.0635m/160≈4×10?4mh=0.0635m/160\approx4\times10^{-4}mh=0.0635m/160≈4×10?4m;
將進出口處與1/4圓弧對應(yīng)的Edge均分為25個單元,這樣單元尺度為h=0.0127m?π/4/25≈4×10?4mh=0.0127m*\pi/4/25\approx4\times10^{-4}mh=0.0127m?π/4/25≈4×10?4m;
將壁面邊界層網(wǎng)格,即O-結(jié)構(gòu)網(wǎng)格的四個徑向Edge設(shè)置為首層網(wǎng)格h0=5×10?6mh_0=5\times10^{-6}mh0?=5×10?6m,等比增長比例為1.1,總共40層網(wǎng)格。
生成網(wǎng)格
Model -> Blocking -> Pre-mesh,Yes。
Model -> Blocking -> Pre-mesh右擊Pre-Mesh Info,可看網(wǎng)格信息:
其節(jié)點總數(shù)約70萬,其六面體單元總數(shù)約70萬。
檢查網(wǎng)格質(zhì)量
Blocking -> Pre-Mesh Quality Histograms,
選擇Criterion為Determinant 2x2x2,Apply。
選擇Criterion為Angle,Apply。
網(wǎng)格質(zhì)量還好的,可以導(dǎo)出網(wǎng)格做計算了。
導(dǎo)出網(wǎng)格文件
Model -> Blocking -> Pre-Mesh,右擊,選擇Convert to Unstruct Mesh,
信息框出現(xiàn)“Current Coordinate system is global”,則表示轉(zhuǎn)化完成。
File -> Mesh -> Save Mesh As,保存當前網(wǎng)格文件為Pipe_Str.uns。
Output Mesh -> Select solver,
保持Output Solver為默認的ANSYS Fluent,Apply。
Output Mesh -> Write input,
Yes保存現(xiàn)有Project,在跳出的窗口中選擇剛才保存的網(wǎng)格文件Pipe_Str.uns,
在隨后彈出的窗口中選擇3D,可將Output file名字改為Pipe_Str,Done完成輸出。
關(guān)閉ICEM軟件。
4 物理問題設(shè)置(Setup設(shè)置)
接下來需要在FLUENT當中進行CFD設(shè)置和計算了。
打開FLUENT軟件:
咱們求解的是三維問題,當然設(shè)置成3D了;
為了提高精度,一般都是用的Double Precision雙精度模式,因為默認的單精度real計算誤差太大了,沒誰會拿它來跑計算的;
Precessing Options下,打開Parallel并行計算,然后把Processes中設(shè)置成你想用的CPU核數(shù),由于LES是非定常和三維的運算,其計算量非常大,如果用Serial串行計算,那么可能要耗費數(shù)月的時間才能算好,所以條件允許的話,盡可能用多核并行計算吧,不然你要等到黃花菜都涼了。。。
讀入/檢查/縮放網(wǎng)格
File -> Read -> Mesh,找到剛才咱們輸出的Pipe_Str.msh文件,把它讀進來。
Setting Up Domain -> Info -> Size,查看下網(wǎng)格單元、面、節(jié)點數(shù)目。
Setting Up Domain -> Display,顯示網(wǎng)格。
Setting Up Domain -> Check,檢查網(wǎng)格,注意看又沒有負體積,計算域的x、y、z上下限是不是跟之前建模的一樣,如果不一致,是不是差了一定倍數(shù),那就需要做Scale縮放處理。這里一切都OK的,表明網(wǎng)格沒問題。
通有設(shè)置
Setup -> General,保持默認Pressure-Based Solver,壓力基求解器適用于不可壓縮流動,將Time下的默認Steady改為Transient,因為LES只能是非定常三維計算。
模型設(shè)置
Setup -> Models -> Viscous - Laminar,Edit…,Model選擇Large Eddy Simulation (LES),Subgrid-Scale Model選擇WMLES,這個模型對網(wǎng)格的要求最低(可查看FLUENT幫助文檔里的詳細說明),OK。
跳出來個信息框,其實是告訴你,要是用LES的話,空間離散和時間離散都要求比較高,必須要用到它推薦的這些離散格式。即,空間用有界中心差分,時間用有界二階隱式格式。當然,系統(tǒng)已經(jīng)自動幫你設(shè)置好了。OK關(guān)掉它就好了。
物性參數(shù)
Setup > Materials > Fluid > Create/Edit… ,
在打開的窗口中,設(shè)置密度ρ=1.331kg/m3\rho=1.331kg/m^3ρ=1.331kg/m3,流體動力粘性系數(shù)μ=2.34×10?5kg/(m?s)\mu=2.34\times10^{-5}kg/(m\cdot s)μ=2.34×10?5kg/(m?s),Change/Create,Close。
邊界條件
Setup > Boundary Conditions > inlet,確認其類型是velocity-inlet,Edit…,
設(shè)置Velocity Specification Method為Magnitude, Normal to Boundary:
設(shè)置Velocity Magnitude (m/s)為6.58 m/s,
設(shè)置Fluctuation Velocity Algorithm為Spectral Synthesizer (SS方法,在進口處生成擾動的湍流速度分布,而非均勻恒定的速度場);
設(shè)置Turbulence Specification Method為Intensity and Hydraulic Diameter,用強度和水力直徑設(shè)置湍流:
設(shè)置Turbulent Intensity (%)為10 %,湍流強度10%,
設(shè)置Hydraulic Diameter (m)為0.0127 m,水力直徑用圓管直徑給定即可;
設(shè)置Reynolds-Stress Specification Method為K or Turbulent Intensity。
OK關(guān)閉進口速度設(shè)置對話框。
確保outlet邊界類型為Pressure-outlet,并保持默認設(shè)置即可。
確保pipewall邊界類型為wall,并保持默認設(shè)置即可。
確保int_fluid邊界類型為interior,這個實際上是計算域內(nèi)部的面。
設(shè)置參考值
Setup > Reference Values,選擇Compute from為inlet,即可。
至此,物理問題描述完畢,可以選擇保存下文件到先前的文件夾…\turbulent_pipe_LES\FLUENT:
File -> Write -> Case,將Case命名為Pipe即可,OK。
5 求解器設(shè)置(Solution設(shè)置和求解)
求解方法
Solution -> Solution Methods,
設(shè)置Scheme為默認的SIMPLE,
設(shè)置Momentum為默認的Bounded Central Differencing,(還記得前面選擇LES模型時的提示么)
設(shè)置Pressure為Second order,
設(shè)置Transient為Second Order Implicit。
初始流場設(shè)置
Solution -> Solution Initialization, 選擇Compute from為inlet,然后點擊initialize來完成初始化流場。
設(shè)置收斂標準
Solution -> Monitors > Residuals - Print, Plot,Edit…,
將連續(xù)方程和x、y、z速度的收斂指標全部設(shè)置為1e-6,OK。
開始計算
先算到統(tǒng)計穩(wěn)態(tài)狀態(tài)
Solution -> Run Calculation,
設(shè)置Time Step Size(s)為1e-05,
設(shè)置Number of Time Steps為10000,
勾選Extrapolate Variables選項,
設(shè)置Max Iterations/Time Step為默認的20次。
為確保萬全,先保存下case和dat,然后點Calculate開始計算,愉快地出現(xiàn)了殘差收斂曲線,雖然剛開始20步不足以收斂到1e-6,但沒過幾部之后,系統(tǒng)就不需要20步就能收斂到1e-6了,所以計算精度還是OK的。
注意:這里的10000步算到0.1s,這樣的計算值僅僅是為了讓流場發(fā)展到統(tǒng)計穩(wěn)定狀態(tài)(也就是類似于振動力學(xué)里面的穩(wěn)態(tài)響應(yīng)狀態(tài)),然后再從0.1s開始后續(xù)的10000步計算到0.2s,邊算邊統(tǒng)計,以獲取0.1-0.2s的統(tǒng)計平均值。
由于我電腦比較爛,所以用的是單核串行計算……只能苦苦等待了……
經(jīng)過多天等待,終于算完了前10000步的0-0.1s,那么接下來就要再次計算10000步從0.1-0.2s,并且再算的過程中開啟統(tǒng)計功能。
再次計算前,可以看下中面上美麗的瞬態(tài)速度場、渦量場、壓力場的云圖,這些細碎凌亂、看似毫無規(guī)律、細思又隱含若干規(guī)則的結(jié)構(gòu),恰是湍流的瞬態(tài)特性。
速度場
渦量場
壓力場
緊接著再次計算并開啟統(tǒng)計量的計算
Solution -> Run Calculation,
設(shè)置和前面是一樣的,即:設(shè)置Time Step Size(s)為1e-05,設(shè)置Number of Time Steps為10000,勾選Extrapolate Variables選項,設(shè)置Max Iterations/Time Step為默認的20次。
唯一的不同是,此時勾選Data Sampling for Time Statistics選項,來開啟對時均統(tǒng)計量的數(shù)據(jù)采樣處理工作,這樣FLUENT在時域計算的同時便開始了對流場變量的時間平均計算(即統(tǒng)計平均量的計算),當時域計算完成后,除了能看到0.2s瞬時的物理量,還得到了0.1-0.2s的時均統(tǒng)計平均物理量。
點擊Calculate開始時域計算和統(tǒng)計量的計算,又是漫長的等待過程……
終于終于算好了,噢耶,可以后處理和結(jié)果分析了。
6 數(shù)值結(jié)果與分析
注意在LES模擬中有兩種速度,瞬時速度和時均速度。瞬時速度是計算域在某個時刻的真實速度,當我們收集統(tǒng)計量的時候,將瞬時速度做時間平均來獲取平均速度。讓咱們在計算域中設(shè)置一個中面來看看在它上面的瞬時軸向速度和平均軸向速度的云圖是個什么樣子吧。
軸向速度云圖
Results -> Graphics and Animation -> Contours -> Set Up…,
在Countours窗口,點擊New Surface -> Plane…,
在Plane Surface窗口,勾選Point and Normal,輸入點坐標為(x0 (m), y0 (m), z0 (m)) = (0,0,0), 輸入法向量為(ix (m), iy (m), iz (m)) = (0,0,1),將New Surface Name命名為. Name the surface as plane-mid,點擊Create完成中面創(chuàng)建。
可以用Graphics and Animations -> Mesh -> Set Up…,選擇和查看這個創(chuàng)建好的面。
回到contour窗口,在Contours of中選擇Velocity… 和X Velocity,在Surfaces中選擇plane-mid,點擊Display,下圖顯示出了瞬態(tài)軸向速度云圖。
若想繪制平均軸向速度,則在Contours of中選擇Unsteady Statistics…和Mean X Velocity,在Surfaces中選擇plane-mid,點擊Display,下圖顯示出了平均軸向速度云圖。
從這倆圖中,可以很清晰地看出瞬態(tài)量和平均量的不同之處,平均量恰如是做了定常假設(shè)的湍流場,而瞬態(tài)量則是瞬時湍流結(jié)構(gòu)的真實體現(xiàn),很多細碎漂亮的小渦結(jié)構(gòu),即傳說中的相干結(jié)構(gòu)(擬序結(jié)構(gòu))。
軸向速度的XY圖
咱們接下來整個橫穿管道中心的線,來看看這根線上平均軸向速度的分布情況,并在下個小節(jié)中將結(jié)果與k??k-\epsilonk??模型的結(jié)果作對比。
點擊Plots -> XY Plot -> Set Up…,
在Solution XY Plot窗口,點擊New Surface -> Line/Rake…,
在the Line/Rake Surface窗口,輸入端部點坐標(x0 (m), y0 (m), z0 (m)) = (0.03175,-0.00635,0)和(x1 (m), y1 (m), z1 (m)) = (0.03175,0.00635,0),將New Surface Name命名為line-mid,點擊Create。完成中線的創(chuàng)建。
這根線可以在Graphics and Animations > Mesh > Set Up…中查看。
回到Solution XY Plot窗口,在Y Axis Function中選擇Unsteady Statistics… 和Mean X Velocity,在Plot Direction輸入(X,Y,Z) = (0,1,0),在Surfaces下選擇line-mid,點擊Plot繪制平均軸向速度沿著中線的分布圖。
若要保存這些離散點,那么在Solution XY Plot窗口中選擇Write to File,點擊Write,選擇位置并將其保存為MeanXVel_LES.xy文件即可。
當然,也可以查看瞬態(tài)x速度的分布情況,與平均值相比,它凌亂多了。
7 驗證和確認(Verification & Validation)
用k??k-\epsilonk??模型來計算湍流管道流動
當然要把原有的cas和dat再復(fù)制一個復(fù)件出來,然后在上面更改設(shè)置了:
Setup -> General,選擇Steady定常問題,點擊OK關(guān)掉跳出來的警告。它提示LES必須用非定航算,不能用定常,咱們后面把LES改成k??k-\epsilonk??模型的。
Setup -> Models -> Viscous,Edit…,
在Viscous Model窗口,Model中選擇k-epsilon (2 eqn), k-epsilon Model中選擇Realizable,Near-Wall Treatment選擇Enhanced Wall Treatment,點擊OK確認。
因為咱們這個網(wǎng)格畫的足夠好(y+<1),所以可以用加強壁面處理方式來做近壁處理。
Setup -> Boundary Conditions,選擇Outlet點擊Edit…,將Turbulence Specification Method改為 Intensity and Hydraulic Diameter,將Backflow Turbulent Intensity設(shè)為10%,將Backflow Hydraulic Diameter設(shè)為0.0127m,OK確認。
接下來,Solution -> Solution Methods,對于Turbulent Kinetic Energy和Turbulent Dissipation Rate選擇Second Order Upwind。
Solution -> Monitors > Residual,點擊Edit…,對k和epsilon設(shè)置為1e-6。為了獲得較快的收斂,將continuity設(shè)為1e-5,OK確認。
Solution -> Solution Initialization,點擊Compute from選擇inlet,點擊initialize完成初始化。
Solutioin -> Run Calculation,將Number of Iterations設(shè)為1000,點擊Calculate,開始計算。
大約600步便收斂了。
中線速度比較
算好后,繪制中線上的x-velocity,并load進之前LES算得的平均x-velocity,將兩者繪制于同一個XY Plot中,白色為k-e,紅色為LES的統(tǒng)計平均量,可見,兩者吻合較好,從而驗證了LES計算結(jié)果的正確性。
總結(jié)
以上是生活随笔為你收集整理的FLUENT算例 —— Turbulent Pipe Flow (LES) 圆管湍流流动(大涡模拟)的全部內(nèi)容,希望文章能夠幫你解決所遇到的問題。
- 上一篇: WRF系列教程1:WRF如何得到更好的模
- 下一篇: MATLAB 格拉布斯准则代码