目錄
第8 章 有限體積法.2
§8.1 有限體積法的基本概念.2
8.1.1 幾個(gè)術(shù)語(yǔ)2
8.1.2 控制體積的選擇.3
8.1.3 結(jié)構(gòu)與非結(jié)構(gòu)網(wǎng)格.5
§8.2 有限體積法構(gòu)造7
8.9-1 方程離散及基本格式7
8.9-2 物理特性要求.11
8.9-3 迎風(fēng)型通量格式14
8.9-3 TVD 格式18
§8.3 非結(jié)構(gòu)網(wǎng)格上的有限體積法23
8.3.1 基本方程.23
8.3.2 離散基本思路24
8.3.3 數(shù)值通量近似.25
第9 章 在水力學(xué)問(wèn)題中的應(yīng)用.29
§9.1 滲流問(wèn)題中的應(yīng)用.29
9.1.1 飽和—非飽和地下水運(yùn)動(dòng)基本控制方程29
9.1.2 方程的離散31
9.1.3 算例【陳揚(yáng) 碩士論文】33
§9.2 二維明渠非恒定流計(jì)算.38
9.2.1 水流基本方程38
9.2.2 控制方程的離散39
9.2.3 計(jì)算實(shí)例.53
§9.3 三維紊動(dòng)分層流計(jì)算.64
9.3.1 紊動(dòng)分層流基本方程64
9.3.2 紊流模型及控制方程離散.65
9.3.3 壓力校正法66
9.3.3 邊界條件69
9.3.4 鹽度引起的負(fù)浮力流動(dòng)71

第8章 有限體積法
有限差分方法是從描述各種物理現(xiàn)象的基本微分方程出發(fā)構(gòu)造離散方程的，前文已經(jīng)對(duì)其作了翔實(shí)、周密的論述。該部分將從基礎(chǔ)算法入手分析介紹在計(jì)算流體力學(xué)界廣為應(yīng)用的有限體積法。基于有限體積法的實(shí)用算法在計(jì)算流體力學(xué)、計(jì)算傳熱學(xué)等領(lǐng)域得到了飛速發(fā)展 [1-3]。在水力學(xué)諸多問(wèn)題，如水流物質(zhì)輸運(yùn)模擬，水工水力學(xué)模擬以及潰壩洪水波演進(jìn)等水流模擬中也得到了廣泛應(yīng)用。
§8.1 有限體積法的基本概念
有限體積法，又稱為有限容積法，它正是從物理量守恒這一基本要求出發(fā)提出的。這也是其受計(jì)算流體力學(xué)界廣為稱道和喜歡之處。其以守恒型的方程為出發(fā)點(diǎn)，通過(guò)對(duì)流體運(yùn)動(dòng)的有限子區(qū)域的積分離散來(lái)構(gòu)造離散方程。有限體積法有兩種導(dǎo)出方式，一是控制容積積分法，另一個(gè)是控制容積平衡法。不管采用哪種方式導(dǎo)出的離散化方程，都描述了有限各控制容積物理量的守恒性，所以有限體積法是守恒定律的一種最自然的表現(xiàn)形式。
該方法適用于任意類型的單元網(wǎng)格，便于應(yīng)用來(lái)模擬具有復(fù)雜邊界形狀區(qū)域的流體運(yùn)動(dòng)；只要單元邊上相鄰單元估計(jì)的通量是一致的，就能保證方法的守恒性；有限體積法各項(xiàng)近似都含有明確的物理意義；同時(shí)，它可以吸收有限元分片近似的思想以及有限差分方法的思想來(lái)發(fā)展高精度算法。由于物理概念清晰，容易編程；有限體積法成為了工程界最流行的數(shù)值計(jì)算手段。
8.1.1 幾個(gè)術(shù)語(yǔ)
在進(jìn)行數(shù)值計(jì)算時(shí)，要把計(jì)算區(qū)域劃分成一系列互不重疊的離散小區(qū)域，然后在該小區(qū)域上離散控制方程求解待求物理量。在有限差分法中只涉及到網(wǎng)格節(jié)點(diǎn)的概念，而有限體積法因?yàn)槲锢斫忉屝枰纬闪艘韵聨讉€(gè)常用幾何要素的相關(guān)名詞。
控制體積（control volume）：圖中陰影部分，方程積分離散時(shí)的小體積單元（二維為面積單元）。
單元（cell）：控制體積的中心，常用形心來(lái)表示，為待求物理量的幾何位置。圖中用空心園來(lái)表示，如點(diǎn)P、W、E、N、S 等。常用單元來(lái)代表整個(gè)控制體積。以下如若不作特殊說(shuō)明，則用“單元”來(lái)代表控制體積。
網(wǎng)格線（grid lines）：用來(lái)分割計(jì)算區(qū)域內(nèi)各控制體積的交錯(cuò)曲線簇，如圖中折線N1N2、N2N3、N3N4、N4N1等
網(wǎng)格節(jié)點(diǎn)（nodes）：網(wǎng)格線之間的交點(diǎn)， 圖中用黑圓點(diǎn)來(lái)表示，如N1、N2、N3、N4等單元界面（cell faces）：相鄰兩個(gè)控制體積間的公共面(二維則可以認(rèn)為是公共邊，這樣看起來(lái)就和網(wǎng)格線一致了，但是要注意這不是同一個(gè)概念)，圖中用小寫字母e、n、w、s 表示。通常定義e、n、w、s 幾何位置位于交界面的形心點(diǎn)，二維則認(rèn)為在公共邊的中心點(diǎn)。
8.1.2 控制體積的選擇
當(dāng)你開(kāi)始用有限體積法模擬流體流動(dòng)時(shí)，而且劃分好網(wǎng)格后，你必須選定控制體積的形成方式。目前，常用的有兩種方法：?jiǎn)卧行姆绞剑╟ell-centered）和頂點(diǎn)中心方式
（vertex-centered）。另外一些學(xué)者還發(fā)展了由兩種方式綜合形成的混合方式。根據(jù)問(wèn)題的特點(diǎn)和要求，不同的變量可以采用不同的控制體積，因此又產(chǎn)生了交錯(cuò)網(wǎng)格和同位網(wǎng)格的稱謂，這里不再深入介紹，讀者根據(jù)需要可以參考相關(guān)文獻(xiàn)[1-3]。
前面我們我們?cè)谟懻撚邢摅w積法的術(shù)語(yǔ)時(shí)，已經(jīng)看到了單元中心控制體積的形成方式。
這里再作一個(gè)說(shuō)明并對(duì)兩者作簡(jiǎn)單的比較。圖8-2 為兩種控制體積選擇方式示意圖。陰影部分表示單元的控制體積。空心圓點(diǎn)表示單元，實(shí)心圓點(diǎn)表示網(wǎng)格節(jié)點(diǎn)，在頂點(diǎn)中心方式中，單元和網(wǎng)格節(jié)點(diǎn)重合。

