OpenFOAM

OpenFOAM 中的边界条件（一）

本系列解读 OpenFOAM 中边界条件的实现。主要关心一些几个问题：

  1. OpenFOAM 中边界条件是怎样与有限体积离散部分交互的？
  2. 怎么从代码看懂一个边界条件具体是怎么计算边界上的值的？
  3. 怎么定制一个边界条件？

本篇先阐述第一个问题。

从有限体积离散的角度来看，离散过程中，可能要用到的边界信息包括两类：一是某个场在边界上的值，另一是某个场在边界上的梯度。前者在对流项的离散中需要用到，举例说
$$
\int_v \nabla \cdot (\rho \mathbf{U} \phi) dV = \sum_f m_f \phi_f
$$
当组成体积元的面中有边界面时，需要用到这些边界面上的值 $\phi_f$。
而在扩散项的离散过程中
$$
\int_v \nabla \cdot (\Gamma \nabla \phi) dV = \sum_f (\Gamma \nabla \phi)_f \cdot \mathbf{S}_f
$$
这时，如果组成体积元的某个面是边界面，就需要该边界面上 $\phi$ 的梯度值 $\nabla \phi_f$ 了。

边界上的某个场值，或者梯度值，其计算方法可以用如下通式表示
$$
\begin{align}
\phi_f & = A_1 \phi_C + B_1 \\
\nabla \phi_f & = A_2 \phi_C + B_2
\end{align}
$$
这里，$ \phi_C$ 表示 $\phi$ 在邻近边界的网格中心的值，$A_1$ ，$B_1$，$A_2$，$B_2$ 是系数。

OpenFOAM 中的边界条件类中，有四个函数分别对应上面四个系数：valueInternalCoeffs 对应 $A_1$，valueBoundaryCoeffs 对应 $B_1$，gradientInternalCoeffs 对应 $A_2$，gradientBoundaryCoeffs 对应 $B_2$。
所以，看懂 OpenFOAM 中的边界条件，很关键的一步就是看懂这四个函数的定义。

此外，还有一个函数，updateCoeffs， 也很重要。这个函数负责对边界条件进行显式地更新。浏览一下 OpenFOAM 边界条件的代码，会发现很多边界条件都是在 updateCoeffs 这个函数中进行边界值的计算的。

另外，还有些边界条件，似乎是在 evaluate 函数中进行边界值的指定的。

至于边界条件是在什么地方调用的，cfd-online 上有一个帖子，Hrvoje Jasak 对这个问题的回答是：
“
Easy:
- on correctBoundaryConditions() for a field
- on updateCoeffs() at matrix creation
correctBoundaryConditions is also called after the linear solver call automatically.
”
其他网友还提供了一些有价值的信息，比如，在 correctBoundaryConditions 函数中，

1
2
3
4
5
6
7
8
template<class Type, template<class> class PatchField, class GeoMesh>
void Foam::GeometricField<Type, PatchField, GeoMesh>::
correctBoundaryConditions()
{
    this->setUpToDate();
    storeOldTimes();
    boundaryField_.evaluate();
}

调用了 evaluate 函数。
为什么有些边界条件用 updateCoeffs()，而有些则用 evaluate() 呢？目前的理解是这样的： updateCoeffs() 主要用来显式地计算并更新变量在边界上的值，当边界上的值是通过某个依赖于外部参数的公式来计算，并且值会随着迭代的进行而不断改变时，则需要用 updateCoeffs()。有些边界条件，比如最基本的 zeroGradient()，不需要外部的参数，只需要每一次将临近网格的值赋给边界就可以了，这时就可以用 evaluate()。

关于边界条件调用的具体过程，需要在看了 fvMatrix 类以后才能更深入地理解，博主目前只能给出一个粗浅的理解。

参考资料：
The Finite Volume Method in Computational Fluid Dynamics: An Advanced Introduction with OpenFOAM® and Matlab®