中科院計算流體力學講義CFD117講有限體積法

有限體積法中科院計算流體力學講義·第117講主講人：李新亮研究員·中科院力學所Contents目錄中科院計算流體力學講義·有限體積法專題01有限體積法基本概念與核心原理02重構與反演：格式構造的技術核心03迎風型有限體積法與Riemann求解器04多維迎風型有限體積法與坐標旋轉CHAPTER01有限體積法基本概念與核心原理從守恒律的積分形式出發，建立離散求解的基本框架ComputationalFluidDynamics有限體積法的基本思想有限體積法以守恒律的積分形式為出發點，將計算域劃分為互不重疊的控制體積，通過對每個控制體積積分建立離散方程。其核心優勢在于天然保證守恒性，且對包含間斷的流場具有堅實的數學基礎。01控制體積劃分：將計算域劃分為一系列互不重疊的控制體積（cell），每個控制體積圍繞一個網格節點或位于單元中心02積分方程轉化：對每個控制體積上的守恒方程進行體積分，將微分方程轉化為積分方程，使間斷解在弱解意義下仍然成立03離散守恒律：控制體積內物理量的時間變化率等于通過其所有界面流入流出的凈通量，離散層面嚴格滿足守恒律04與有限差分的本質區別：有限差分基于微分形式的點態離散，有限體積基于積分形式的區域平均離散有限體積法控制體積劃分示意CFD·NumericalMethods有限體積法與有限差分法的本質對比有限體積法與有限差分法在離散基礎、守恒性保證、網格適應性和物理直觀性四個維度存在根本差異。有限差分法FiniteDifferenceMethodTaylor展開近似—基于微分形式在網格節點處展開，要求網格光滑且拓撲規則DifferentialForm守恒性需特殊設計—普通中心差分在間斷附近可能產生非物理解Non-Conservative高精度格式成熟—WENO、緊致差分等，適合基礎研究與規則幾何域WENO/Compact有限體積法FiniteVolumeMethod天然離散守恒—基于積分形式對控制體積內平均值演化，無需額外設計IntegralForm非結構網格適應—三角形、四面體等，對復雜幾何外形適應力優異UnstructuredMesh物理圖像直觀—界面通量對應真實輸運過程，便于理解和實現邊界條件PhysicalFluxFINITEVOLUMEMETHOD一維守恒律的有限體積離散推導一維守恒律方程的有限體積離散將微分方程轉化為控制體積平均值的時間演化方程。半離散格式的核心在于界面數值通量的計算，這正是重構與反演技術要解決的問題。01定義控制體積與單元平均值——在計算域上劃分控制體積[xj-?,xj+?]，引入ūj=(1/Δx)∫udx，對守恒律方程在整個控制體積上積分。02Newton-Leibniz轉化為邊界通量——dūj/dt=?(1/Δx)[fj+??fj-?]，體積積分化為界面通量差，即半離散有限體積格式。03核心信息缺口——已知的是單元平均值ūj，但格式需要的是界面處的數值通量fj+?，從平均值到界面值之間存在信息缺口。04時間與空間離散——時間推進可采用多種Runge-Kutta方法，空間精度完全取決于界面通量fj+?的重構精度。FINITEVOLUMEMETHOD控制體積的定義方式與選擇有限體積法的控制體積有節點中心型和單元中心型兩種基本定義方式。在一維均勻網格上兩者等價，但在多維非結構網格上各有優劣，選擇需綜合考慮精度、穩定性和計算效率。節點中心型Vertex-Centered控制體積圍繞網格節點構建，相鄰單元中點連線形成對偶網格作為實際控制體積存儲量較少（節點數通常少于單元數），適合大規模三維計算，但非結構網格上重構較復雜對偶網格單元中心型Cell-Centered未知量存儲在網格單元中心，控制體積即為網格單元本身，概念更直觀、實現更簡單在激波捕捉和間斷分辨率方面表現更穩健，高階重構（如MUSCL、WENO）更容易在單元中心框架下實現高階重構FINITEVOLUMEMETHOD離散守恒性：有限體積法的理論基石有限體積法在離散層面嚴格保持守恒性，相鄰控制體積共享界面通量使得內部通量精確抵消。這一性質保證了激波的Rankine-Hugoniot條件在離散層面自動滿足，是激波捕捉格式的數學基礎。01相鄰控制體積共享界面通量：fj+1/2對第j個單元是流出，對第j+1個單元是流入，符號相反、數值相同界面通量共享02全域求和時內部通量精確抵消（telescopingproperty），總物理量變化僅由邊界通量決定，離散守恒嚴格成立TelescopingProperty03Rankine-Hugoniot跳躍條件在離散層面自動滿足，保證激波傳播速度的正確性，無需人工干預激波位置R-H條件04守恒性是獲得正確弱解的必要條件：Lax-Wendroff定理表明，守恒型格式的收斂極限必為原方程的弱解Lax-Wendroff定理CHAPTER02重構與反演：格式構造的技術核心從單元平均值到界面數值通量的兩步關鍵技術FINITEVOLUMEMETHOD重構與反演：有限體積法的兩步核心技術有限體積法的空間離散分為重構和反演兩個串聯步驟。重構從單元平均值恢復界面處的物理量分布，決定空間精度；反演從界面兩側狀態計算數值通量，決定波傳播和間斷處理能力。重構Reconstruction空間分布構造從單元平均值ūj構造單元內物理量的空間分布uj(x)，可采用分段常數、分段線性或高階多項式逼近，是有限體積法實現高階精度的基礎環節精度決定階數k階多項式重構配合適當的數值通量可獲得k+1階空間精度，重構階數直接決定了格式的理論收斂階間斷處理高階重構需引入限制器或WENO權值，防止間斷附近產生虛假振蕩，保證解的TVD/TVB有界性反演FluxEvaluation數值通量計算利用界面左右狀態uL和uR計算數值通量fj+?，本質為求解局部Riemann問題，是有限體積法的核心計算步驟耗散與分辨率從簡單迎風格式到精確Riemann求解器，方法選擇決定格式的數值耗散特性與激波分辨率，影響計算結果的銳度相容性條件當uL=uR=u時，數值通量應退化為物理通量f(u)，確保格式在光滑區域具有一致性FINITEVOLUMEMETHOD重構階數：從分段常數到分段線性重構的階數直接決定有限體積格式的空間精度。分段常數重構僅一階精度但絕對穩定，分段線性重構（MUSCL）達到二階精度但需配合限制器防止振蕩。重構階數的選擇需要在精度和魯棒性之間權衡。分段常數重構一階精度u_j(x)=ū_j，界面左右值即為相鄰單元平均值。實現簡單、計算量小，耗散較大但絕對穩定，適用于對穩定性要求極高的場景。絕對穩定·無振蕩分段線性重構MUSCL·二階u_j(x)=ū_j+σ_j(x?x_j)，斜率由相鄰單元差分計算。界面值精度顯著提升，是工程計算中最常用的重構策略。需配合限制器使用斜率計算方式精度?穩定中心差分精度最高但可能引入振蕩，單側差分更穩定但精度降低。實際應用中需根據流場特征和網格質量進行權衡選擇。Minmod·Superbee·MC高階重構≥三階·WENO三階及以上需使用二次或更高次多項式，配合WENO等自適應權值技術。在光滑區實現高分辨率，間斷處自動降階保持穩定性。高分辨率·計算成本高LIMITERFUNCTIONS斜率限制器：高階重構的振蕩控制斜率限制器是保證高階重構在間斷附近不產生虛假振蕩的關鍵技術。通過在光滑區保留完整斜率、在間斷附近自動壓縮斜率，實現TVD性質。經典限制器函數MinMod限制器φ(r)=max(0,min(1,r))，取左右差分中絕對值較小者，最保守但保證無振蕩，激波附近耗散偏大VanLeer限制器φ(r)=(r+|r|)/(1+r)，用連續函數實現平滑過渡，兼顧穩定性和分辨率，工程中應用最廣Superbee限制器φ(r)=max(0,min(2r,1),min(r,2))，在接觸間斷附近允許更陡斜率，分辨率高但對強激波可能產生輕微振蕩理論基礎Godunov定理線性單調格式最高一階精度，高階無振蕩格式必須是非線性的——限制器正是引入非線性的機制TVD條件TotalVariationDiminishing要求限制器函數滿足特定約束區域，保證總變差不增大從而杜絕虛假振蕩FINITEVOLUMEMETHODWENO重構：自適應高階精度方法WENO重構通過多候選模板的自適應加權組合，在光滑區域實現最優階數精度，在間斷附近自動降階并保持無振蕩特性。相比傳統限制器，WENO在復雜波系結構中展現出更優異的分辨率和魯棒性。構造r個候選子模板，每個子模板生成一個k階重構多項式，整體最優模板可達2k?1階精度2k?1階光滑度指示子βk衡量各子模板上解的光滑程度，含間斷的模板權重被自動壓低至接近零βk最終重構值為各子模板多項式的非線性加權平均：光滑區恢復最優精度，間斷區自適應降階保證無振蕩自適應降階WENO-Z改進版通過全局光滑度指標τ重新定義權重，在臨界點處恢復設計精度，減少經典WENO的過度耗散τFINITEVOLUMEMETHOD數值通量計算：基礎反演方法數值通量的計算（反演）是有限體積法的第二步核心技術。從中心通量到Lax-Friedrichs再到Rusanov通量，通過引入與波速相關的數值耗散實現格式穩定。01中心通量f=[f(uL)+f(uR)]/2，無耗散但絕對不穩定，不能單獨用于雙曲型方程的數值求解。不穩定02Lax-Friedrichs通量加入全局最大波速α控制的耗散項?α(uR?uL)/2實現穩定，但耗散過大。全局α03Rusanov通量使用界面局部最大特征速度αj+?代替全局α，顯著減少光滑區域的數值耗散。局部α04共同特點不需要求解Riemann問題，實現簡單且通用，適合多物理場耦合中的初步計算。簡單通用CHAPTER03迎風型有限體積法與Riemann求解器基于波傳播物理圖像的數值通量構造方法ComputationalFluidDynamicsRiemann問題：迎風型方法的數學基礎Riemann問題是初始條件為分段常數雙曲型方程的初值問題。間斷分解產生的波系結構構成了迎風型數值通量的物理基礎。01Riemann問題定義：雙曲守恒律在初始時刻具有分段常數初值的初值問題，界面處02解的波系結構：初始間斷分解為若干基本波——激波、稀疏波、接觸間斷，波之間由恒定03Euler方程的Riemann解：通常包含三個波——左行波（激波或膨脹波）、接觸間斷、右行波，中間有兩個星區狀態u*L和u*R04面通量，保證格式的物理一致性激波管實驗設備——Riemann問題的物理實現FINITEVOLUMEMETHOD·CFD117Godunov格式：迎風型方法的奠基之作Godunov格式（1959）通過在每個界面求解精確Riemann問題來計算數值通量，是迎風型有限體積法的開創性工作。格式天然滿足熵條件但僅一階精度，高精度擴展需結合MUSCL重構，計算效率受限于精確Riemann求解器的迭代開銷。01核心步驟將單元平均值視為分段常數初值→在每個界面構造局部Riemann問題→求精確解在x/t=0處的值→計算界面通量02天然滿足熵條件精確Riemann求解器自動排除膨脹激波等非物理解，無需額外添加熵修正03一階精度限制空間精度僅為一階（分段常數假設），需結合MUSCL重構升級為二階Godunov型格式（即MUSCL-Hancock方法）04計算開銷精確Riemann求解器需要對非線性方程組進行Newton迭代，計算量大，尤其對含復雜狀態方程的實際氣體問題NumericalFlux·FiniteVolumeMethodRoe近似Riemann求解器Roe求解器（1981）通過構造滿足特定性質的Roe平均矩陣，將非線性Riemann問題線性化后用特征分解精確求解。格式具有高分辨率和良好的計算效率，是工程CFD中應用最廣泛的數值通量方法之一。Roe性質—構造Roe平均矩陣?(uL,uR)，滿足相容性（uL=uR時退化為Jacobian）、守恒性（通量差等于矩陣乘狀態差）、雙曲性三個條件通量公式—F=?[f(uL)+f(uR)]??|?|(uR?uL)，第一項為中心通量，第二項為基于特征分解的迎風耗散項特征分解—|?|=R|Λ|R?1，R為右特征向量矩陣，Λ為特征值對角陣，耗散按各特征波分量獨立施加分辨率優勢—對接觸間斷的分辨率顯著優于L-F和Rusanov格式，因接觸間斷對應的特征值耗散項較小CFD·有限體積法·LECTURE117Roe格式的熵修正與魯棒性分析Roe格式在跨聲速膨脹區可能因特征值過零而產生不滿足熵條件的膨脹激波。通過Harten-Yee等熵修正方法在特征值接近零時引入適當耗散，可有效消除這一缺陷。極端條件下的正定性保持是保證格式魯棒性的關鍵。熵修正EntropyFix01Roe格式中特征值|λ|在跨聲速點附近趨近于零，耗散項消失，可能產生違反熵條件的膨脹激波02Harten-Yee修正：當|λ|<ε時用ε2+λ2/(2ε)代替|λ|，在零點附近保持最小耗散量ε，消除膨脹激波03修正參數ε選取需謹慎：過小無法完全消除膨脹激波，過大則增加不必要的數值耗散降低接觸間斷分辨率魯棒性增強策略01極端流動條件下（高馬赫數、強激波相互作用）Roe平均可能產生負壓力/密度，需配合positivity-preserving限制器02HLLC等替代求解器在保持接觸間斷分辨率的同時天然保證正定性，成為近年工程應用的熱門選擇Harten-YeeEntropyFixλ?=|λ||λ|≥ε/ε2+λ2─────────2ε|λ|<εMinimumDissipationatSonicPointεWaveDecompositionRoe格式的新理解：特征波分解視角從特征波分解的角度重新理解Roe格式：通量差被分解為各特征波分量的疊加，耗散沿每個波的傳播方向獨立施加。01通量差的特征分解f(uR)?f(uL)=RΛR?1(uR?uL)=Σkαkλkrk每個分量對應一個特征波的貢獻，將非線性通量差線性化為特征模態的疊加02波強度與迎風耗散αk=lk·(uR?uL)，Roe耗散項對每個波分量乘以|λk|實現迎風處理，保證數值耗散的方向與信息傳播方向一致03物理意義：波傳播方向正特征值（右行波）耗散來自左側狀態，負特征值（左行波）耗散來自右側狀態符合雙曲型方程的信息傳播特性，確保數值解滿足熵條件04統一框架：近似Riemann求解器不同格式的區別在于如何處理特征波分解和耗散施加方式為設計新型求解器提供理論基礎，可擴展至高分辨率格式和隱式方法APPROXIMATERIEMANNSOLVERSHLL與HLLC近似Riemann求解器HLL求解器通過雙波假設實現極強的魯棒性但犧牲了接觸間斷分辨率；HLLC在此基礎上增加中間波恢復接觸間斷分辨率并保持正定性，成為現代工程CFD中最受歡迎的數值通量方法之一。HLL求解器01假設Riemann解由兩個波（速度SL和SR）分隔為三個區域，中間僅一個恒定狀態，極大簡化了求解過程02通量公式簡潔且無條件保證正定性，適合高馬赫數流動，但接觸間斷和剪切波被完全耗散，分辨率低于Roe雙波結構HLLC求解器01在HLL的雙波結構中增加一個中間波S*（接觸間斷速度），恢復三個中間狀態，顯著提升接觸間斷分辨率02天然保證正定性（密度和壓力非負），在強激波相互作用等極端條件下比Roe更魯棒，已成為主流工程CFD軟件的默認選擇三波恢復FiniteVolumeMethod·FluxSchemes數值通量方法特性系統對比不同數值通量方法在精度、分辨率、魯棒性和計算效率四個維度各有優劣。工程實踐中需根據流動特征選擇合適的方法：低速不可壓縮流適合中心格式+人工耗散，高速可壓縮流推薦HLLC或Roe+MUSCL組合。常見數值通量方法特性對比通量方法基礎精度接觸間斷分辨率魯棒性熵條件計算開銷Lax-Friedrichs一階差（過度耗散）極強自動滿足最低Rusanov一階較差極強自動滿足低Roe一階*優中等（需熵修正）需修正中等HLL一階差（單中間態）極強自動滿足低HLLC一階*優強自動滿足中等ExactRiemann一階*最優強自動滿足最高*配合MUSCL/WENO重構可提升至二階或更高階精度；HLLC在魯棒性和分辨率之間取得了最佳平衡FluxSplittingEuler方程的通量分裂方法通量分裂通過將Euler方程通量分解為正負兩部分實現迎風離散，是構造迎風格式的另一重要途徑。01Steger-Warming分裂基于Jacobian矩陣特征值符號將f=A?u+A?u，A?對應正特征值（右行波），A?對應負特征值（左行波）A?u+A?u02vanLeer分裂構造關于馬赫數的光滑分裂函數，在聲速點（M=1）處連續可微，消除了Steger-Warming分裂在聲速點的非光滑缺陷M=103Roe特征投影分裂基于Roe平均矩陣的R|Λ|R?1分解實現通量差分裂，物理意義最清晰，推薦使用Roe平均R|Λ|R?104與Riemann求解器通量分裂本質上是Riemann求解器的另一種表述形式，兩者在數學上可以互相轉化等價轉化CFD·有限體積法隱式時間推進：LU-SGS方法LU-SGS方法通過近似LU分解將隱式方程組的求解轉化為前后兩次掃描推進，在保持隱式方法大時間步長優勢的同時避免了求解大型聯立方程組。LU-SGS方法原理線性化：引入差量q=Δu將隱式方程線性化為(I+Δt·?A/?x)q=?Δt·?f/?x≡RHS，避免直接求解非線性方程組近似分解：對系數矩陣進行LU分解(L+D+U)≈(L+D)D?1(D+U)，L步前向掃描q*=D?1(RHS?Lq*)，U步后向修正效率提升：單側差分使每次掃描只需代入已知值即可推進，無需解聯立方程組，計算量從O(N3)降至O(N)工程應用要點大CFL收斂：CFL數可取到10–100量級，相比顯式方法（CFL<1）定常問題收斂速度提升一個數量級以上通量分裂：A矩陣有正有負時強行單側差分會不穩定，需利用通量分裂A=A?+A?分別進行迎風處理CHAPTER04多維迎風型有限體積法與坐標旋轉從一維到多維的關鍵技術跨越與工程實現FINITEVOLUMEMETHOD多維有限體積法的基本框架多維有限體積法將一維方法自然推廣到二維和三維：對每個控制體積的所有界面分別計算法向數值通量。界面通量計算可轉化為沿法線方向的一維Riemann問題，使得一維Riemann求解器可直接應用于多維場景。01二維積分形式—dū/dt=-(1/|Ω|)ΣkFk·nk·lk，對控制體積Ω的每條邊k計算法向通量F·n并乘以邊長lk。該形式將面積分轉化為沿各邊的線積分，是有限體積法的核心離散框架。02旋轉不變性—Euler方程在坐標旋轉下保持不變，界面法向通量可轉化為沿法線方向的一維Riemann問題求解。這一性質保證了多維問題求解的數學嚴謹性。03局部坐標系分解—每個界面定義法向n和切向τ，將速度分解為vn和vτ，法向分量參與一維Riemann求解。通過局部坐標變換，將多維問題降維處理。04切向速度被動輸運—切向速度在Riemann求解中作為被動標量輸運，不影響法向波系結構，求解完畢后再轉回全局坐標系。這種處理方式大幅簡化了計算復雜度。FINITEVOLUMEMETHOD坐標旋轉技術的數學實現坐標旋轉通過正交變換矩陣將多維界面通量問題轉化為法向一維Riemann問題。旋轉矩陣的正交性保證了變換的高效性和數值精確性。01旋轉矩陣與速度分解旋轉矩陣T將全局速度(u,v)分解為法向分量v?=un?+vn?和切向分量vτ=?un?+vn?，密度和壓力保持不變。v?·vτ02通量計算三步流程U′L=TUL，U′R=TUR→F′?=Riemann(U′L,U′R)→F·n=T?1F′?，三步完成通量求解。3Steps03正交矩陣高效性T為正交矩陣（T?1=T?），逆變換僅需轉置運算，計算開銷極低，且保證變換過程的數值精確性。T?1=T?04三維空間推廣需要兩個切向分量，旋轉矩陣擴展為3×3速度子塊，但坐標旋轉與法向通量求解的基本思路完全一致。3×3CFD·有限體積法多維網格類型與工程實踐多維有限體積法支持結構化、非結構和混合三類網格。結構化網格實現簡單但幾何適應性差，非結構網格適應復雜外形但數據結構復雜，混合網格結合兩者優勢。網格類型選擇01結構化網格（四邊形/六面體）：拓撲關系明確、內存訪問連續、計算效率高，但復雜幾何外形下網格生成困難HighEfficiency02非結構網格（三角形/四面體）：自動生成容易、幾何適應性強，但鄰居關系需額外存儲、內存訪問不連續FlexibleGeometry03混合網格：邊界層用棱柱層/六面體保證壁面法向分辨率，遠場用四面體適應外形，工程中應用最廣BestPractice梯度計算與網格質量01非結構網格上梯度計算采用最小二乘法或Green-Gauss方法，精度和魯棒性受網格質量直接影響LSQ/Green-Gauss02網格質量關鍵指標：縱橫比（aspectratio）、偏斜度（skewness）、最小角度，差質量網格導致重構精度嚴重退化Aspect·Skewness·AngleFINITEVOLUMEMETHOD高階多維有限體積方法高階多維有限體積法通過k-exact重構或DG方法實現三階及以上空間精度。k-exact重構基于擴展模板的最小二乘擬合，DG方法在單元內獨立存儲高階多項式。兩者在精度、效率和實現復雜度上各有取舍。01k-exact重構構造k階多項式p(x,y)使得對每個鄰居單元的積分平均值等于已知值，超定系統用最小二乘求解。最小二乘02三階重構（k=2）需要確定6個系數，至少6個約束條件（自身+5個鄰居），模板選取和條件數控制是關鍵難點。6個系數03DG方法每個單元獨立存儲k階多項式的(k+1)(k+2)/2個系數，無需跨單元重構，但存儲量和計算量成倍增加。(k+1)(k+2)/204工程現狀二階MUSCL+限制器仍是主流選擇，高階方法主要用于氣動聲學、大渦模擬等對精度要求極高的領域。二階MUSCLBOUNDARYCONDITIONS邊界條件處理與GhostCell方法有限體積法通過GhostCell（虛