（1）單元中心方式
單元位于控制體積的中心，即將單一的網(wǎng)格單元作為控制體積，網(wǎng)格單元互不重疊。目前用非結(jié)構(gòu)網(wǎng)格做流動(dòng)計(jì)算多使用單元中心格式,它有利于處理陸地和給定流量邊界條件。
（2）單元頂點(diǎn)方式
以網(wǎng)格節(jié)點(diǎn)為中心，它的形成有多種方式。其常用的構(gòu)成方式是由以該節(jié)點(diǎn)為頂點(diǎn)的格子的形心以及各共頂點(diǎn)的網(wǎng)格線中心點(diǎn)的一系列連接線段所構(gòu)成一個(gè)多邊形控制體積。也可以由環(huán)繞該節(jié)點(diǎn)的一組格子組成控制體。
因此計(jì)算同樣多的變量，單元中心方式變量布置簡(jiǎn)單直觀，易于處理邊界條件和保持離散的守恒性，而且需要的網(wǎng)格數(shù)要比單元頂點(diǎn)方式少得多，可節(jié)省計(jì)算時(shí)間。
8.1.3 結(jié)構(gòu)與非結(jié)構(gòu)網(wǎng)格
1 結(jié)構(gòu)網(wǎng)格
結(jié)構(gòu)網(wǎng)格，指的是一類具有一定的分布特征，可以用相應(yīng)的行列關(guān)系來(lái)順序描述的網(wǎng)格。
常用的矩形網(wǎng)格、曲線網(wǎng)格以及塊結(jié)構(gòu)網(wǎng)格等。
圖 8-3 結(jié)構(gòu)網(wǎng)格
矩形網(wǎng)格是FDM 或FVM 最為常用的網(wǎng)格系統(tǒng)。它不對(duì)原始的控制方程（直角坐標(biāo)）作任何坐標(biāo)變換，在邊界處采用簡(jiǎn)單的階梯形網(wǎng)格近似復(fù)雜的邊界。或者，用規(guī)則矩形網(wǎng)格覆蓋整個(gè)包括陸邊界在內(nèi)的計(jì)算區(qū)域，對(duì)陸邊界及陸域處的網(wǎng)格作特殊的技術(shù)處理，如所謂的凍結(jié)法，即令某些控制體積“凍結(jié)”，或者說(shuō)，“堵塞”某些控制體積，則剩下的“活動(dòng)”的控制體積可形成所需要的不規(guī)則區(qū)域[1]，“凍結(jié)”亦即令該控制體積中的因變量在運(yùn)算過(guò)程中始終保持不變。常采用大數(shù)值源項(xiàng)和大粘性系數(shù)法實(shí)現(xiàn)凍結(jié)區(qū)內(nèi)流速恒為零。盡管這類對(duì)規(guī)則網(wǎng)格作特殊的技術(shù)處理的方法對(duì)邊界概化過(guò)于嚴(yán)重，邊界計(jì)算值過(guò)于粗糙，陸邊界附近會(huì)形成虛假的曲折水流。但是，其網(wǎng)格生成方便，計(jì)算方法簡(jiǎn)單等優(yōu)點(diǎn)使該方法得到了大量的應(yīng)用。
在計(jì)算天然水域時(shí)，利用規(guī)則矩形網(wǎng)格階梯近似邊界不僅改變了原始物理邊界附近的流態(tài)，而且增加了邊界條件設(shè)置的復(fù)雜性；計(jì)算精度因此也大受影響，尤其在近邊界高梯度區(qū)域。自然地，研究人員引入廣義曲線坐標(biāo)變換的思想。即在物理計(jì)算區(qū)域內(nèi)構(gòu)筑曲線網(wǎng)格，使得網(wǎng)格曲線與所求解區(qū)域的邊界相重合，而后利用坐標(biāo)變換將復(fù)雜的物理域變換到規(guī)則的計(jì)算空間內(nèi)；在規(guī)則計(jì)算區(qū)域上離散求解變換后的控制方程，將得到的解投影至原物理區(qū)域得到數(shù)值解；或者在曲線網(wǎng)格上直接應(yīng)用原始因變量求解。除將復(fù)雜的物理空間變換成為易于求解的規(guī)則網(wǎng)格之外，曲線坐標(biāo)變換方法還能式計(jì)算網(wǎng)格點(diǎn)與物理域的活動(dòng)網(wǎng)格節(jié)點(diǎn)相對(duì)應(yīng)，從而可適應(yīng)非恒定流動(dòng)邊界的模擬。目前常用的生成曲線網(wǎng)格的方法有微分方法和代數(shù)方法。微分方法常用求解橢圓型方程或雙曲型方程實(shí)現(xiàn)。其中，尤以Thompson 為代表提出
的利用Poisson 方程生成貼體曲線網(wǎng)格的方法最為著名和流行。該類型方法可通過(guò)調(diào)整源項(xiàng)控制網(wǎng)格的分布，在梯度大的區(qū)域人工或自動(dòng)加密網(wǎng)格。因此，該方法提出后，得到了廣泛的發(fā)展和應(yīng)用。

