附錄F-有限體積法計(jì)算方腔流(F).doc

《附錄F-有限體積法計(jì)算方腔流(F).doc》由會(huì)員上傳分享，免費(fèi)在線閱讀，更多相關(guān)內(nèi)容在教育資源-天天文庫(kù)。

1、附錄F二維不可壓縮黏性流體方腔流動(dòng)問(wèn)題的有限體積算法與計(jì)算程序二維方腔流動(dòng)問(wèn)題是一個(gè)不可壓縮黏性流動(dòng)中典型流動(dòng)。雖然目前尚不能求得它的解析解，但是它常被用來(lái)作為檢驗(yàn)各種數(shù)值算法計(jì)算精度和可靠性的算例。文獻(xiàn)中幾乎大多數(shù)算法都對(duì)它進(jìn)行過(guò)計(jì)算。在本算例中采用有限體積算法三階迎風(fēng)型離散格式對(duì)它進(jìn)行數(shù)值求解。同時(shí)，為了初學(xué)者入門和練習(xí)方便,這里給出了用語(yǔ)言和語(yǔ)言編寫的計(jì)算二維不可壓縮黏性方腔流動(dòng)問(wèn)題計(jì)算程序，供大家學(xué)習(xí)參考。F-1利用有限體積算法三階迎風(fēng)型離散格式求解二維不可壓縮黏性流體方腔流動(dòng)問(wèn)題1.二維不可壓縮黏性流體方腔流動(dòng)問(wèn)題二維不可壓縮黏性流體方腔流動(dòng)（cavit

2、yflow圖F.1二維不可壓縮黏性方腔流問(wèn)題示意圖）：有一正方形腔室，其量綱為一的寬度為，里面充滿靜止的不可壓縮黏性流體，方腔內(nèi)初始時(shí)刻壓力和密度為它周圍壁面（左右壁面和底面）固定不動(dòng)，上壁面以量綱為一的速度沿著上壁面方向自左向右運(yùn)動(dòng)（圖F.1）。2.基本方程組、初始條件和邊界條件圖F.1二維不可壓縮黏性方腔流動(dòng)問(wèn)題示意圖設(shè)流體是黏性流體。二維方腔流動(dòng)問(wèn)題在數(shù)學(xué)上可以由二維不可壓縮黏性流動(dòng)N-S方程組來(lái)表示，把它寫成通用變量的微分方程組形式，有：（F.1）其中為變量在水平方向的流速，為在垂直方向的流速，為黏度，為源項(xiàng)。源項(xiàng)中不僅包含壓力梯度項(xiàng)，也包含時(shí)間導(dǎo)數(shù)項(xiàng)。初

3、始條件：方腔上壁面以量綱為一的速度沿著上壁面方向自左向右運(yùn)動(dòng)。邊界條件：流動(dòng)速度均可采用無(wú)滑移邊界條件，壓力采用自由輸出邊界條件。3.計(jì)算網(wǎng)格劃分和控制體單元與節(jié)點(diǎn)定義采用交錯(cuò)網(wǎng)格，圖F.2和圖F.3是計(jì)算網(wǎng)格、控制體單元和節(jié)點(diǎn)示意圖。圖F.2方腔流動(dòng)計(jì)算網(wǎng)格、控制體單元和節(jié)點(diǎn)示意圖圖F.3計(jì)算采用的交錯(cuò)網(wǎng)格示意圖節(jié)點(diǎn)所在主控制單元如圖F.2中有陰影部分所示。在方向與節(jié)點(diǎn)相鄰的節(jié)點(diǎn)為和，在方向與節(jié)點(diǎn)相鄰的節(jié)點(diǎn)為和，主控制單元界面分別為。壓力和速度分別在三套不同網(wǎng)格中如圖F.3中有陰影部分所示。4.有限體積算法三階迎風(fēng)型離散格式對(duì)方程（F.1）在圖F.2所示節(jié)點(diǎn)所在

