有限体积法介绍

有限体积法1 有限体积法基本原理上一章讲到的有限差分法将数值网格的节点上定义为计算节点，并在网格节点上对微分形式的流体基本方程进行离散，用网格节点上的物理量的代数方程作为原PDE 的近似。

在本章所要学习的有限体积法则采用了不同的离散形式。

首先，有限体积法离散的是积分形式的流体力学基本方程：•d q ds ds SS⎰⎰⎰ΩΩ+∇⋅Γ=⋅φφρφn n v(1)计算域用数值网格划分成若干小控制体。

和有限差分法不同的是，有限体积法的网格定义了控制体的边界，而不是计算节点。

有限体积法的计算节点定义在小控制体内部。

一般有限体积法的计算节点有两种定义方法，一种是将网格节点定义在控制体的中心，另一种方法中，相邻两个控制体的计算节点到公共边界的距离相等。

第一种方法的优点在于用计算节点的值作为控制体上物理量的平均值具有二阶的精度；第二种方法的好处是在控制体边界上的中心差分格式具有较高的精度。

积分形式的守恒方程在小控制体和计算域上都是成立的。

为了获得每一个控制体上的代数方程，面积分和体积分需要用求面积公式来近似。

2 面积分的近似采用结构化网格，在二维情况下，每一个控制体有4个面，二维情况，每一个控制体有6个表面。

计算节点用大写字母表示，控制体边界和节点用小写字母表示。

为了保证守恒性，控制体不能重叠，每一个面都是相邻两个控制体的唯一公共边界。

控制体边界上的积分等于控制体个表面的积分的和：∑⎰⎰=kS Skfds fdS(2)上式中，f 可以表示n u ρφ或n∂∂Γφ。

显然，为了获得边界上的积分，必须知道f 在边界上的详细分布情况，这是不可能实现的，由于只是计算节点上的函数值，因此必须采用近似的方法来计算积分。

整个近似过程分成两步第一步：用边界上几个点的近似积分公式第二步：边界点上的函数值用计算节点函数值的插值函数近似 面积分可采用以下不同精度的积分公式： 二阶精度积分：e e e e S e Sf S f fds F e≈==⎰(3)上式中e f 为边界中点出的函数值。

近似为方格中心点的值乘以方格的面积。

三阶精度积分：e sene S e S f f fds F e2+≈=⎰ (4)四阶精度积分：e see ne S e Sf f f fds F e64++≈=⎰(5)应该注意的是，采用不同精度的积分公式，在相应的边界点的插值时也应采用相应精度的插值函数。

积分公式的精度越高，近似公式就越复杂。