曲線網(wǎng)格應(yīng)用遺留眾多需要解決的問(wèn)題。首先，大多數(shù)網(wǎng)格的生成方法只提供了離散點(diǎn)的變換，而不給出解析函數(shù)形式的變換關(guān)系。使用不光滑的網(wǎng)格時(shí)，對(duì)變換關(guān)系的差分近似會(huì)造成了很大的數(shù)值誤差，甚至?xí)?dǎo)致不切實(shí)際的值。其次，如果網(wǎng)格嚴(yán)重偏離正交性，就會(huì)極大損壞原有的迭代方法的收斂速率。再次，因變量的選擇也須謹(jǐn)慎考慮。在曲線網(wǎng)格中，可取原始笛卡爾坐標(biāo)系變量或曲線坐標(biāo)系中沿網(wǎng)格方向的協(xié)變量?jī)煞N作為因變量。
滿足正交性的網(wǎng)格具有優(yōu)良的特性。在正交曲線坐標(biāo)系下，變換后的某些項(xiàng)消失；控制方程更為簡(jiǎn)潔；數(shù)值穩(wěn)定性和解的精度高；計(jì)算量和收斂速度也有了很大的改善。
但是，正交曲線網(wǎng)格極大的受制于邊界的幾何形狀，對(duì)于天然水域，由于邊界極其復(fù)雜，
很難生成完全正交的曲線網(wǎng)格，特別對(duì)于三維問(wèn)題，正交曲線網(wǎng)格幾乎是不可能實(shí)現(xiàn)的。因此，在實(shí)際工程水流計(jì)算上較少采用。
正是由于復(fù)雜邊界下生成正交坐標(biāo)系的困難，而非正交曲線網(wǎng)格不僅相對(duì)較容易生成，而且可隨意調(diào)節(jié)網(wǎng)格的疏密，在高梯度區(qū)域人工或自動(dòng)加密網(wǎng)格。因此，研究人員退而求其次，放棄完全正交性的要求，轉(zhuǎn)而尋求非正交曲線網(wǎng)格生成和數(shù)值處理。如生成邊界正交曲
線網(wǎng)格，或生成盡量接近正交的曲線網(wǎng)格，以及探求非正交的光滑或不光滑或網(wǎng)格的求解方法[2]等。為充分利用網(wǎng)格正交性的優(yōu)點(diǎn)，常希望盡可能生成接近正交的網(wǎng)格，至少在接近
邊界高梯度區(qū)域，以恢復(fù)數(shù)值近似損失的精度。為此，Sorenson[11]開(kāi)發(fā)了一種能控制網(wǎng)格疏密和角點(diǎn)的邊界正交曲線網(wǎng)格生成的微分方法。
2 非結(jié)構(gòu)網(wǎng)格
八十年代以來(lái)，為了更好地?cái)M合自然邊界，非正交網(wǎng)格上的數(shù)值計(jì)算開(kāi)發(fā)和應(yīng)用得到了飛速的發(fā)展。
為了更好的說(shuō)明非結(jié)構(gòu)網(wǎng)格，可以先看圖8-4。
圖 8-4 非結(jié)構(gòu)網(wǎng)格
從圖中可以看出，非結(jié)構(gòu)網(wǎng)格中單元格分布不再規(guī)則一致，其位置很難再憑借行列索引關(guān)系確定。非結(jié)構(gòu)網(wǎng)格可以采用任意形狀的單元格，單元邊的數(shù)目也無(wú)限制，彌補(bǔ)了結(jié)構(gòu)化網(wǎng)格不能夠解決任意形狀和任意連通區(qū)域的網(wǎng)格剖分的缺欠。實(shí)用上，為簡(jiǎn)化編程以及貼合邊界要求，一般應(yīng)用四面體、六面體（二維為三角形、四邊形）等已經(jīng)足夠擬合天然水域邊界。非結(jié)構(gòu)網(wǎng)格最重要的一個(gè)特征是控制方程離散得到的代數(shù)方程的系數(shù)矩陣不再是結(jié)構(gòu)網(wǎng)格下有規(guī)律的對(duì)角結(jié)構(gòu)；若用對(duì)角形式存儲(chǔ)，其帶寬只能通過(guò)適當(dāng)?shù)牟贾脝卧幪?hào)順序來(lái)減少。非結(jié)構(gòu)網(wǎng)格原則上可應(yīng)用于任何類型的數(shù)值方法，包括FDM，但FEM 和非結(jié)構(gòu)網(wǎng)格上的FVM 算法更成熟，應(yīng)用更廣。