4、控制體單元內(nèi)積分，有：（F.2）由于二維不可壓縮黏性流體方腔流動(dòng)是二維問(wèn)題，因此控制體單元體積僅是面積，而它的邊界是長(zhǎng)度。設(shè)，利用定理，可將方程（F.2）改寫成如下有限體積算法離散格式：（F.3）對(duì)上式中采用一階向前差分近似，則有：（F.4）同時(shí)記：（F.5）（F.6）則可由式（F.2）寫成：（F.7）式中都是控制體單元內(nèi)節(jié)點(diǎn)上的已知量，如果利用差分計(jì)算得到控制體單元邊界上的流通量，就可以求出節(jié)點(diǎn)上未知量。圖F.4三階迎風(fēng)型離散格式示意圖為了便于討論，現(xiàn)對(duì)一維對(duì)流擴(kuò)散方程的三階迎風(fēng)型離散格式進(jìn)行分析：在三階迎風(fēng)型離散格式中，計(jì)算主控制單元界面上流動(dòng)量需要取主控制單

5、元界面兩側(cè)3個(gè)節(jié)點(diǎn)處的流動(dòng)量值進(jìn)行插值計(jì)算得到，其中兩個(gè)節(jié)點(diǎn)位于界面緊鄰的兩側(cè)，第三個(gè)節(jié)點(diǎn)位于迎風(fēng)一側(cè)較遠(yuǎn)鄰點(diǎn)，如圖F.4所示。當(dāng)時(shí)，通過(guò)、和三個(gè)節(jié)點(diǎn)值擬合曲線來(lái)計(jì)算主控制單元左側(cè)界面參數(shù)。通過(guò)節(jié)點(diǎn)、和三個(gè)節(jié)點(diǎn)值擬合曲線來(lái)計(jì)算主控制單元右側(cè)界面參數(shù)。當(dāng)，則分別通過(guò)節(jié)點(diǎn)、、和、、三個(gè)節(jié)點(diǎn)值計(jì)算主控制單元左、右兩側(cè)界面參數(shù)和。根據(jù)上述計(jì)算原則，可以得到界面參數(shù)計(jì)算公式如下：當(dāng)時(shí)，界面參數(shù)計(jì)算公式為：（F.8a）當(dāng)時(shí)，界面參數(shù)計(jì)算公式為：（F.8b）對(duì)于一維無(wú)源項(xiàng)一維對(duì)流擴(kuò)散方程三階迎風(fēng)型離散格式：當(dāng)時(shí)，三階迎風(fēng)型離散格式為：（F.9）其中（F.9）同理，若，三階迎風(fēng)型

6、離散格式為：（F.10）其中（F.10a）將兩種流動(dòng)方向離散方程（F.9）和（F.10）合并后，可得到統(tǒng)一的一維對(duì)流擴(kuò)散方程三階迎風(fēng)型離散格式：（F.11）其中（F.11a）式中（F.11b）同理，可以得到帶有源項(xiàng)的二維對(duì)流擴(kuò)散方程三階迎風(fēng)型離散格式為：（F.12）其中為有限體積算法中源項(xiàng)平均值。式中各個(gè)系數(shù)為：（F.12a）式中（F.12b）源項(xiàng)為：（F.13）若把表示時(shí)刻動(dòng)量，表示時(shí)刻動(dòng)量，則可以得到源項(xiàng)離散格式為：（F.14）最后，得到有限體積算法二維對(duì)流擴(kuò)散方程三階迎風(fēng)型離散格式：（F.15）式中系數(shù)為一階迎風(fēng)格式中各對(duì)應(yīng)系數(shù)。5.計(jì)算結(jié)果分析利用三階迎風(fēng)

7、型離散格式和相應(yīng)的初始條件和邊界條件，求解二維不可壓縮黏性流體方腔流動(dòng)問(wèn)題。圖F.5是不同雷諾數(shù)條件下采用三階迎風(fēng)型離散格式得到的二維不可壓縮黏性流體方腔流動(dòng)的計(jì)算結(jié)果。計(jì)算結(jié)果和文獻(xiàn)中其他高精度算法得到的計(jì)算結(jié)果進(jìn)行了比較，兩者計(jì)算結(jié)果十分吻合，能把方腔下壁面兩個(gè)底角附近二次小渦清晰地計(jì)算出來(lái)。這表明有限體積算法三階迎風(fēng)型離散格式具有相當(dāng)高的計(jì)算精度。=1000=100=10000=5000圖F.5不同雷諾數(shù)條件下采用三階迎風(fēng)型離散格式計(jì)算二維不可壓縮黏性方腔流動(dòng)的計(jì)算結(jié)果從圖F.5中可以看出：二維不可壓縮黏性流體方腔流動(dòng)的中心大渦并不在中心位置，方腔內(nèi)流動(dòng)也并

8、不對(duì)稱。這