FVM in CFD 学习笔记_第6章_有限体积网格

最新推荐文章于 2023-02-13 08:52:29 发布

梅冠华

最新推荐文章于 2023-02-13 08:52:29 发布

阅读量6.7k

点赞数 18

分类专栏： FVM in CFD 学习笔记文章标签： FVM CFD 非结构网格拓扑关系

本文链接：https://blog.csdn.net/meiguanhua/article/details/104745225

版权

本文详细介绍了有限体积法在计算流体力学（CFD）中的网格划分，包括结构网格和非结构网格的特点。结构网格具有规整的拓扑关系，适合简单几何形状，而非结构网格则适用于复杂几何区域，具有更高的适应性。在计算过程中，需要存储几何信息（如单元体积、面积矢量）和拓扑关系（如单元-面、面-角点的对应关系）。通过面形心处的变量值和面积矢量，可以计算单元梯度。在OpenFOAM中，使用非结构网格，但同时也支持结构网格。

摘要由CSDN通过智能技术生成

学习自F. Moukalled, L. Mangani, M. Darwish所著The Finite Volume Method in Computational Fluid Dynamics - An Advanced Introduction with OpenFOAM and Matlab
Chapter 6 The Finite Volume Mesh

划分网格是有限体积法的首要工作，即，把原本连续的计算域划分成一系列不重叠的单元集合，这些单元由一系列面构成（对于每个面来说，其分属于两个单元），同时还要标记边界面来定义物理边界。随后，还要计算网格相关的几何信息（比如单元体积、面积矢量等），以及拓扑关系（单元-面，面-角点的编号关系）。这些是如何处理的呢？本章便着重阐明有限体积法网格中的几何信息和拓扑信息。

1 区域离散

在这里插入图片描述如图所示，对于（a）所示的计算域，可以用结构网格剖分（b），也可用非结构网格剖分（c），将原本控制方程在连续区域上的解转化为控制方程在离散区域上的解，最终获得的是离散点上的变量值（即，有限个离散点处的值）。

结构网格的拓扑关系是规整的，所以无需额外存储，其编程较为简单、计算效率较高，然而其适用的几何计算域较为苛刻，并不是任何计算域都能很容易地划分成结构网格的，当然，可以用分块结构网格的剖分方式来提高其几何适应性，但是工作量很大，且难度较大。

非结构网格的拓扑关系比较复杂，需要额外存储（单元-面，面-角点的关系表），且编程较为复杂（需要考虑网格的非正交修正），计算效率也比较低，然而非结构网格的优点也很明显，其对于复杂几何区域的适应性很强，其划分网格的过程非常方便快捷。

实际上，非结构网格是更一般的情形，而结构网格则是非结构网格的一种特例而已，因此非结构网格的所有计算是兼容结构网格的。

在OpenFOAM中，使用的是非结构网格，但是也能使用分块结构网格和分块嵌套结构网格。下面就先讲讲简单的结构网格，然后再讲讲复杂得非结构网格。

2 有限体积网格

以梯度计算为例，看看结构网格和非结构网格都需要存储什么几何量和拓扑关系。

梯度计算所需的网格信息
在这里插入图片描述
如图，变量一般存储在单元形心 $C$ 处，若采用Green-Gauss方法来计算有限体积单元形心处的变量梯度值 $\nabla\phi_C$ ，有
$\nabla\phi_C=\frac{1}{V_C}\sum_{f=1}^{nb(C)}\phi_f \vec S_f$
即，用面形心处的变量值 $\phi_f$ 与面积矢量 $\vec S_f$ 复合后，加和，再除以单元体积 $V_C$ ，便可得到单元形心处的变量梯度值了。

而面形心处的变量值 $\phi_f$ 则一般用该面两侧单元形心处变量值的加权平均来获取，即
$\phi_f=g_F\phi_F+g_C\phi_C$
一种计算权重系数的方法是
$g_F=\frac{V_C}{V_C+V_F} \quad\quad g_C=\frac{V_F}{V_C+V_F}=1-g_F$
这里只是做简要叙述，梯度的具体计算方法（含非正交修正方法），权系数的其它插值方法将在后续章节中详解。

可见，要存储的几何信息和拓扑关系还真不少呢，有

单元形心坐标值、面形心坐标值、面积矢量值、单元标识与组成单元的面标识的对应关系、面标识与面左右两侧单元标识的对应关系，等

3 结构网格

结构网格的几何信息和拓扑信息要简单明了得狠，以下图为例
在这里插入图片描述
在结构网格中， $i, j, k$ 方向分别划分 $N_i,N_j,N_k$ 个单元，则单元总数为 $N_i*N_j*N_k$ ，且直接可用标识 $(i, j, k)$ 索引到特定的单元，而特定单元的邻近单元可由标识 $(i, j, k)$ 加或减1来寻访到，至于邻近单元的数目，3维情况单元为六面体单元，每个单元有6个面、8个角点、有6个邻近单元，2维情况单元为四边形单元，每个单元有4个面、4个角点、有4个邻近单元。由于拓扑关系明晰，所以无需额外存储拓扑信息。

3.1 拓扑信息

在这里插入图片描述
在计算的时候，最终还是要生成整体矩阵来做计算的，那么局部标识 $(i, j)$ 与整体标识 $n$ 是如何关联的呢？其实很简单，如果整体标识是按照一行一行来标记变量的话，有
$n=i+(j-1)N_i\quad\quad 1\le i \le N_i\quad\quad 1 \le j \le N_j$
那么，单元 $(i, j)$ 邻近单元的局部标识和整体标识为
$\rightarrow n \\ (i+1, j) \rightarrow n+1 \\ (i-1, j) \rightarrow n-1 \\ (i, j+1) \rightarrow n+N_j \\ (i, j-1) \rightarrow n-N_j$
如果是3维问题，则
$n=i+(j-1)N_i+(k-1)N_iN_j \quad\quad 1\le i \le N_i \quad\quad 1 \le j \le N_j \quad\quad 1 \le k \le N_k$
单元 $(i, j, k)$ 邻近单元的局部标识和整体标识为
$\rightarrow n \\ (i+1, j, k) \rightarrow n+1 \\ (i-1, j, k) \rightarrow n-1 \\ (i, j+1, k) \rightarrow n+N_j \\ (i, j-1, k) \rightarrow n-N_j \\ (i, j, k+1) \rightarrow n+N_iN_j \\ (i, j, k-1) \rightarrow n-N_iN_j$
非常方便，组装整体矩阵的时候也很快捷，因为在整体矩阵中的存储位置一目了然，不烦取寻找对应关系了。

3.2 几何信息

在这里插入图片描述访问某特定单元的局部几何信息较为简便，对单元 $(i, j)$ 来说，其构成面为 $\vec {S1}(i,j)$ ， $\vec {S2}(i,j)$ ， $\vec {S1}(i+1,j)$ ， $\vec {S2}(i,j+1)$ ，由于构成单元的面矢量是向外为正的，所以有
$\vec S_{i-1/2,j}=-\vec {S1}(i,j) \quad\quad \vec S_{i,j-1/2}=-\vec {S2}(i,j)$