非結(jié)構(gòu)網(wǎng)格最早用于FEM，但流體流動(dòng)是高度非線性問(wèn)題，而且FEM 計(jì)算量較大，這些問(wèn)題使得基于FEM 的非結(jié)構(gòu)網(wǎng)格技術(shù)未能在對(duì)流問(wèn)題為主的地面水流（如淺水流動(dòng)，水波運(yùn)動(dòng)等）計(jì)算上得到重視。八十年代以來(lái)，基于FVM 的非結(jié)構(gòu)網(wǎng)格技術(shù)在空氣動(dòng)力學(xué)得到了廣泛的發(fā)展和應(yīng)用。九十年代開(kāi)始一些專家學(xué)者根據(jù)淺水流動(dòng)特征，將這些算法引入到計(jì)算淺水動(dòng)力學(xué)中，并在模擬涌潮，潰壩等水力計(jì)算難題上取得了成功[5]。
與此同時(shí)，粘性流動(dòng)的非結(jié)構(gòu)網(wǎng)格FVM 模擬也開(kāi)始出現(xiàn)。并在20 世紀(jì)90 年代中后期掀起了研究高潮[6]。作為全球計(jì)算流體力學(xué)軟件供應(yīng)商和技術(shù)服務(wù)商的Fluent 公司已經(jīng)將最新的非結(jié)構(gòu)網(wǎng)格研究成果集成，實(shí)現(xiàn)了研究成果的商業(yè)化。
非結(jié)構(gòu)網(wǎng)格能很好地模擬自然邊界及水下地形，利于邊界調(diào)節(jié)的實(shí)現(xiàn)；便于控制網(wǎng)格密
度，易作修改和適應(yīng)性調(diào)整；網(wǎng)格生成有眾多富有成效的方法和自適應(yīng)技術(shù)[9]，比曲線網(wǎng)格更易得到高質(zhì)量的單元格。但是非結(jié)構(gòu)網(wǎng)格應(yīng)用也帶來(lái)一系列需要解決的問(wèn)題：?jiǎn)卧衽帕胁灰?guī)則，須建立相應(yīng)的數(shù)據(jù)結(jié)構(gòu)存儲(chǔ)單元格信息；控制方程離散得到的代數(shù)方程的系數(shù)矩陣是高度稀疏的非對(duì)角型矩陣，需要尋求合適的存儲(chǔ)方式及解法；隱格式較難實(shí)現(xiàn)，粘性項(xiàng)處理困難；數(shù)值解后處理工作量大；二階非結(jié)構(gòu)FVM 較易實(shí)現(xiàn)，若要擴(kuò)展到高階格式，則
需花費(fèi)較大的代價(jià)。
§8.2 有限體積法構(gòu)造
8.9-1 方程離散及基本格式
1 基本思路
與有限差分法不同，有限體積法是基于守恒型的控制方程。這里以一維對(duì)流方程為例，
說(shuō)明有限體積法離散的基本思路。

