计算不可压缩流体 -- 差分/有限体积离散格式

不可压缩NS方程求解选用速度，压力为自变量。两者通过守恒方程耦合。

1 自然格式: 速度分量， 压力均在网格节点。 边界条件直接作用在节点上，完美自然边界（Dirichilet B. C)， 问题是满足速度无散度，将在计算中引入非物理的速度。

2 交错格式(staggered grid): 速度分量在单元边界中点， 压力值在单元中心；无滑移条件不能精确满足。如果将因变量通量（质量，动量）移动到边界中点，又导致因变量通量不精确满足。

3 局部交错格式： 速度分量在网格节点，压力值在单元中心； 满足无滑移条件，通量在单元边界守恒。压力Dirichilet边界不精确满足，计算复杂，举例： ALE

4 co-locate格式： 速度分量，压力都存在单元中心；速度，压力边界均不精确满足，但是实现容易，商用软件用的多；

5 multi-momentum交错格式： 速度分量在单元两方向都计算，压力值位于单元中心。每个速度分量在每个单元边界面上都计算，所以精确满足无滑移/通量守恒，单色计算量较（2）大一倍。



一 自然差分格式举例 （2D NS, no body force）

u/t + u du/dx + v du/dy = - grad_p_x + miu * laplace(u)

时间差分： backward Euler ; 空间差分： 中心差分

-- > 非紧致5对角矩阵 A u^(n+1) + B u^n = f^n <--- 时间分裂推进求解速度

压力波松方程(PPE) laplace(p) = -2 ( du/dy * dv/dx - du/dx * dv/dy)

考虑齐次Newman边界 及边界不可穿透条件，上式可解。

关于压力边界： Newman问题，不具有唯一解。一般采用给任意初始压力边界，即PPE实际上求解了 Dirichilet - Newman 混合问题，至于这个初始任意压力，对后续计算是没啥影响的。

自然差分格式算法缺点：

1 没有保证速度无散度条件，PPE简化中使用了 div_v = 0 , 但 动量方程本身没有求解 div_v = 0. 对于初值满足 div_v = 0 在较短的计算时间内可以满足 div_v =0，但是长时间还算没保证啊。

2 齐次压力边界只有在 miu -> 0 才满足

3 中心差分格式诱发速度，压力值的棋盘分布，满足div_v = 0 的伪速度。

4 div-stability 条件 ??

差分CFD基本不用自然格式。

二 交错网格格式（FV)

使用交错网格，2D速度分量和压力共采用三套单元计算，分别即作u单元，v单元 和压力单元。 u，v 单元中心定义单元平均速度， 边线中点定义压力。物理意义上，压力作为一类面力，理应定义在面上。

考虑x分量

LHS 时间导数项 = d ( u_(i, j-1/2) * hx * hy ) /dt % u_(i, j - 1/2) 即单元平均速度

LHS 对流项 = ( u^2_(i+1/2, j-1/2) - u^2_(i-1/2, j-1/2) ) * hy + ( v_(i,j)*u_(i,j) - v_(i,j-1) * u(i, j-1) ) * hx

其中 u_(i+1/2, j-1/2) = ( u_(i, j-1/2) + u_(i+1, j-1/2) ) / 2 ， 注意 u_(i, j-1/2)即单前单元平均速度, 同理 u_(i-1/2, j-1/2)

u_(i, j) = ( u(i, j-1/2) + u(i, j+ 1/2)) / 2 %使用y方向平均, v_(i, j) = ( v(i - 1/2, j) + v( i+1/2, j) ) / 2 % 使用x方向平均

压力梯度项 = ( p(i + 1/2, j-1/2) - p(i-1/2, j-1/2) * hy %单元边界

扩散项 = miu * ( ux_(i+1/2, j - 1/2) - ux_(i - 1/2, j - 1/2)) * hy + ( vy_(i, j) - vy_(i, j-1)) * hx )

同理y分量. 整理各项，得到半隐式交错格式有限体积法的NS离散方程。

上述各项， 只有对流项使用了平均（线性插值）, 截断误差分析表明， 该格式整体精度为网格间距的偶数次幂，采用Richardson 外插可进一步提高精度。

另外，交错格式使用了大量半点的值，实际CFD后处理一般都只计算网格点的值，所以还需要回插网格点值，形成后处理数据。

最后，u, v速度边界处理， 对于u边界速度 u(i, 1) = 2 U_0 - u(i, 2); v边界速度 v(1, j) = V_0

最最后， 交错网格下的PPE求解，同自然格式，注意Dirichilet边界。

三 co-locate 差分格式

所有因变量都在单元中心计算。 诸如， Rhie & Chow算法， Zang算法。后续投影法，SIMPLE算法详细之。