3 体积分的近似和面积分相似，体积分也有不同精度的近似公式 二阶精度积分公式∆Ω≈==⎰P e S q S q qds Q e(6)采用双二次样条函数228272652423210),(y x a xy a y x a xy a y a x a y a x a a y x q ++++++++=(7)可以得到四阶精度的积分公式：()nw ne sw se s n w s P S q q q q q q q q q qds Q e444444441636++++++++∆Ω≈=⎰ (8) 4 函数的插值在上节讲到的积分的近似公式中用到了非计算节点上的函数值，被积函数f 中包含了多个物理量及其偏微分，如对流项n v ⋅=ρφcf,扩散项φ∇⋅Γ=n df，在源项中也有类似情况，这里假定流场和流体的物性参数是已知的，物理量φ及其偏导数在控制面上的值需要通过计算节点上物理量的插值得到。

下面已e 面为例进行讨论。

4.1 迎风插值（UDS ）e φ用上游计算节点的函数值近似相当于对一阶偏导数采用迎风格式，因此用UDS 来表示这种近似方法，在UDS 中：()()⎩⎨⎧<⋅>⋅=00e Ee P e if if n v n v φφφ (9)UDS 是唯一无条件满足有界性要求的近似格式，在数值过程中不会产生数值振荡。

UDS 存在数值粘性。

根据Taylor 公式，该格式具有一阶精度，并具有数值粘性：()2/x u e num e ∆=Γρ(10)在多维问题中，如果流动方向和网格是斜交的，截断误差会在垂直于流动方向以及流线方向产生扩散，这是一种非常严重的误差，函数的峰值或函数值的快速变化会被抹平，为了得到高精度结果需要采用非常精细的网格。

4.2 线性插值（CDS ）P E E E e φλλφφ)1(-+=(11)PE Pe E x x x x --=λ(12)线性插值具有二阶精度，线性插值相当于FDM 中的CDS 格式，因此用CDS 表示。

CDS 格式会产生数值振荡。

对于扩散项PE P Ee x x x --≈⎪⎭⎫⎝⎛∂∂φφφ (13)4.3 三阶迎风格式（QUICK ）和UDS 类似，QUICK 格式也和流动方向有关()()⎩⎨⎧<⋅+-+->⋅+-+-=0)1(0)1(43432121e E EE Pe P W E e ifg g g g if g g g g n v n v φφφφφφφ(14)其中：()We P e Pe W e g ,,2,,111λλλλ-+-=；()()We P e We Pe g ,,2,,2111λλλλ-+--=(15a)()()Pe E e Pe We g ,,2,,3111λλλλ-+-+=；Pe E e Pe E e g ,,,2,41λλλλ-+=(15b)4.4 高阶格式（4阶精度CDS ）采用三次曲线可拟合出四阶精度的中心插值公式，在均匀网格中，四阶公式为：48332727EEW E P e φφφφφ--+=(16)x x EE W P E e∆-+-=⎪⎭⎫⎝⎛∂∂242727φφφφφ (17)5 边界的处理对于对流项，在入口处一般给出了流量或函数值，在边界和对称面上流量为零，在出口处假设和出口的法向坐标无关，因此可采用迎风格式。

对于扩散项则可能需要采用偏心格式。

6 有限体积法应用举例例：考虑一标量在已知流场中的输运过程（如图4.4所示），输运方程为：⎰⎰⋅∇Γ=⋅SSdS dS n n v φρφ(18)边界条件：0=φ；北部入口边界y -=1φ；西部壁面边界对称条件；南部边界 梯度为0；东部出口条件x u x =，y u y -=，流线方程c xy =对流项：e e S c e m dS F eφρφ&≈⋅=⎰n v(19)()y u dS m e x S e e∆=⋅=⎰ρρn v &为质量通量。

⎩⎨⎧+-+=CDSfor )1UDS for )0,min()0,max(E E e P e e Ee P e c e m （mm m F φλφλφφ&&&& (20)若采用UDS 格式，代数方程组中各项系数为：)0,min(e cE m A &=；)0,min(w c W m A &= )0,min(n cN m A &=；)0,min(s c S m A &= (21))(cS c N c W c E c P A A A A A +++-=若采用CDS 格式，代数方程组中各项系数为e e c E m A λ&=；w w cW m A λ&= n n c N m A λ&=；s s cS m A λ&= (22))(cS c N c W c E c P A A A A A +++-=根据连续性方程：0=+++s n w e m m m m &&&&(23)相邻CV 之间的关系：W e P w m m ,,&&-=；W e P w ,,1λλ-=(24)其余相邻CV 有类似关系φ0=φ，扩散项采用CDS 格式()P E P E eS c e x x y y x dS F eφφφφ--Γ∆≈∆⎪⎭⎫⎝⎛∂∂Γ≈⋅∇Γ=⎰n (25)代数方程组中扩散项系数为：P E dE x x y A -Γ∆-=；WP dW x x y A -Γ∆-=P N dN x x x A -Γ∆-=；SP dW x x x A -Γ∆-=(26))(dS d N d W d E d P A A A A A +++-=对于任意控制体P E E N N P P S S W W Q A A A A A =++++φφφφφ (27) d l c l l A A A +=，l 为任意指标P ，E ，W ，S ，N 。

(28)边界条件的处理：对于西部和北部边界，由于给定了函数值，对流项可直接代入函数值而无需插值，扩散项则采用一侧差分WP W PW x x x --=⎪⎭⎫⎝⎛∂∂φφφ (29)这里，W 点和P 的w 边中点重合。

南边和西边的梯度为零，以南边为例，由于梯度为零，S P φφ=，代数方程变为：P E E N N P P S W W Q A A A A A =++++φφφφ)((30)6 SIMPLE 方法考虑定常不可压流动问题，控制方程为： 连续性方程：0=⋅⎰SVdS n v ρ(31)动量方程：⎰⎰⎰⎰Ω+-∇⋅=⋅CVSVSVSVd dS p dS dS b n v n n vv ρμρ(32)不可压缩问题求解的困难在于压力场的求解。

主要原因在于压力p 没有独立的方程组。

先考虑一维问题： 对于动量方程：()()w e we w e p p x u x u uu uu +-⎪⎭⎫⎝⎛∂∂-⎪⎭⎫⎝⎛∂∂=-μμρρ (33)若采用CDS 格式()()()()2222W P E P W P PE WP E P p p p p x u u x u u uu uu uu uu +++-∆--∆-=+-+μμρρρρ 简化后得：()()222W E P w E WE p p x u u u uu uu +-∆-+=-μρρ(34)根据连续性方程，c u u u i i i ===-+11，则有11+-=i i p p ，由于相邻节点之间的压力没有联系方程，容易造成压力交错现象。

为了解决这一问题，可采用交错网格技术，即速度场和压力场采用不同的网格。

以二维问题为例，交错网格的布置如下图所示：主控制体为压力控制体（黑色实线网格），u 的控制体（红色虚线网格）的计算节点在主控制体的e 边，控制体的e ，w 边界通过主控制体的计算节点，v 控制体（蓝色双点划线网格）的计算节点在主控制体的n 边，该控制体的n ，s 面经过主控制体的计算节点。

在u 的控制体中，采用有限体积法离散可得u 的代数方程：∑∑-++=-++=eN P nb nb n n e E P nb nb e e A p p Q v a v a A p p Q u a u a )()( (35)压力场的求解采用压力校正方法。