為推導(dǎo)離散化的方程，我們必須事先對(duì)計(jì)算區(qū)域進(jìn)行網(wǎng)格劃分，然后在各個(gè)單元上進(jìn)行積分平均。假設(shè)我們對(duì)一維計(jì)算域劃分得到一個(gè)網(wǎng)格系統(tǒng)，如圖 所示的為該網(wǎng)格系統(tǒng)中某單元鄰接關(guān)系示意圖。重點(diǎn)考察單元P，單元E 和W 分別為它的兩個(gè)相鄰單元。虛線為單元界面，圖中用小寫字母表示，界面e、w 包圍形成的陰影部分即為單元P 的控制體積。控
制體積實(shí)際是三維體的的概念，只不過(guò)在考慮一維控制體積時(shí)，實(shí)際上把y、z 兩個(gè)坐標(biāo)方向假設(shè)成單位厚度，這樣在一維問(wèn)題中看起來(lái)成為了線的概念。

2 狀態(tài)變量分布近似
從前面可以知道，用有限體積法推導(dǎo)離散化方程時(shí)，必須確定物理量的局部分布，這是離散過(guò)程極為重要的一步，不僅關(guān)系到守恒性是否保持，而且關(guān)系到算法精度，這是由于不同的狀態(tài)變量分布對(duì)應(yīng)了不同的離散格式。

值得指出的是，在用有限體積法計(jì)算時(shí)，和有限差分法一樣，方程的解是用單元節(jié)點(diǎn)上離散點(diǎn)值構(gòu)成的，而不關(guān)心單元間的狀態(tài)變量是怎么變化的，也就是不關(guān)心解的分布。這和有限單元法中一旦選定了分布曲線，就確定了狀態(tài)變量的分布函數(shù)是不同的。盡管我們?cè)陔x散方程的時(shí)候要假定單元分布曲線，但是這只是為了推導(dǎo)公式時(shí)計(jì)算積分近似而采用的一些輔助關(guān)系式。一旦離散化方程推導(dǎo)出來(lái)了，就可以不用再管這些分布曲線近似。這種觀點(diǎn)使我們?cè)诓捎煤畏N分布曲線近似方法有完全的自由。在積分離散時(shí)，根據(jù)數(shù)值模擬的需要，對(duì)控制方程中的每一項(xiàng)都可以采用不同的分布曲線來(lái)近似單元界面上的狀態(tài)變量，而不必要追求近似假設(shè)的一致性。例如，我們用一階迎風(fēng)格式離散對(duì)流擴(kuò)散方程時(shí)，對(duì)流項(xiàng)中因變量是用階梯型分布近似，而擴(kuò)散項(xiàng)中也取為階梯分布的話，則根本導(dǎo)不出離散式，至少要有二階的分段線性分布才行。上面已經(jīng)提到，分布曲線近似關(guān)系到算法的精度，實(shí)際上，有限體積法中，每一種不同的狀態(tài)變量分布近似，體現(xiàn)了不同的離散格式的幾何構(gòu)造方法。如對(duì)流問(wèn)題中，一階迎風(fēng)格式為階梯型分布，中心格式采用的分段線性分布，更為高階的格式則需要更多的單元加入，形成復(fù)雜的狀態(tài)變量分布近似關(guān)系，例如對(duì)流項(xiàng)二次迎風(fēng)插值QUICK，分段拋物插值PPM 等。
3 基本格式
利用這些交界面表達(dá)式可以得到具體的離散化的代數(shù)方程，這里不在敘述。在這里，我們還可以看出，有限體積法離散的思路和有限差分法的有著本質(zhì)的區(qū)別，有限差分法是應(yīng)用Taylor 展開(kāi)式近似微分控制方程中的導(dǎo)數(shù)，而有限體積法則是通過(guò)積分將導(dǎo)數(shù)項(xiàng)化成交界面的狀態(tài)變量的表達(dá)式，然后，代入根據(jù)局部分布近似導(dǎo)出的交界面值得到完整的離散化的代數(shù)方程。
對(duì)于中心格式和一階迎風(fēng)格式的特點(diǎn)和有限差分法的一致。中心格式精度較高，具有二階精度，但是不能反映對(duì)流傳輸機(jī)制。對(duì)非線性問(wèn)題，常會(huì)出現(xiàn)非線性不穩(wěn)定，為了消除該影響，需要添加人工粘性。一階迎風(fēng)格式穩(wěn)定性好，但是數(shù)值耗散大。值得強(qiáng)調(diào)的是，一階迎風(fēng)格式離散思想中蘊(yùn)涵的物理本質(zhì)，即波動(dòng)傳波的信息主要來(lái)自上游方向，是我們構(gòu)造高分辨率格式的基礎(chǔ)。

8.9-2 物理特性要求
1 守恒性
如果對(duì)一個(gè)離散方程在定義域的任一有限空間內(nèi)作求和的運(yùn)算（相當(dāng)于連續(xù)問(wèn)題中對(duì)微分方程作積分），所得的表達(dá)式滿足該區(qū)域上物理量守恒的關(guān)系時(shí)，則稱該離散格式具有守恒特征。
我們知道，根據(jù)流體運(yùn)動(dòng)的不同特征，可以采用拉格朗日（Lagrange）法和歐拉（Euler）法來(lái)描述流體的運(yùn)動(dòng)。其中歐拉法著眼于研究固定空間位置上各流體質(zhì)點(diǎn)運(yùn)動(dòng)要素的空間變化情況，在推導(dǎo)流體運(yùn)動(dòng)的控制方程時(shí)，歐拉法被廣為采用。讓我們回顧一下連續(xù)方程的推導(dǎo)過(guò)程，首先在流場(chǎng)中劃定一個(gè)固定的空間區(qū)域作為流體運(yùn)動(dòng)的控制體積，考察流體流入、流出該控制體積的情況，利用質(zhì)量守恒原理，可導(dǎo)出流體運(yùn)動(dòng)的連續(xù)性方程。其他描述流體中物質(zhì)（能量）輸運(yùn)的控制方程也可根據(jù)各自的守恒律推導(dǎo)出來(lái)。在數(shù)值計(jì)算中，物理量的守恒性也必須得到保證，否則會(huì)導(dǎo)致非物理解的產(chǎn)生。因此守恒性成了數(shù)值求解過(guò)程中非常重要的一個(gè)方面。
有限體積法正是從物理量守恒這一基本要求出發(fā)提出的。它有兩種導(dǎo)出方式[3]，一是控制容積積分法，另一個(gè)是控制容積平衡法。不管采用哪種方式導(dǎo)出的離散化方程，都描述了有限各控制容積物理量的守恒性。就像微分方程是表示關(guān)于無(wú)窮小控制體積內(nèi)物理量的守恒一樣，帕坦卡[1]（S.V. Patankar）非常生動(dòng)的對(duì)此作了注解：如果把我們想象成處于微積分以前的時(shí)期，控制容積方程就會(huì)成為我們表達(dá)守恒原理的唯一途徑。有限體積法的離散化方程滿足了單個(gè)控制體積的平衡，當(dāng)然在整個(gè)計(jì)算區(qū)域內(nèi)，諸如質(zhì)量、動(dòng)量等物理量的積分守恒也就都能精確得到滿足。無(wú)論在數(shù)值計(jì)算中采用巨大數(shù)目的細(xì)網(wǎng)格和少數(shù)的粗網(wǎng)格，數(shù)值解也照樣顯示準(zhǔn)確的積分平衡。有限體積法的離散思想自動(dòng)滿足守恒定律，如質(zhì)量守恒，動(dòng)量守恒，能量守恒等等。所以有限體積法是守恒定律的一種最自然的表現(xiàn)形式。
為了得到守恒型的離散方程，需要滿足一定的條件。這里仍以前文導(dǎo)出的全離散的守恒型顯格式離散方程（8-12）為例。一般對(duì)l + k +1個(gè)單元（稱為節(jié)點(diǎn)模板），守恒型數(shù)值通量

守恒格式可以自動(dòng)算出激波并給出其正確的位置，因此這種數(shù)值方法稱為激波捕獲法（shock-capturing method）。如果使用非守恒格式，則需要在有關(guān)位置使用局部Rankine—Hugoniot 條件，以獲得準(zhǔn)確的激波位置，這種方法稱為激波擬合法(shock-fitting method)。
2 迎風(fēng)性
在有限差分法部分中的基本理論部分，我們已經(jīng)深入的討論了雙曲型方程（對(duì)流）的特征線概念。特征值體現(xiàn)了微分方程的解（或者說(shuō)擾動(dòng)、信息）的傳播方向。它既有很強(qiáng)的數(shù)學(xué)意義，更反映了問(wèn)題的物理機(jī)制，表達(dá)了波動(dòng)、能量等的傳播方向。這啟示我們?cè)诶碚摲?br /> 析和數(shù)值方法的設(shè)計(jì)上，應(yīng)該充分考慮這一物理特性，與之相適應(yīng)，而不可反其道而行之。根據(jù)特征（信息）傳播方向，設(shè)計(jì)得到的數(shù)值方法可統(tǒng)稱為特征型或迎風(fēng)型格式。例如前文
提到的傳統(tǒng)的一階迎風(fēng)格式，盡管非常簡(jiǎn)單，卻是蘊(yùn)含了最原始的的迎風(fēng)設(shè)計(jì)思想，即在計(jì)算單元界面處的數(shù)值通量時(shí)，應(yīng)根據(jù)特征傳播的方向性，對(duì)于向右傳播的分量，應(yīng)該使用左
邊的單元值來(lái)計(jì)算(因?yàn)橛绊憗?lái)自左邊)；對(duì)于向左傳播的分量，應(yīng)該使用右邊的值來(lái)計(jì)算(因?yàn)橛绊憗?lái)自右邊)。
當(dāng)前，這種根據(jù)信息傳播方向設(shè)計(jì)數(shù)值格式的迎風(fēng)思想已經(jīng)成為了流體數(shù)值計(jì)算發(fā)展的重要基石，它不斷地啟發(fā)研究人員構(gòu)造出各種各樣的離散格式，如特征迎風(fēng)格式，矢通量分裂格式，通量差分裂格式等。并在這些基于迎風(fēng)思想的格式基礎(chǔ)上，發(fā)展了分門別類的高階、高分辨率格式，如MUSCL 格式、QUICK 格式、ENO（WENO）格式等等。
3 TVD 性

既要利用格式的高階精度，又能使格式具有TVD 性質(zhì)，從而求得高精度的物理上有意義的解，是研究人員的目標(biāo)，為此，過(guò)去的二十幾年中，各種二階（甚至更高階）的TVD 格式得到了大量的發(fā)展。這類高階TVD 格式開(kāi)發(fā)應(yīng)用成為了計(jì)算流體力學(xué)界最前沿的研究。
8.9-3 迎風(fēng)型通量格式
對(duì)于雙曲型問(wèn)題，不同特征分量傳播的方向不同。考慮傳播方向的影響從物理上分析更為合理。我們對(duì)傳統(tǒng)的一階迎風(fēng)格式已經(jīng)作了介紹。這里我們將介紹20 世紀(jì)80 年代發(fā)展起來(lái)的通量向量分裂格式（Flux Vector Splitting），通量差分裂格式（Flux Difference Splitting）以及Godunov 法和近似黎曼解算子法（Approximate Riemann Solvers）。
1 矩陣分裂
實(shí)際流動(dòng)的控制方程是一個(gè)非線性系統(tǒng)，各個(gè)變量之間相互耦合。我們這里考慮一維的雙曲型方程組
2 通量向量分裂格式（FVS）
Steger 和Warming 根據(jù)矩陣分裂方法，于1981 年首次提出了一種通量向量分裂格式。
最早該格式應(yīng)用于空氣動(dòng)力學(xué)計(jì)算上，20 世紀(jì)90 年代，開(kāi)始擴(kuò)展到計(jì)算水力學(xué)領(lǐng)域。

在把通量向量格式應(yīng)用到水力學(xué)中時(shí)，仍需要注意方程組的特性。因?yàn)樗鬟\(yùn)動(dòng)的控制方程——圣維南（Saint Venant）方程組的通量項(xiàng)不具有齊次性，需要補(bǔ)充一個(gè)冗余的能量方程得到具有齊次通量的方程組，然后就可以按照前面的方法建立離散化方程。
3 通量差分裂格式（FDS）
通量向量分裂格式是將通量函數(shù)按不同特征方向進(jìn)行分裂，然后積分離散。與通量向量分裂格式不同，通量差分裂格式是對(duì)任意區(qū)間上通量函數(shù)的差進(jìn)行分裂的。
4 Godunov 格式和近似黎曼解法

8.9-3 TVD 格式
在前文討論中，我們已經(jīng)知道，滿足TVD 特性的格式求出的數(shù)值解才是真正符合實(shí)際的物理解。然而，一階迎風(fēng)格式盡管具有TVD 性，但是數(shù)值耗散過(guò)大，將本來(lái)應(yīng)該很陡峭的激波給抹得很平．為了在解光滑的區(qū)域獲得較高精度的解，并且在激波等間斷區(qū)域獲得沒(méi)有數(shù)值振蕩或者是數(shù)值振蕩可以令人接受的較陡峭的數(shù)值解，人們開(kāi)始研究高分辨格式（High Resolution schemes）。1983 年Harten 在他的論文中提出了高分辨率和總變差減少的概念，首先證明了計(jì)算格式具有TVD 性質(zhì)的充分條件，具體構(gòu)造了修正通量的高分辨率（二階）TVD 格式。他的成果在數(shù)值算法的研究上具有十分重要的意義，，為構(gòu)造和分析高分辨率格式提供了有力的工具。隨后，各種各樣的TVD 格式如雨后春筍般的涌現(xiàn)出來(lái)，因其具有對(duì)間斷的自動(dòng)高分辯識(shí)別能力，而得到了廣泛的應(yīng)用。這種格式具有兩個(gè)特點(diǎn)：（1）對(duì)標(biāo)量非線性方程及常系數(shù)雙曲型方程組，格式的解是總變差遞減的（Total Variation Diminishing）；（2）與守恒律和熵不等式是相容的。TVD 格式提高了對(duì)激波的捕捉能力，但其不足之處在極值點(diǎn)上格式只有一階逼近精度。

§8.3 非結(jié)構(gòu)網(wǎng)格上的有限體積法
前面主要對(duì)有限體積法基本概念和離散格式作了介紹。在這節(jié)中，我們將介紹二維非結(jié)構(gòu)網(wǎng)格上的有限體積法，以便于應(yīng)用它來(lái)模擬自然中復(fù)雜區(qū)域內(nèi)的流動(dòng)及物質(zhì)輸運(yùn)現(xiàn)象。本節(jié)只對(duì)算法的空間離散進(jìn)行討論，因?yàn)闀r(shí)間的離散和有限差分法一致，因此，不在專門介紹。

標(biāo)簽： 點(diǎn)擊： 評(píng)論: