从二十多年前的公开发行开始，OpenFOAM从根本上颠覆了计算流体力学的世界，其在全世界的学术界和研究中心广泛地传播。 它的成功很大程度上与项目的开源性质和软件的革命性设计有关。 最重要的是，OpenFOAM降低了CFD软件的成本。 由于这个原因，学术界以其丰富的学生劳动，引导了早期OpenFOAM的大部分采用。 开放源码的灵活性允许研究人员以无与伦比的自由进行创新。该行业还对OpenFOAM的定制能力和针对内部问题专门设计的求解器感兴趣。 如今OpenFOAM作为CFD的一个标准为该领域的每个人所熟悉。然而OpenFOAM的所有优势都有一个弱点：缺乏文档及相对陡峭的学习曲线。

从理论上讲，在开源的情况下代码就是文档。一个有经验的用户可以深入到C++的各个层面并理解代码的基本功能。然而对于新用户来说，这种期望是开始使用该软件的巨大障碍。 C++语言和模板化OpenFOAM代码的复杂性加剧了这一困难。即便是像Doxygen这样的工具（它试图对OpenFOAM类结构进行编目）也不能充分减少这些障碍。在用户论坛上，一个沮丧的人哀叹说，对一个C++初学者来说，挖掘Doxygen的输出就像倒着阅读中文。OpenFOAM社区非常需要一种全面的且易于访问的文档形式，为世界各地的OpenFOAM新手提供一个入口。因为OpenFOAM是免费的，学习曲线是OpenFOAM用户进入的根本障碍。

人们可以以合适的价格从OpenFOAM专家那里购买培训，但一流的参考书依然是无可替代的。 多年来，没有任何参考文档为新的OpenFOAM用户提供深度和广度的支持。 这一切在2014年发生了变化，当时Tomislav Marić, Jens Höpken, 及Kyle Mooney出版了《OpenFOAM Technology Primer》。 这本书对那些想要从一个单一的来源可以理解OpenFOAM的人来说是一笔意外之财。 该书的目标受众是那些对CFD有基本背景，但希望更深入地研究OpenFOAM的工作原理和使用方法的人。 尽管我有近十年的OpenFOAM经验，但那本书成了我研究实验室的圣经，就在我的桌子上触手可及。 但在取得巨大成功的几年后，这本书完全消失了。 社区渴望这个可爱的参考资料的回归--直到2021年。

我们现在很幸运，作者们投入了更多的时间来更新和重新发布它。 编写OpenFOAM文档的一个根本挑战是，它是一个不断发展的软件，不断地改变和改进。保持这样一个文档是一个永远做不完的任务，需要一个无休止的承诺。出于这些原因，我很激动地看到作者选择在知识共享许可证下发布文档。就像开源原始OpenFOAM代码的决定一样，这本有价值的书将引领OpenFOAM使用和文档编制的新时代。

Professor David P. Schmidt

Dept. of Mechanical and Industrial Engineering

University of Massachusetts Amherst

Tomislav Marić

Tomislav在克罗地亚萨格勒布大学学习机械工程，并于2017年在德国达姆施塔特工业大学数学系Dieter Bothe教授领导的数学建模与分析研究所（MMA）获得博士学位。在德国达姆施塔特工业大学数学系Dieter Bothe教授的领导下，于2017年获得博士学位。Tomislav目前在达姆施塔特工业大学作为雅典青年研究员工作（2020年10月）。Tomislav自2008年以来一直在开发非结构化的拉格朗日/欧拉界面近似（LEIA）方法，用于OpenFOAM开源软件中的两相流模拟。他与Jens共同创立了sourceflux，自2016年以来，作为达姆施塔特工业大学合作研究中心（CRC）1194的成员，Tomislav在科学软件开发和研究数据管理方面支持CRC-1194的研究人员。

Jens Höpken

Jens在杜伊斯堡-埃森大学学习海军建筑学，毕业于德国杜伊斯堡-埃森大学的船舶技术、海洋工程和运输系统研究所（ISMT）。Jens从2007年开始与OpenFOAM合作，他是一位OpenFOAM专家，在为海军水动力学开发OpenFOAM方面有超过十年的经验。Jens与Tomislav共同创立了sourceflux。

Kyle G. Mooney

Kyle于2016年在马萨诸塞大学阿默斯特分校获得机械工程博士学位。他的研究涉及粘弹性流体和喷雾液滴动力学的数值模拟。研究生毕业后，他加入了ICON技术流程与咨询公司，帮助福特汽车公司、菲亚特克莱斯勒汽车公司和通用汽车公司创新汽车空气动力学模拟流程。在搬到旧金山后，他转向领导流体机械研发工作，同时在几个初创公司进行硬件和软件产品开发。他目前是Geminus.AI的高级应用工程师，为工业应用创建高保真度流动系统仿真模型。在他的空闲时间，你可能会发现他在打拳击、玩滑板、在高山上徒步旅行，或者表演音乐。

计算流体动力学(CFD)开源软件OpenFOAM被广泛应用于工业界和学术界。相较于使用专有CFD软件，使用OpenFOAM的优势在于开源通用公共许可证(GPL)，其允许用户自由使用和修改现代高端CFD代码。开源许可消除了产品优化周期中的许可成本，实现了参数变化的直接自动化，并加速了新的数值方法和模型的开发。加快了新方法的开发和实现，因为它们是从OpenFOAM中的现有功能开始的，而不是从头开始的。

除了上面提到的所有优点之外，使用OpenFOAM还有一个缺点。OpenFOAM使用了C++编程语言和现代软件设计模式，因此学习如何以模块化和可持续的方式在OpenFOAM中编写新方法需要付出大量的努力。使用OpenFOAM需要结合不同背景的知识，包括应用数学、物理、软件开发、C++编程语言和高性能计算(并行编程和性能测量)等。

这本书致力于在一个地方描述OpenFOAM的不同方面，帮助初学OpenFOAM的用户发展成为中级OpenFOAM程序员。为了实现这一目标，强烈建议读者通读所涵盖的例子。本书没有介绍OpenFOAM的一些核心部分所必须的C++编程语言的更高级部分。但是这些知识可以在有关C++编程语言和软件设计模式方面的书籍中找到。

本书涵盖了使用OpenFOAM的两个主要方面：使用应用程序以及开发和扩展OpenFOAM应用程序和库。本书第一部分使用几个OpenFOAM实用程序和应用程序描述OpenFOAM工作流程，第二部分介绍了OpenFOAM新求解器和库的开发。

1 目标读者

本书面向对开源计算流体动力学(CFD)感兴趣的任何人。

但是不可能在一本书中提供有效开发OpenFOAM所需的C++编程语言、软件设计、计算流体动力学(CFD)和高性能计算(HPC)的所有背景信息，本书将重点放在OpenFOAM解算器和库的使用、设计和开发上，并引导读者了解那些没有详细介绍的主题所涉及的其他深入信息来源。

因此，假设C++编程语言中的一些面向对象编程知识，涉及类（封装、继承和组合）、虚函数（动态多态性）和运算符重载，这些主题的背景信息在书的第二部分提供。然而，读者也应该使用每章末尾所引用的文献独立地了解这些主题。本书所提供的示例故意避开了验证和确认，因为这会让读者偏离学习OpenFOAM的轨道。然而如果没有计算流体动力学的知识，以及非结构有限体积法（FVM）的知识，就不可能学习和理解OpenFOAM，这本书中简要介绍了该方法，在其他地方也有更详细的介绍。

2 本书涵盖的内容

• 第1章 概述了CFD模拟的工作流程，以及OpenFOAM中的非结构有限体积法（FVM）。
• 第2章介绍了区域离散化（网格生成和转换）及区域分解。
• 第3章描述了仿真实例的结构和设置：设置初始条件与边界条件，配置模拟控制参数及数值参数。
• 第4章概述了前后处理实用程序和数据可视化。
• 第5章对OpenFOAM库进行了深入的概述。
• 第6章描述了如何以高效且可持续的方式编程OpenFOAM：开发和使用库、使用git版本控制系统、调试和分析
• 第7章简要概述了湍流建模：将湍流引入模拟案例并配置湍流模型。
• 第8章介绍了OpenFOAM预处理和后处理应用程序的编程。
• 第9章介绍了OpenFOAM中解算器设计的背景，并展示了如何使用新功能扩展现有求解器。
• 第10章介绍了OpenFOAM中边界条件的数值背景和软件设计。提供了一个自定义边界条件的实现示例，该示例使用了第6章中描述的原则。
• 第11章以温度相关粘度模型为例，介绍了OpenFOAM中输运模型的数值背景、设计和实现。
• 第12章介绍了OpenFOAM中的函数对象的设计与实现，并与C++函数对象进行了比较。
• 第13章介绍OpenFOAM中的动态网格处理。介绍OpenFOAM中动态网格引擎的设计和使用，以及使用动态网格处理扩展解算器。
• 第14章对本书进行了总结。

3 如何阅读此书

对于刚开始使用OpenFOAM的用户，建议从头到尾阅读这本书，并对示例进行独立的研究。有经验的OpenFOAM用户可以从第二部分中选择一章，其中有关于如何对OpenFOAM的某个特定部分进行编程的相关信息。

4 OpenFOAM版本

OpenFOAM有不同的版本可供选择，例如OpenFOAM Foundation、Foam Extended和OpenFOAM。本书没有介绍这些forks之间的区别和相似之处。书中的内容和实例库与OpenFOAM-v2012相匹配，并且该书将只遵循此OpenFOAM版本。

有关如何安装此版本OpenFOAM的信息，请访问其官方网站。可以选择将OpenFOAM作为Linux包安装，编译源代码的快照，或者编译克隆的git存储库。由于本书的目标是解决OpenFOAM编程问题，因此建议编译源代码快照或OpenFOAM git存储库的克隆。

5 命名和排版约定

命令行使用带有前缀?> 的输入。下面是一个示例：

?>  ls $FOAM_TUTORIALS
Allclean basic electromagnetics lagrangian
Allrun combustion financial mesh
Alltest compressible heatTransfer multiphase
DNS discreteMethods incompressible stressAnalysis

C++代码排版形式如下所示：

template<class GeoMesh> 
tmp<DimensionedField<scalar, GeoMesh>  >  stabilise
(
    const DimensionedField<scalar, GeoMesh> &,
    const dimensioned<scalar> &
);

OpenFOAM模拟的配置依赖于所谓的字典文件。字典文件是一种文本文件，以OpenFOAM特有的格式存储键值对的列表：

ddtSchemes
{
    default Euler;
}

在方程式中，标量(如通量$\varphi$)的排版没有强调；矢量的排版使用粗体(例如，速度$\mathbf{U}$)；张量使用粗体和下划线(如单位矩阵$\underline{\mathbf{E}}$)。

6 随书示例

书中包含的示例可以在GitLab上找到: https://gitlab.com/ofbook-/ofprimer。

最新版本合并到主分支中，所有版本都与OpenFOAM git标签版本相匹配，如OpenFOAM-v2012。

7 贡献

错误报告和主题建议使用GitLab服务台处理，并且可以投票表决。支持率最高的话题将有更大的机会在下一版中得到解决。

在提交错误报告或功能请求之前，请搜索existing issues以确保该问题尚未被报告。

本章概述了在OpenFOAM中求解计算流体力学（CFD）问题的工作流程。

此外本章还包含了OpenFOAM所使用的非结构化网格上的有限体积法（FVM）的背景信息，以及对OpenFOAM平台顶层结构的描述。

CFD 分析的目标是更深入地理解所考虑的问题。由于模拟结果通常伴随着实验数据，因此需要考虑结果的有效性和模拟的整体目标。此外，模拟的物理过程需要得到适当的数学建模。为此，CFD 工程师需要在 OpenFOAM 框架中选择合适的求解器。在 CFD 分析过程中，假设可能会使分析变得更加复杂或简化。通常会对模拟域的几何形状进行分析和简化，考虑到可用的计算资源、所需的模拟保真度和期望的运行

任何 CFD 分析的目的都是为了更深入地了解所考虑的问题。由于仿真结果往往伴随着实验数据，因此需要考虑结果的有效性和仿真的总体目标。此外需要对模拟的物理过程进行适当的数学建模，这又要求 CFD 工程师在 OpenFOAM 框架内能够正确地选择合适的求解器。在 CFD 分析过程中，任何假设都可能会使分析变得复杂或简化。考虑到可用的计算资源、所需的仿真保真度以及仿真所需的周转时间，通常会对仿真区域的几何形状进行分析和简化。

下面概述了在进行 CFD 分析项目之前可能提出的实用问题。

1、一般考虑

• CFD 分析得出的结论应该是什么？
• 如何定义计算结果的精度程度？
• 用什么方法来验证计算结果？
• 这个项目有多少时间可用？

2、热物理问题

• 流动是层流、湍流还是转捩流动？
• 流动是可压缩的还是不可压缩的？
• 流动是否包含多种流体相或化学组分？
• 传热在问题中是否起着重要作用？
• 材料属性是否为自变量的函数？例如剪切稀化流体。
• 边界条件和初始条件方面是否有足够的信息？它们是否有适当的模型，是否可以准确地进行近似？

3、几何与网格

• 能否建立流体域的精确离散形式?
• 在模拟过程中，计算域是否会变形或移动?
• 是否可以在不影响求解精度的情况下降低计算区域的复杂性?

4、计算资源

• 有多少可用的模拟计算时间?
• 可用的分布式计算资源有哪些?
• 一次计算就足够，还是需要多次计算?

这些问题有助于对任何流动问题进行精确的 CFD 分析。使用 OpenFOAM 或任何其他 CFD 模拟软件都需要对物理学、数值方法以及可用的计算资源有适当的理解。CFD 的跨学科性质极大地增加了它的复杂性。

CFD 分析通常可细分为 5 个主要步骤。其中一些步骤有时需要进行多次以获取所需的高质量结果。

1.2.1 问题定义

从工程的角度出发，数值模型应尽可能简单，并能准确地描述实际工程系统。忽略模拟问题的无关方面可以提高 CFD 分析的效率，因为它简化了物理过程，从而简化了描述它的数学模型。例如，尽管空气是可压缩的，但在某些流动状态下，模拟机翼上的流动时可以将空气视为不可压缩的流体。

1.2.2 数学建模

一旦物理过程的相关方面被独立出来，问题就需要以数学模型的形式进行描述，这在 CFD 中通常是一组偏微分方程 (PDE)。CFD 工程师必须了解用来描述不同物理现象的模型。在 OpenFOAM 框架中，用户可以在几十种求解器之间进行选择。每个求解器实现一个特定的数学模型，选择正确的数学模型对于获得模拟问题的有效解通常是至关重要的。例如在翼型绕流的问题中，采用不可压缩假设会忽略对能量方程的求解。作为另一个例子，势流只受拉普拉斯方程控制（详情可见 Ferziger 及 Peric（2002）的教材）。若需要考虑更加复杂的物理输运现象，则数学模型的复杂性也会随之增加，这通常导致更复杂的数学模型。例如用于模拟湍流使用雷诺平均 NS 方程(RANSE)。

数学模型描述了流动的细节，这意味着数值模拟最多只能近似模型的解，无法产生比数学模型本身所能描述的更多的关于流动的信息。有关 OpenFOAM 中湍流建模的更多信息可以在第 7 章中找到。关于特定数学模型的更多细节可以在流体力学教科书中找到。

1.2.3 前处理及网格生成

数学模型将物理场定义为模型方程的因变量。在 CFD 中，方程通常描述的是一个边值和初值问题。因此，在开始模拟之前需要对物理场进行初始设置 (预处理)。如果物理场在空间上变化，则可以使用不同的应用程序 (utility)来计算和预处理。有些应用程序是随 OpenFOAM 一起发布的 (例如 setFields 实用程序)，或者是其它项目的一部分 (例如 swak 4 Foam 项目的 funkySetFields 实用程序)。

一些可用的预处理程序的使用将在第 8 章进行描述。

注：有关 swak 4 Foam 项目的更多信息可以在 http://openfoamwiki.net/index.php/Contrib/swk4Foam上的OpenFOAM wiki 上找到。

为了在数值上近似模型解，必须对模拟区域进行离散化。模拟域的空间离散化包括将区域分割成由不同形状的控制体积 (单元)组成的计算网格。所有这些控制体一起被称为“网格”或计算网格。通常情况下，在感兴趣的区域需要对网格进行细化：例如在出现较大的流场梯度的那些区域。此外，还必须注意数学模型的准确性和正确的选择。以空间的方式解析流动特征并不能补偿一个最初没有考虑这些特征的模型。

另一方面，瞬态模拟中增加网格分辨率可能会极大地降低模拟速度。这是因为当使用显式离散方式时，为了获得稳定的解，离散时间步长通常需要设置成较小的值。如果数值模拟不收敛，则网格可能是模拟中最有可能需要改变的部分。失败的模拟经常是由网格质量不足引起的。OpenFOAM 附带两个不同的网格生成器，即 blockMesh 和 snappyHexMesh。这两种网格生成器的用法都将在第 2 章中介绍。

此外，预处理还包括其他各种任务，例如，如果模拟是在多台计算机或 CPU 核心上并行运行，则要对计算域进行分解。

1.2.4 求解

除了网格生成外，计算求解通常是 CFD 分析中最耗时的部分。

计算求解所需的时间在很大程度上取决于数学模型、用于近似求解的数值格式以及计算网格的几何和拓扑性质。在这一步中，微分数学模型被 (线性化的)代数方程组所取代。在 CFD 中，这样的代数线性方程组通常很大，导致矩阵具有数百万或数十亿的系数。这些代数方程组是用专门为此开发的算法——迭代线性求解器来求解的。

OpenFOAM 框架支持选择多种线性求解器，尽管求解器应用程序通常有预设或理想的线性求解器和参数选择。熟练的 CFD 工程师有机会修改求解器和相应的参数。线性求解器和离散化策略的选择是影响计算速度和稳定性的重要因素。

1.2.5 后处理

在模拟成功完成之后，用户通常会有大量的数据需要分析和讨论。为了检查流动细节，必须适当地提取数据、绘制数据和/或可视化数据。通过使用诸如 paraView 之类的专用工具，可以相当容易地讨论这些数据。为了分析模拟结果，OpenFOAM 提供了广泛的后处理应用程序。

OpenFQAM 自带标准后处理工具 paraView。它是一个开源工具，可以从以下网站免费获得： www.paraview.org 。

1.2.6 验证和确认

用户必须自己决定是否信任计算结果。通常 CFD 软件非常复杂，其依赖于可配置的参数，这给错误留下了很大的空间。如果在前面的步骤中出现了错误，很有可能在验证和确认过程中被发现。

验证（Verification）保证了数值方法能正确地求解它所近似的数学模型。换句话说，验证检查数学模型的解是否被适当逼近。

确认（Validation）将仿真结果与实验数据进行比较：当涉及到仿真结果的置信度时，引入了更严格的安全系数。当与实验进行比较时，验证确保选择了正确的数学模型，其解充分反映了现实。如果模拟结果不满足要求，则必须重新进行之前的 CFD 分析步骤。

本节简要概述OpenFOAM中的有限体积法（Finite Volume Method，FVM）。关于有限体积法的更详细描述，读者可以参考一下相关文献：Ferziger和Peric(2002)、Versteeg and Malalasekra，Weller,Tabor, Jasak及Fureby、Jasak, Jemcov及Tuković、Moukalled, Mangani, Darwish, et al等。

OpenFOAM中非结构有限体积法离散步骤与1.2节中描述的CFD分析步骤有一定的关联。用于定义流体流动的物理属性(如压力、速度或温度/焓等)是数学模型（描述流体流动的偏微分方程组）中的相关变量。

看似不同的物理过程有时可以使用相似的数学模型进行描述，例如热传导过程与水中糖分浓度扩散过程都可以使用扩散模型进行描述。如Ferziger和Peric 2002所描述的通用标量输运方程中包含了用于模拟不同物理过程的项（微分算子），如：流体粒子随流体速度的输运(对流项)、热源(源项)等。由于通用方程中包含了经常遇到的项，因此常用其描述FVM的离散过程。通用输运方程为：

$$\frac{\partial \varphi }{\partial t}+\nabla \cdot (\mathbf{U}\varphi )+\nabla \cdot (D\nabla \varphi )=S_{\varphi }\text{\hspace{1em}(1.1)}$$

式中，$\varphi$ 为待求标量；$\mathbf{U}$ 为指定的速度向量；$D$ 为扩散系数。方程(1.1)中各项从左至右分别为：瞬态项、对流项、扩散项及源项。方程中的每一项都描述了一个以不同方式改变场变量 $\varphi$ 的物理过程。

有时可以根据物理过程的特性忽略一些项：如对于无粘流体流动，可忽略动量的扩散，因此可以将其从动量方程中去掉。此外出现在某些项中的系数可以是恒定值，也可以是空间和/或时间上变化的物理场，或与介质的物性参数有关。

数值方法的目的为了获取数学模型的近似解。对复杂物理过程进行近似求解是有必要的，因为精确求解通常只能在一些非常特殊的情况下才能得到，而这些特殊情况往往难以应用在工程实际中。

数学模型的近似解通过求解控制方程组的离散近似而得到。有限体积法（FVM）的近似过程涉及使用相应的线性代数方程组替代PDE系统，然后使用计算机进行求解。非结构有限体积法（FVM）生成代数方程组包括两个主要步骤：解域离散和方程离散。

1.3.1 解域离散

数学模型使用连续变量。为了逼近数学模型的解，空间被离散成有限数量的控制体积（网格单元）。这些有限体积（单元）构成了有限体积网格。

图1.1显示了从连续流体域$\Omega$到离散流体域${\Omega }_D$的转变。

连续流体域$\Omega$近似为网格单元（有限体积）的集合：

$$\Omega \approx {\Omega }_D={\cup }_{c\in C}{\Omega }_c\text{\hspace{1em}(1.2)}$$

其中，${\Omega }_D$为有限体积网格，$C$是离散区域${\Omega }_D$中所有网格索引的集合。在图1.1a中流体填充的空间的每一点上定义的连续物理场在有限体积${\Omega }_c$内被线性近似，如图1.1b所示。

每个有限体积都存储与其中心${\mathbf{x}}_c$相关的物理属性（如温度）的体积平均值。将该值与网格中心相关联可使区域离散达到二阶精度。为了解为何会出现这种情况，可以假设物理场$\varphi$可以使用泰勒级数展开为：

$$\varphi (\mathbf{x})=\varphi \left({\mathbf{x}}_c\right)+\nabla \varphi \left({\mathbf{x}}_c\right)\cdot \left(\mathbf{x}-{\mathbf{x}}_c\right)+\nabla \nabla \varphi \left({\mathbf{x}}_c\right):\left(\mathbf{x}-{\mathbf{x}}_c\right)\otimes \left(\mathbf{x}-{\mathbf{x}}_c\right)+\dots \text{\hspace{1em}(1.3)}$$

引入在点 ${\mathbf{x}}_c$ 处定义的量（这些量在 ${\Omega }_c$ 上是恒定的），用式(1.3)中的泰勒级数展开表示网格 ${\Omega }_c$ 内部 $\varphi (\mathbf{x})$ 的体积平均值，可得到：

$$\begin{array}{rl}{\varphi }_c& =\varphi \left({\mathbf{x}}_c\right)+\nabla \varphi \left({\mathbf{x}}_c\right)\cdot \frac{1}{\left|{\Omega }_c\right|}{\int }_{{\Omega }_c}\left(\mathbf{x}-{\mathbf{x}}_c\right)dV\\ & +\nabla \nabla \varphi \left({\mathbf{x}}_c\right):{\int }_{{\Omega }_c}\left(\mathbf{x}-{\mathbf{x}}_c\right)\otimes \left(\mathbf{x}-{\mathbf{x}}_c\right)+\dots dV.\end{array}\text{\hspace{1em}(1.4)}$$

体积域${\Omega }_c$的中心${\mathbf{x}}_c$定义为：

$${\int }_{{\Omega }_c}\left(\mathbf{x}-{\mathbf{x}}_c\right)dV=0\text{\hspace{1em}(1.5)}$$

将式 $(1.5)$ 到式 $(1.4)$ 中，可得到：

$${\varphi }_c=\varphi \left({\mathbf{x}}_c\right)+{\int }_{{\Omega }_c}\nabla \nabla \varphi \left({\mathbf{x}}_c\right):\left(\mathbf{x}-{\mathbf{x}}_c\right)\otimes \left(\mathbf{x}-{\mathbf{x}}_c\right)+\dots dV\text{\hspace{1em}(1.6)}$$

由式 $(1.6)$ 可知，对于线性分布的 $\varphi$，由于其高阶导数为零，因此 $\varphi$ 在有限体积 ${\Omega }_c$ 上的平均值正好等于 $\varphi$ 在 ${\Omega }_c$ 的中心 ${\mathbf{x}}_c$ 处的值。换句话说，有限体积中心处的单元平均值(以网格为中心)可以准确地恢复线性物理场的值。能精确恢复线性函数值的方法至少具有二阶精度。

注：在网格${\{{\Omega }_c\}}_{c\in C}$的中心$\{{\mathbf{x}}_c{\}}_{c\in C}$处得到$\varphi$的单元平均值${\varphi }_c$的非结构FVM的区域离散具有二阶精度。

基于式 $(1.6)$ 的变量 $\varphi$ 的插值：

$${\varphi }_c\approx \frac{1}{\left|{\Omega }_c\right|}{\int }_{{\Omega }_c}{\varphi }_l\left({\mathbf{x}}_c\right)dV+O\left({\Vert \mathbf{x}-{\mathbf{x}}_c\Vert }_2^2\right)\text{\hspace{1em}(1.7)}$$

式 $(1.6)$ 具有二阶精度，因为式中泰勒级数截断部分的最大项与 ${\Vert \mathbf{x}-{\mathbf{x}}_c\Vert }_2^2$ 成比例。

在式 $(1.3)$ 中，网格中心处的梯度 $\nabla \varphi ({\mathbf{x}}_c)$ 也必须估计，计算方法在1.3.2节方程离散中会介绍。在理解方程如何离散之前，应当更详细地定义网格 ${\Omega }_D$，特别是网格单元 ${\{{\Omega }_c\}}_{c\in C}$ 的相互连接方式，因为这种连接性决定了哪些区域可以离散化。

我们在此区分三种类型的网格：结构网格、块结构网格及非结构网格。这三种网格对区域和方程的离散有不同的要求，并且网格单元相互连接和寻址方式也不同。网格单元之间的连接(网格连通性)决定了相邻单元的访问方式，这对于方程离散非常重要。这进一步影响了在访问网格单元时可能进行的优化，这些优化反映了数值运算的效率以及如何将这些运算并行化。

例如，网格单元的非结构化寻址使得很难在任何特定方向上访问网格单元。这使得在高阶插值格式所需的非结构化网格上构建更大的模板变得复杂，因为它们依赖于来自更宽邻域的单元平均(单元中心)值。使用域分解和消息传递方法对这种依赖单元中心值的高阶插值进行并行化也很复杂，而且可能效率不高，因为必须跨进程边界传递长度可变的大消息。因此，网格的连通性决定了在网格上可以有效地计算什么，有时会产生这样的效果：即为特定的网格选择错误的数值方法所造成的效率低下，使得该方法无法用于实际问题，而这些问题通常需要更多的网格单元。在OpenFOAM中，网格连通性在数值算法的实现方式中也起着非常重要的作用，OpenFOAM仅支持非结构化网格。

结构网格支持对任意网格单元的邻居进行直接寻址以及对网格单元的直接遍历：单元被标记为沿坐标轴方向递增的索引(参见图1.2a)。另一方面，非结构化网格没有明显的方向(见图1.3a)。

结构网格连通性提高了FVM中涉及的插值的绝对精度，但是当将其用于复杂几何区域的网格生成时，会使网格的灵活性降低。

在网格生成过程中，用户通常希望在解发生较大变化的地方生成密集的网格，而不希望在变量梯度较小的流动区域中浪费网格。所有网格生成步骤应在最短的时间内完成，如果没有实质性的改善，使用结构化网格是无法实现此目的的。对结构网格进行局部网格细化是不可能实现的，因为网格细化必须通过整个网格传播到各自的方向。图1.3b给出了结构网格和非结构网格之间拓扑差异的示意图。

图1.2a显示了一个二维结构网格。这种网格也被称之为笛卡尔网格，因为其由体积中心沿坐标轴方向分布的控制体积组成的。通过改变索引$i,j$在网格中移动具有一个显著的优势：使用这种网格的数值方法可以通过简单地将索引$i$，$j$递增或递减1来访问任何相邻的网格。例如当通过单元面的通量${\varphi }_f$从单元中心插值到面中心时，可以轻松地应用高阶插值格式以提高求解精度。高阶插值格式意味着将更多的网格单元（不一定与当前网格相邻）包括在插值计算中。

但存在一个问题：在结构网格中进行网格局部细化是很难实现的。

当结构网格上可用的插值阶次仍然不够精确时，在极端变化的区域中常会出现这种情况。例如，在激波模拟或两种不相溶流体的两相流模拟中，会出现物性参数值的显著变化。将它们分开的两种流体之间形成了一个界面，物性参数值可能会存在几个数量级的变化（如气-水界面，密度比约为1000），如图1.4中的示意图所示。

为求解这种数值的急剧变化，需要对网格进行局部细化。这种网格细化既可以在前处理期间完成，也可以在求解计算时进行。如前所述，结构网格的细化不能在局部完成：结构化网格的拓扑结构迫使网格在整个方向上进行细化（请参见图1.2b），其中在两个方向上对单个网格的细化会在整个方向上产生细化。符合弯曲几何形状的结构化网格尤其难以生成，为保持结构网格的连通性，必须对弯曲边界进行参数化，这仅适用于相对简单的几何边界。

然而，结构网格离散在实际应用中有一些扩展，其中一些扩展允许对网格进行局部和动态细化。为了提高局部精度，可以使用块结构细化，这是一个构建由多个结构块组成的网格的过程。当组装这样的块结构网格时，不同的块可以具有不同的网格密度。不幸的是，这引入了新的复杂性，因为数值方法必须能够处理不一致的悬挂节点。或者通过精心地处理结构块，以使不同密度的相邻块上的网格节点完全匹配。即使对于简单的流体域，构建块结构网格也是一个复杂的问题，这往往使得在许多涉及复杂几何形状或边界形状的技术应用中，块结构的网格都不是一个好的选择。细化块结构的网格会导致细化区域扩展到整个块中。使用依赖于边界一致性的块结构网格的标准求解器，细化会使得网格生成更加复杂。

OpenFOAM中的动态自适应局部细化是通过引入额外的数据结构来实现的，这些数据结构产生并存储与网格细化过程相关的信息。这种方法的其中一个例子是基于八叉树的细化，其利用八叉树数据结构将每个结构笛卡尔网格元拆分为八个。这种方法要求网格（或至少是正在细化的网格子集）仅由六面体单元组成。数值插值程序（离散微分算子）使用八叉树数据结构存储信息，并考虑局部网格细化导致的网格拓扑变化。然后可以使用切割单元方法将处理更复杂几何域的可能性添加到八叉树精细结构网格中。在这种情况下，通过边界的分段线性近似来切割弯曲域边界的网格。基于八叉树的自适应网格细化在效率方面具有优势，这取决于在底层结构化网格上执行拓扑操作的方式。然而基于八叉树的细化逻辑要求细化域必须是盒形的。更多关于局部自适应网格细化程序的信息可以在[文献](8 “Tomasz Plewa, Timur Linde, and V Gregory Weirs. Adaptive mesh refinement-theory and applications. Vol. 41. Springer, 2005, pp. 3–5.”)中找到。

OpenFOAM实现了一个支持任意非结构网格的二阶收敛的FVM。除了非结构化网格的连通性之外，任意非结构网格的单元可以是任何形状。这允许用户对几何复杂度非常高的流体区域进行离散。非结构网格允许非常快速、有时甚至是自动的生成，这对于工业应用来说非常重要，因为在工业应用过程中需要尽可能快速地获得结果。

图1.3a显示了一个二维非结构网格。由于网格寻址不是结构化的，所以仅出于解释网格连通性的目的对单元进行了标记。无序的网格单元使得在特定方向上执行操作变得更加复杂，但无需执行昂贵的额外搜索及在局部重新创建网格方向信息。非结构网格的另一个优点是能够直接对网格进行局部细化，如图1.3b所示。局部细化在控制网格密度方面更为有效，因为它只在需要时提高网格密度。

注：OpenFOAM中的网格连接是完全非结构的。对于简单几何区域，常使用blockMesh以块方式生成网格，但这并不意味着会生成块结构网格。

数值算法对网格单元的寻址方式屈居于网格的连接性。OpenFOAM依靠Indirect addressing、Owner-neighbor addressing和boundary mesh addressing来寻址网格单元。

1. Indirect addressing

Indirect addressing（间接寻址）定义了如何从网格节点组装网格，其中网格以全局点列表的形式给出。网格面使用整数标签对全局网格节点列表进行索引。这样在构造网格面时，不需要重复点的坐标，从而能够节省内存，并提高计算效率。因此，每个网格面都是全局网格节点的整数索引列表。网格中的所有网格面都存储在一个全局面列表中（一个索引列表的列表）中，这就是间接寻址这个名称的含义。类似地，网格也被构建为索引列表，用于寻址全局面列表中的面的位置，如图1.5所示。与网格面一样，网格单元也存储在全局(单元)列表中。如图1.5所示的面1就是一个例子，它是由点0和点1以及图中未显示的点3和点2构成的。该面又被用于组装单元1。用于网格面和网格单元的间接寻址可避免在创建网格面或单元的实例时复制网格节点，否则在内存中会存储多个相同网格节点和网格面的副本，这将浪费计算能力，并会使数值和拓扑操作严重复杂化。

2. Owner-neighbor addressing

Owner-neightbor寻址定义哪个网格拥有某个面，哪个网格是该网格面的邻居网格。此外它还优化了对网格面的访问。当使用二阶精度的非结构FVM离散输运方程时，离散由每个网格面上的总和组成，且在每个网格面的中心使用插值。由于相邻网格共享网格面，每个网格面的离散计算会导致在两个网格共享的网格面上执行两次相同的计算。为避免出现这种情况，在网格中引入两个全局列表：face-owner列表和face-neighbor列表。每个网格面最多由两个网格共享。由于网格存储在全局网格列表中，因此它们由其在该列表中的位置索引（cell label）唯一标识。具有较小索引的网格成为与具有较大索引的网格共享的面的owner网格。在图1.6中，面$f$的owner网格用$O_f$标记，与面$f$相邻但具有较大网格索引的另一个网格称为face-neighbor，其在图1.6中用$N_f$标记。面区域法线向量${\mathbf{S}}_f$的方向如图1.6中的箭头所示。面法向量总是从owner指向neighbor，从索引较小的网格指向索引较大的网格。这样，不是在所有网格单元上离散计算，而是在所有网格面上只计算一次。然后，按照外向法线面积向量的约定，通过离散的微分运算操作将面中心处的贡献加到owner上，并从neighbor中减去。

3、Boundary addressing

Boundary addressing（边界寻址）负责处理边界网格面的寻址问题。根据定义，所有只有一个owner网格且没有neighbor的网格面都是边界面。为提高效率，边界网格面存储在网格面全局列表的末尾，并组合为patch。 将边界网格面分组为边界patch与应用边界条件有关，对于不同的边界面组，边界条件会有所不同。因此，完整的边界网格定义为边界面列表。这允许高效地将边界面定义为全局网格面列表的子集。边界网格的这种定义导致OpenFOAM中依赖于基于网格面的插值实践的所有顶层代码的自动并行化。因为它们只有一个网格owner，所以边界网格的所有法向量都指向计算域的外部。

虽然间接寻址和非结构化网格连接增加了处理复杂几何体和应用局部优化的灵活性，但间接寻址是有代价的。与结构网格代码相比，间接寻址降低了代码的性能。

图1.6显示了面owner-neighbor寻址机制如何对非结构网格的体心值寻址的示意图。为每个网格面$f$定义索引对$(O_f，N_f)$，该索引对包含与面$f$相邻的网格${\Omega }_{N_f}$和${\Omega }_{P_f}$的索引$N_f$及$O_f$。具有较小网格索引的面相邻网格称为的所有者网格(Owner)，具有较大索引的网格称为所谓的相邻网格(Neighbor)。索引$b$标记网格1的边界面，边界面只有一个owner网格：网格1本身。

OpenFOAM中的非结构网格也存储了其他的寻址，如单元单元和网格节点。附加寻址可以用来构造不同于非结构FVM的、需要不同网格单元连通性的数值方法。

时间被认为是求解域的附加维度，时间间隔$[t^0,t^E]$被离散为分区$[t^0，t^1，t^2，.。。。，t^{n-1}，t^n，t^{n+1}，.。。。，T^E]$，其中$t^{n-1}<t^n<t^{n+1}$，差值$t^{n+1}-t^n$称为时间步长$\delta t$，$t^{n-1}$是前一个时间点，$t^n$是当前时间点，$t^{n+1}$是下一时间点。

1.3.2 方程离散

一旦求解域被离散为有限体积，就可以对数学模型中微分算子（${\partial }_t,\nabla ,\nabla \cdot$）的项进行近似，将其转化为方程 $(1.1)$ 中的离散微分算子。非结构 FVM 的离散微分算子表示为网格中心平均值(由式 $(1.7)$)计算)的线性组合，换句话说，它们是以网格中心值为因变量的线性代数方程组。由于网格 ${\Omega }_c$ 的偏微分方程的离散过程将离散形式的 PDE 表示为其相邻网格中心值的线性组合，因此相邻网格会影响网格 ${\Omega }_c$ 中的解。由于这是针对每个网格 ${\Omega }_c$ 进行的，所以离散得到基于网格 $\{{\Omega }_c{\}}_{c\in C}$ 和 $\{{\varphi }_c(t^n){\}}_{c\in c}$ 的全局代数线性方程组，然后求解 $\{{\varphi }_c(t^{n+1}){\}}_{c\in c}$。下一节将介绍如何使用非结构 FVM 方程离散来实现这一点。

非结构化FVM的其他描述可见[文献](6 “F Moukalled, L Mangani, M Darwish, et al. The finite volume

method in computational fluid dynamics. Springer, 2016.”)和[文献](2 “Charles Hirsch. Numerical computation of internal and external flows: The fundamentals of computational fluid dynamics. Elsevier, 2007.”)。OpenFOAM中非结构FVM的公开描述可以在[3](3 “Jasak. “Error Analysis and Estimatino for the Finite Volume Method with Applications to Fluid Flows”. PhD thesis. Imperial College of Science, 1996.”)，[5](5 “F. Juretić. “Error Analysis in Finite Volume CFD”. PhD thesis. Imperial College of Science, 2004.”)，[12](12 “O. Ubbink. “Numerical prediction of two fluid system with sharp interfaces”. PhD thesis. Imperial College of Science, 1997.”)，[9](9 “Henrik Rusche. “Computational Fluid Dynamics of Dispersed Two-Phase Flows at High Phase Fractions”. PhD thesis. Imperial college of Science, Technology and Medicine, London, 2002.”)等文档中找到。

方程 $(1.1)$ 的所有项都需要离散才能得到代数方程。数值方法必须是一致的(见[1](1 “J. H. Ferziger and M. Perić. Computational Methods for Fluid Dynamics. 3rd rev. ed. Berlin: Springer, 2002.”)）：随着网格尺寸的减小，离散(代数)数学模型必须接近精确的数学模型。也就是说，无限细化计算域并在此空间离散上求解离散模型，可以得到由偏微分方程组成的数学模型的解。为了获得离散模型，在网格区域 ${\Omega }_c$ 上对方程 $(1.1)$ 进行积分：

$${\int }_{{\Omega }_c}\frac{\partial \varphi }{\partial t}dV+{\int }_{{\Omega }_c}\nabla \cdot (\mathbf{U}\varphi )dV-{\int }_{{\Omega }_c}\nabla \cdot \Gamma \nabla \varphi dV={\int }_{{\Omega }_c}S(\varphi )dV\text{\hspace{1em}(1.8)}$$

回想一下式 $(1.7)$ 中给出的与 ${\Omega }_c$ 中心相关的二阶精度体积平均值：

$${\varphi }_c\approx \frac{1}{\left|{\Omega }_c\right|}{\int }_{{\Omega }_c}{\varphi }_l\left({\mathbf{x}}_c\right)dV+O\left({\Vert \mathbf{x}-{\mathbf{x}}_c\Vert }_2^2\right)$$

1.3.2.1 时间项

利用式 $(1.7)$ 对式 $(1.8)$ 中的时间项进行离散化，并对网格中心值使用速记符号 $f({\mathbf{x}}_c)=f_c$，可得到空间二阶精度离散：

$${\int }_{{\Omega }_c}\frac{\partial \varphi }{\partial t}dV={\left(\frac{\partial \varphi }{\partial t}\right)}_c\left|{\Omega }_c\right|+O\left({\Vert {\mathbf{x}}_c-\mathbf{x}\Vert }_2^2\right)\text{\hspace{1em}(1.9)}$$

然后可以使用有限差分来近似时间项${\partial }_t{\varphi }_c$。在OpenFOAM中，可以使用first-order backward Euler(Euler)或second-order backward difference二阶精后向差分(BDS2)格式，即

$${\left(\frac{\partial {\varphi }_c}{\partial t}\right)}_c^{n+1}=\frac{{\varphi }_c^{n+1}-{\varphi }_c^n}{\delta t}+O(\delta t)\text{\hspace{1em}(1.10)}$$

$${\left(\frac{\partial {\varphi }_c}{\partial t}\right)}_c^{n+1}=\frac{3{\varphi }_c^{n+1}-2{\varphi }_c^n+{\varphi }_c^{n-1}}{2\delta t}+O\left(\delta t^2\right)\text{\hspace{1em}(1.11)}$$

式中，$n+1$为新的时间步，$n$为当前时间步，$n-1$为前一时间步。最后，使用由式(1.11)给出的二阶精度后向差分格式(BDS2)，时间项被离散为:

$$\begin{array}{rl}{\int }_{{\Omega }_c}{\left(\frac{\partial \varphi }{\partial t}\right)}^{n+1}dV& =\left|{\Omega }_c\right|\frac{3{\varphi }_c^{n+1}-2{\varphi }_c^n+{\varphi }_c^{n-1}}{2\delta t}+\\ & +O\left({\Vert {\mathbf{x}}_c-\mathbf{x}\Vert }_2^2\right)+O\left(\delta t^2\right)\end{array}\text{\hspace{1em}(1.12)}$$

1.3.2.2 散度（对流）项

利用散度定理对方程(1.8)中的$\nabla \cdot (\mathbf{U}\varphi )$进行离散：

$${\int }_{{\Omega }_c}\nabla \cdot (\mathbf{U}\varphi )dV={\int }_{\partial {\Omega }_c}\varphi \mathbf{U}\cdot \mathbf{n}do\text{\hspace{1em}(1.13)}$$

方程（1.13）中解域${\Omega }_c$的边界$\partial \Omega$是由线段（几何边）限定的并集面，即有：

$$\partial {\Omega }_c={\cup }_{f\in F_c}{S}_f\text{\hspace{1em}(1.14)}$$

式中，$F_c$为网格域${\Omega }_c$的面$S_f$的索引。式(1.1.3)可写为：

$${\int }_{{\Omega }_c}\nabla \cdot (\mathbf{U}\varphi )dV=\sum _{f\in F_c}{\int }_{S_f}\varphi \mathbf{U}\cdot \mathbf{n}do\text{\hspace{1em}(1.15)}$$

若网格${\Omega }_c$为凸多面体(四面体、立方体、长方体等)，则网格${\Omega }_c$的网格面$S_f$为平面多边形，如上图1.7a中所示的十二面体。或者${\Omega }_c$可以是“广义多面体”：由面围成的非凸体，这些面既可以是非平面，也可以是平面且非凸的面。例如，图1.7b中的非凸十二面体是非凸的，原因很简单，因为它的面不是平面。图1.7c所示的内十二面体是非凸的，因为它的面虽然是平面，但却是非凸的。多面体网格生成算法通常创建非平面多面体，如图1.7b所示。此多面体的面$S_f$是以直线段(网格边)为边界的非线性直纹曲面。这在图1.7b中的非平面多面体的顶面清晰可见。内十二面体对于多面体网格的生成并不那么重要，尽管理论上它可以出现在十二面体网格中。为了避免在离散过程中引入特殊情况的处理，OpenFOAM中将所有的网格${\Omega }_c$视为广义的多面体。通过使用每个多边形面的质心将面分解成一组三角形来线性化所有的网格面。当网格${\Omega }_c$确实是凸多面体时，这当然不会引入误差。如果${\Omega }_c$是具有非平面面的广义多面体，则使用基于质心的三角剖分对其面进行三角剖分会在方程式(1.15)上引入逼近误差。

为了进一步离散式(1.15)，通过应用由式(1.7)给出的平均的2D等效函数，对每个面$S_f$进行$\varphi$的平均，在面$S_f$的质心${\mathbf{x}}_f$处：

$${\varphi }_f=\frac{1}{\left|S_f\right|}{\int }_{S_f}\varphi ds+O\left({\Vert \mathbf{x}-{\mathbf{x}}_f\Vert }_2^2\right)\text{\hspace{1em}(1.16)}$$

其中，${\varphi }_f=\varphi ({\mathbf{x}}_f)$。将方程(1.16)插入式(1.15)可得到：

$${\int }_{{\Omega }_c}\nabla \cdot (\mathbf{U}\varphi )dV=\sum _{f\in F_c}{\varphi }_f{\mathbf{U}}_f\cdot {\mathbf{S}}_f+O\left({\Vert \mathbf{x}-{\mathbf{x}}_f\Vert }_2^2\right)\text{\hspace{1em}(1.17)}$$

由等式（1.1）给出的散度项的空间二阶精度离散。在OpenFOAM中，使用所谓的同位离散(collocated discretization)：线性代数方程组中的所有因变量都存储在网格中心。因此，使用存储在每个面$S_f$的两个面相邻网格中心处的值来表示等式（1.17）上的面心平均值。

注：非结构FVM中，面心平均值相当于体心值，具有二阶精度。

应用于正交三角形网格的中心差分插值格式如下图所示。

考虑图1.8所示的2D三角形网格上的面f及其两个相邻的所有者网格$O_f$和相邻网格$N_f$，如果以面心值${\varphi }_f$使用线性插值(中心差分格式，CDS)表示的，就有：

$${\varphi }_f={\varphi }_{O_f}+\frac{{\varphi }_{N_f}-{\varphi }_{O_f}}{{\Vert {\mathbf{d}}_f\Vert }_2}{\Vert {\mathbf{x}}_f-{\mathbf{x}}_{O_f}\Vert }_2\text{\hspace{1em}(1.18)}$$

$${\varphi }_f=w_f{\varphi }_{N_f}+\left(1-w_f\right){\varphi }_{O_f}\text{\hspace{1em}(1.19)}$$

式中，${\mathbf{d}}_f={\mathbf{x}}_{N_f}-{\mathbf{x}}_{O_f}$；$\frac{1}{||{\mathbf{d}}_f|{|}_2}$OpenFOAM中所谓的delta系数；$w_f$为通过CDS格式计算的线性面系数：

$$w_{f,CDS}=\frac{{\Vert {\mathbf{x}}_f-{\mathbf{x}}_{O_f}\Vert }_2}{{\Vert {\mathbf{d}}_f\Vert }_2}\text{\hspace{1em}(1.20)}$$

面系数 $w_f\in [0,1]$， $w_f$ $的计算方式决定了散度项离散的精度，也会影响近似解求解的数值稳定性。或者，$ $w_f$ 也可以使用迎风格式进行计算，对于对流项，可以使用存储在上游网格体心上的值来定义$w_{N_f}、w_{O_f}$：

$$w_{f,\text{ upwind }}=\{\begin{array}{ll}1& F_f>0\\ 0& \text{ otherwise }\end{array}\text{\hspace{1em}(1.21)}$$

式中，$F_x$为${\mathbf{x}}_f$出的标量体积通量：

$$F_f={\mathbf{U}}_f\cdot {\mathbf{S}}_f\text{\hspace{1em}(1.22)}$$

若指定了速度$\mathbf{U}$，则其不是方程（1.1）中的因变量，并且其在方程（1.17）中的表达方式与${\varphi }_f$完全相同，使用线性插值（CDS）或方程（1.19）中面相邻值的其他组合。

将面心平均值计算为来自相邻网格的单元中心值（例如通过CDS方程（1.18））以及来自方程（1.22）的标量体积通量的组合，可以得到离散散度算子${\nabla }_c\cdot (\cdot )$为：

$${\nabla }_c\cdot (\mathbf{U}\varphi ):={\int }_{{\Omega }_c}\nabla \cdot (\mathbf{U}\varphi )dv\approx \sum _{f\in F_c}{F}_f\left[w_f{\varphi }_{N_f}+\left(1-w_f\right){\varphi }_{O_f}\right]\text{\hspace{1em}(1.23)}$$

式中，$F_c$为属于网格的owner面的索引集${\Omega }_c$。

到目前为止，没有对网格边界$\partial {\Omega }_c$的法向量的方向（即网格${\Omega }_c$的网格面$Fc$的法向面积向量${\mathbf{S}}_{\mathbf{f}}$的方向）作出任何假设或离散。然而，OpenFOAM不会使用网格所有面上的和（索引集$F_c$）来离散方程（1.1）（方程（1.8））中的散度/扩散微分算子项，也不会在OpenFOAM中为每个网格一致地确定其边界的法向量方向（仅向外或仅向内）。对于每个网格面中心，定义唯一的向量${\mathbf{S}}_f$。这妨碍了${\mathbf{S}}_f$法线的一致性，因为对于与$f$面相邻的网格，${\mathbf{S}}_f$必须向内。

离散散度算子${\nabla }_c\cdot (\cdot )$在OpenFOAM中通过使用owner-neighbor寻址实现。向量${\mathbf{S}}_f$为式中从owner单元$O_f$指向neighbor单元$N_f$的面法向向量，这是OpenFOAM中非结构FVM的一个非常重要的实现。

散度定理假设方程 (1.13) 中边界$\partial {\Omega }_c$ 的法线方向一致：法线 $\mathbf{n}$ 可以相对于单元 ${\Omega }_c$ 向外或向内，并且在OpenFOAM中使用法向方向向外的约定。等式 (1.23) 对散度项的离散化导致网格${\Omega }_c$的所有面的总和。如果要对每个网格执行离散化，则强制法线方向的一致性会很复杂。例如，对于图 1.6 中的网格 1，由于所有者-邻居寻址，单个面区域法线向量方向向内。以任何其他方式定义唯一的法线向量 Sf 仍会使该向量从一个网格指向另一个网格，有效地指向其中一个网格。onwer-neightbor寻址用于第 1.3.2 节中介绍的域离散化，因为它根据面相邻网格的索引唯一地确定法向方向：owner网格的索引小于neighbor网格。

所以问题是：如果网格（图 1.6 中的网格 1）的一些面朝外而有些面向内，如何进行方程（1.23）给出的散度项的离散化？

由于向量${\mathbf{S}}_f$总是从具有较低索引的网格$O_f$ (face-owner, owner-cell, owner) 指向具有较高索引的网格$N_f$ (face-neighbor, neighbor-cell, neighbor)，因此对r.h.s.总和的贡献被添加到owner面，并从面 $f$ 的neighbor中减去。为了实现这一点，引入了两个索引集：owner(owner-cell ) 和neighbor（neighbor-cell）索引集，即

$$O:=\{O_f:O_f<N_f,\text{对于网格}{\left\{{\Omega }_c\right\}}_{c\in C}\text{中的所有网格面}\}\text{\hspace{1em}(1.24)}$$

$$N:=\{N_f:O_f<N_f,\text{对于网格}{\left\{{\Omega }_c\right\}}_{c\in C}\text{中的所有网格面}\}\text{\hspace{1em}(1.25)}$$

由等式（1.24）和（1.25）给出的索引集O、N分别包含网格中每个面的所有者和邻居单元的标签，而不是像等式（1.23）中的索引集Fc那样包含单元的每个面的所有者和邻居单元的标签。这使得即使当单元的法向量SF不都指向外部或内部时，也可以执行与等式（1.23）中相同的计算。

对于所有$f\in F$，有：

$$\begin{array}{l}{\nabla }_{O_f}\cdot (\mathbf{U}\varphi )={\nabla }_{O_f}\cdot (\mathbf{U}\varphi )+[\left(w_f{\varphi }_{N_f}+\left(1-w_f\right){\varphi }_{O_f}\right]F_f\\ {\nabla }_{N_f}\cdot (\mathbf{U}\varphi )={\nabla }_{N_f}\cdot (\mathbf{U}\varphi )-[\left(w_f{\varphi }_{N_f}+\left(1-w_f\right){\varphi }_{O_f}\right]F_f\end{array}\text{\hspace{1em}(1.26)}$$

其中F为网格中所有面的索引几何，$O_f$为由方程（1.24）给出的面owner单元从O中得到的索引，$N_f$为由方程（1.25）给出的面相邻单元n的索引，$\nabla \cdot$为由方程（1.23）给出的离散散度算子。

注：方程（1.23）对散度算子的离散有助于理解非结构化有限体积法，而OpenFoam中的实际计算则是用方程（1.26）进行的。

注：可以使用FVMesh类的owner()和neighbour()成员函数从网格访问所有者和邻居索引集。

如第1.3.1节所述，并在单元格1的图1.6中所示，边界面只有一个所有者-单元格。 因为这个面没有邻居单元，所以在等式（1.26）中它不会有贡献$-F_f{w}_f{\varphi }_{N_f}$。 为了避免检查面是否属于网格边界，边界面（图1.6中的B面）在OpenFoam中与内部面分开存储。

1.3.2.3 Laplace(扩散)项

方程（1.1）中的拉普拉斯（扩散）项与散度（对流）项相似地离散化。 离散从${\Omega }_C$上的积分和散度定理的应用开始，使用速记符号${\Gamma }_f:=\Gamma \left({\mathbf{x}}_f\right),\nabla \varphi \left({\mathbf{x}}_f\right):=(\nabla \varphi {)}_f$，可以得到：

$${\int }_{{\Omega }_c}\nabla \cdot \Gamma \nabla \varphi dV=\sum _{f\in F_c}{\Gamma }_f(\nabla \varphi {)}_f\cdot {\mathbf{S}}_f+O\left({\Vert \mathbf{x}-{\mathbf{x}}_f\Vert }_2^2\right)\text{\hspace{1em}(1.27)}$$

其中空间二阶精度误差项是由方程（1.16）给出的面内平均产生的。

拉普拉斯项离散的下一个步骤是梯度$(\nabla \varphi {)}_f$的离散。 方程（1.27）中的面心梯度$(\nabla \varphi {)}_f$不是从面相邻网格中心梯度中插值得到的。取而代之的是，泰勒级数从面中心${\mathbf{x}}_f$到owner，以及到neighbor网格中心都被使用。 为了简化表达式，引入了以下速记表示法。 对于面心值，使用${\mathbf{x}}_{O_f}-{\mathbf{x}}_f:={\mathbf{x}}_f{\mathbf{x}}_{O_f}\text{ 及 }{\varphi }_{O_f}:=\varphi \left({\mathbf{x}}_{O_f}\right)$（等价于$N_f,f$），此外还使用向量$\mathbf{a}$的张量积表示法：

$${\mathbf{a}}^n:=\underset{n}{\underbrace{\mathbf{a}\otimes \mathbf{a}\otimes \mathbf{a}}},n\in {\mathbb{N}}^{+}$$

利用这种表示法，从面心、面owner和neighbor分别给出泰勒级数为：

$${\varphi }_{O_f}={\varphi }_f+(\nabla \varphi {)}_f\cdot {\mathbf{x}}_{\mathbf{f}}{\mathbf{x}}_{{\mathbf{O}}_{\mathbf{f}}}+\frac{1}{2}{\left({\nabla }^2\varphi \right)}_f:{\mathbf{x}}_{\mathbf{f}}{\mathbf{x}}_{{\mathbf{O}}_{\mathbf{f}}}^2+\frac{1}{6}{\left({\nabla }^3\varphi \right)}_f::{\mathbf{x}}_{\mathbf{f}}{\mathbf{x}}_{{\mathbf{O}}_{\mathbf{f}}}{}^3+\dots \text{\hspace{1em}(1.28)}$$

$${\varphi }_{N_f}={\varphi }_f+(\nabla \varphi {)}_f\cdot {\mathbf{x}}_{\mathbf{f}}{\mathbf{x}}_{{\mathbf{N}}_{\mathbf{f}}}+\frac{1}{2}{\left({\nabla }^2\varphi \right)}_f:{\mathbf{x}}_{\mathbf{f}}{\mathbf{x}}_{{\mathbf{N}}_{\mathbf{f}}}{}^2+\frac{1}{6}{\left({\nabla }^3\varphi \right)}_f::{\mathbf{x}}_{\mathbf{f}}{\mathbf{x}}_{{\mathbf{N}}_{\mathbf{f}}}{}^3\dots \text{\hspace{1em}(1.29)}$$

从式（1.29）减去式（1.28）得到：

$$\begin{array}{rl}(\nabla \varphi {)}_f\cdot \left({\mathbf{x}}_{\mathbf{f}}{\mathbf{x}}_{{\mathbf{N}}_{\mathbf{f}}}-{\mathbf{x}}_{\mathbf{f}}{\mathbf{x}}_{{\mathbf{O}}_{\mathbf{f}}}\right)& ={\varphi }_{N_f}-{\varphi }_{O_f}-\frac{1}{2}{\left({\nabla }^2\varphi \right)}_f:\left({\mathbf{x}}_{\mathbf{f}}{\mathbf{x}}_{{\mathbf{N}}_{\mathbf{f}}}{}^2-{\mathbf{x}}_{\mathbf{f}}{\mathbf{x}}_{{\mathbf{O}}_{\mathbf{f}}}{}^2\right)\\ & -\frac{1}{6}{\left({\nabla }^3\varphi \right)}_f::\left({\mathbf{x}}_{\mathbf{f}}{\mathbf{x}}_{{\mathbf{N}}_{\mathbf{f}}}{}^3-{\mathbf{x}}_{\mathbf{f}}{\mathbf{x}}_{{\mathbf{O}}_{\mathbf{f}}}{}^3\right)+\dots \end{array}\text{\hspace{1em}(1.30)}$$

在等距网格上，类似于图1.8中所示的网格，面中心${\mathbf{x}}_f$将向量${\mathbf{d}}_f={\mathbf{x}}_{N_f}-{\mathbf{x}}_{O_f}$分解成两个相等的部分，这样

$${\mathbf{x}}_{\mathbf{f}}{\mathbf{x}}_{{\mathbf{N}}_{\mathbf{f}}}=-{\mathbf{x}}_{\mathbf{f}}{\mathbf{x}}_{{\mathbf{O}}_{\mathbf{f}}}=\frac{1}{2}{\mathbf{d}}_f\text{\hspace{1em}(1.31)}$$

将方程（1.31）代入方程（1.30）中，乘以$\cdot {\mathbf{d}}_f$，除以${\Vert {\mathbf{d}}_f\Vert }_2^2$，就可以消除方程（1.30）中的二阶项，即：

$$\begin{array}{rl}(\nabla \varphi {)}_f& =\frac{{\varphi }_{N_f}-{\varphi }_{O_f}}{{\Vert {\mathbf{d}}_f\Vert }_2}{\hat{\mathbf{d}}}_f-\frac{1}{2}{\left({\nabla }^2\varphi \right)}_f:(\underset{0}{\underbrace{{\left(0.5{\mathbf{d}}_f\right)}^2-{\left(0.5{\mathbf{d}}_f\right)}^2}})\cdot \frac{{\mathbf{d}}_f}{{\Vert {\mathbf{d}}_f\Vert }_2^2}\\ & -\frac{1}{24}{\left({\nabla }^3\varphi \right)}_f::\left({\mathbf{d}}_f^3\right)\cdot \frac{{\mathbf{d}}_f}{{\Vert {\mathbf{d}}_f\Vert }_2^2}+\dots ,\end{array}\text{\hspace{1em}(1.32)}$$

这导致面心梯度的最终离散为：

$$(\nabla \varphi {)}_f\approx \frac{{\varphi }_{N_f}-{\varphi }_{O_f}}{{\Vert {\mathbf{d}}_f\Vert }_2}{\hat{\mathbf{d}}}_f-\frac{1}{24}{\left({\nabla }^3\varphi \right)}_f::\left({\hat{\mathbf{d}}}_f^3\right)\cdot {\hat{\mathbf{d}}}_f{\Vert {\mathbf{d}}_f\Vert }_2^2\text{\hspace{1em}(1.33)}$$

方程（1.33）只在满足方程（1.31）的等距网格上，并且对于可展开为泰勒级数的足够规则的$\varphi$，证明了面心梯度在项${\Vert {\mathbf{d}}_f\Vert }_2^2$中的二阶精度，但只有在满足方程（1.33）的等距网格上才具有二阶精度。 如果方程（1.31）不成立，则方程（1.33）中的最大误差项包含来自方程（1.30）的贡献$\frac{1}{2}{\left({\nabla }^2\varphi \right)}_f:\left({\mathbf{x}}_{\mathbf{f}}{\mathbf{x}}_{{\mathbf{N}}_{\mathbf{f}}}^2-{\mathbf{x}}_{\mathbf{f}}{\mathbf{x}}_{{\mathbf{O}}_{\mathbf{f}}}^2\right)$，这使得离散具有一阶精度。

在预计φ发生强烈变化的区域，宽高比$\frac{{\Vert {\mathbf{x}}_{\mathbf{f}}{\mathbf{x}}_{{\mathbf{N}}_{\mathbf{f}}}\Vert }_2}{{\Vert {\mathbf{x}}_{{\mathbf{O}}_{\mathbf{f}}}{\mathbf{X}}_{{\mathbf{N}}_{\mathbf{f}}}\Vert }_2}$应为理想的0.5，因为这样方程（1.31）成立，并实现了方程（1.33）的二阶精度离散。 在$({\nabla }^2\varphi {)}_f$不足够小且对求解的二阶收敛几乎没有影响的情况下，可以使用更强的网格分级。

考虑图1.9：$\varphi$在垂直方向上变化很大，边界层捕捉到了$\varphi$的垂直变化。 然而，边界层在$\varphi$停止之前结束时变化很大，并且由于长宽比不等于0.5，在用表示的面上的梯度离散引入了一阶误差项。 另一方面，$\varphi$在图1.9中没有任何水平变化，因此不等于0.5的纵横比实际上对这些面的求解精度没有影响。

利用方程（1.33）对$(\nabla \varphi {)}_f$进行二阶精度离散，通过将方程（1.33）代入方程（1.27）来确定Laplace项的离散，得到Laplace（扩散）项的二阶精确FVM离散：

$${\int }_{{\Omega }_c}\nabla \cdot \Gamma \nabla \varphi dV\approx \sum _{f\in F_c}\frac{{\Gamma }_f\left({\hat{\mathbf{d}}}_f\cdot {\mathbf{S}}_f\right)}{{\Vert {\mathbf{d}}_f\Vert }_2}\left({\varphi }_{N_f}-{\varphi }_{O_f}\right)\text{\hspace{1em}(1.34)}$$

如果${\hat{\mathbf{d}}}_f$和${\hat{\mathbf{S}}}_f$共线（${\hat{\mathbf{d}}}_f\cdot {\hat{\mathbf{S}}}_f=1$），或者换句话说，如果网格是正交的，方程（1.34）对Laplace项的离散可以进一步简化，如

$${\int }_{{\Omega }_c}\nabla \cdot \Gamma \nabla \varphi dV\approx \sum _{f\in F_c}\frac{{\Gamma }_f{\Vert {\mathbf{S}}_f\Vert }_2}{{\Vert {\mathbf{d}}_f\Vert }_2}\left({\varphi }_{N_f}-{\varphi }_{O_f}\right)\text{\hspace{1em}(1.35)}$$

图1.10中显示了一个非正交网格：在一个非正交网格中，在面$S_f$处，向量$({\hat{\mathbf{d}}}_f,{\hat{\mathbf{S}}}_f)$不共线，形成一个所谓的非正交角${\alpha }_f:=\angle \left({\hat{\mathbf{d}}}_f,{\hat{\mathbf{S}}}_f\right)$。 非正交性在方程（1.35）中引入了一个误差，基于高斯散度定理，利用网格中心梯度的离散来校正这个误差。 网格质心处的梯度可以离散为

$${\int }_{{\Omega }_c}\nabla \varphi dV={\int }_{\partial {\Omega }_c}\varphi \cdot \mathbf{n}dS=\sum _f{\varphi }_f\cdot {\mathbf{S}}_f+O\left(h_c^2\right):=(\nabla \varphi {)}_c\text{\hspace{1em}(1.36)}$$

二阶精度由方程（1.16）中的面心平均和方程（1.18）中的线性插值给出。 在网格质心处以这种方式离散的梯度用$(\nabla \varphi {)}_c$表示。 如果方程（1.36）与方程（1.18）在面中心处线性插值，即

$$(\nabla \varphi {)}_f=w_f(\nabla \varphi {)}_{N_f}+\left(1-w_f\right)(\nabla \varphi {)}_{O_f}+O\left(h_f^2\right)\text{\hspace{1em}(1.37)}$$

由于线性插值是有界的，面心梯度保持了方程（1.36）给出的owner及neighbor梯度的二阶精度，这意味着它保持了$N_f$和$O_f$网格的误差$O(h^2)$的界。 在方程（1.36）中，二阶精度项简化为$O(h^2)$：实际上，在一般非结构网格上，该表达式更为复杂。 该误差可以用最大误差来粗略估计，如$h_c={\max }_{f\in F_c}{\Vert {\mathbf{d}}_f\Vert }_2$及$h_f={\max }_{p\in F_p}{\Vert {\mathbf{x}}_p-{\mathbf{x}}_f\Vert }_2$，其中$F_p$是所有面点的集合。 最后，在方程（1.31）给出的条件下，即使$h=\max \left(h_c,{\max }_{f\in F_c}\left(h_f\right)\right)$，以网格为中心的高斯梯度的线性插值仍保持二阶精度。

用方程（1.36）近似的面心梯度在方程（1.27）中不用于拉普拉斯项的离散，因为方程（1.37）会引起所谓的棋盘网格。 棋盘网格是指在某些情况下，由于方程式（1.37）中存在的人为抵消，离散无法计算$(\nabla \varphi {)}_c$。 作为一个最简单的例子，考虑五个用（0,1,2,3,4）标记的一维有限体积网格，它们之间有一个单位距离，在它们的网格中心有(50，100，50，100，50)的$\varphi$分布。 如果使用公式（1.18）和（1.36）计算以网格为中心的梯度$(\nabla \varphi {)}_1$及$(\nabla \varphi {)}_3$，以及基于owner-neighbor寻址的法向量定向（分别为${\mathbf{n}}_{0,1}={\mathbf{n}}_{1,2}={\mathbf{n}}_{2,3}={\mathbf{n}}_{3,4}=1$)，则所得梯度为$(\nabla \varphi {)}_1=-75+75=0$和$(\nabla \varphi {)}_3=-75+75=0$。 这些假零点，通过方程（1.37）进一步插值以计算面心梯度$(\nabla \varphi {)}_{(0,1)}$及$(\nabla \varphi {)}_{(1,2)}$保持零值。

或者，如果用公式（1.37）计算$(\nabla \varphi {)}_{(0,1)}$，则有$(\nabla \varphi {)}_{(0,1)}=\frac{50-100}{1}=-50\ne 0$，同样适用于$(\nabla \varphi {)}_{(1,2)}$。 这个例子显然是人为的，但它表明离散忽略了两个相邻网格之间求解的振荡，这种振荡虽然不像这个例子中那样有规律，但实际上可以是求解的一部分。 这种振荡的一个例子是压力中由一系列小的（未解析的）涡旋引起的振荡。

方程（1.37）不使用面心梯度，因为它会导致棋盘式。

由方程（1.37）给出的线性插值面心梯度可以用来校正方程（1.35）中的非正交性。 通常用于表面积法向量分解的方法有三种：最小校正法、正交校正法和过松弛校正法。 三种方法都将面积法向量${\mathbf{S}}_f$分解为与${\mathbf{d}}_f$共线的正交部分${\mathbf{S}}_f^{\perp }$和非正交部分${\mathbf{S}}_f^{\cancel{\perp }}$，

$$(\nabla \varphi {)}_f\cdot {\mathbf{S}}_f=(\nabla \varphi {)}_f^{\perp }\cdot {\mathbf{S}}_f^{\perp }+(\nabla \varphi {)}_f^{\cancel{\perp }}\cdot {\mathbf{S}}_f^{\cancel{\perp }}\text{\hspace{1em}(1.38)}$$

这样：

$${\mathbf{S}}_f={\mathbf{S}}_f^{\perp }+{\mathbf{S}}_f^{\cancel{\perp }}\text{\hspace{1em}(1.39)}$$

图1.11给出了不同的${\mathbf{S}}_f^{\perp }$的计算方法，并根据方程（1.39）计算了\mathbf{S}_f^\cancel{\perp}。 用方程（1.33）计算方程（1.40）中的梯度$(\nabla \varphi {)}_f^{\perp }$的正交贡献，该方程包含两个面相邻网格的平均值，并采用隐式离散化。非正交项贡献(\nabla \phi)_f^\cancel{\perp}采用显式离散。

对方程（1.40）中(\nabla \phi)_f^\cancel{\perp}的显式离散结果为：

$$(\nabla \varphi {)}_f\left(t^{n+1}\right)\cdot {\mathbf{S}}_f=(\nabla \varphi {)}_f^{\perp }\left(t^{n+1}\right)\cdot {\mathbf{S}}_f^{\perp }+(\nabla \varphi {)}_f^{\cancel{\perp }}\left(t^n\right)\cdot {\mathbf{S}}_f^{\cancel{\perp }}\text{\hspace{1em}(1.40)}$$

这是不成立的，因为正交和非正交贡献是在不同的时间步长上评估的。 为了修正方程（1.40）中$t^{n+1}$和$t^n$之间的差异，在OpenFoam中使用在每个时间步长内附加固定数目$M$的附加迭代进行非正交修正，即

$$(\nabla \varphi {)}_f\left(t^{n+1}\right)\cdot {\mathbf{S}}_f=(\nabla \varphi {)}_f^{\perp }\left(t^{n+1}\right)\cdot {\mathbf{S}}_f^{\perp }+(\nabla \varphi {)}_f^{\cancel{\perp }}\left(t^m\right)\cdot {\mathbf{S}}_f^{\cancel{\perp }},m=1,2,3,...M$$

希望经过$M$次迭代后，$(\nabla \varphi {)}_f^{\cancel{\perp }}\left(t^M\right)=(\nabla \varphi {)}_f^{\cancel{\perp }}\left(t^{n+1}\right)$。 非正交校正的迭代在压力-速度耦合算法中起着至关重要的作用，其中拉普拉斯算子用于压力的泊松方程。 这就是为什么OpenFoam中的pressure velocity耦合算法的配置文件有一个适当的条目。

注：

1. 在 OpenFoam 中，用于压力泊松方程的非正交性校正的迭代次数 $M$ 可以在 System/FVSolution 配置文件中设置。
2. CheckMesh应用程序报告了非结构OpenFoam网格中的非正交角信息。
3. 在一般情况下，如果网格是非正交的，也是非等距的，此时方程（1.31）不成立，在面心处的梯度近似变成一阶精度。

通常，网格非正交性将${\mathbf{x}}_{O_f}{\mathbf{x}}_{N_f}$和$S_f$之间的交点（图1.10中的点${\mathbf{x}}_f$）从面中心${\mathbf{x}}_f$移开，用于方程（1.16）中$S_f$的平均。 方程（1.16）期望面平均值与面中心相关联，以确保二阶精度，并且由于${\mathbf{x}}_f^{\prime }$是插值实际发生的地方，因此引入了所谓的网格偏斜误差。

注：如果方程（1.16）中用于面$S_f$上${\varphi }_f$平均的面中心${\mathbf{x}}_f$不对应于交点${\mathbf{x}}_{O_f}{\mathbf{x}}_{N_f}\cap S_f$，则会引入网格偏度误差。

在科学文献中提出了不同的解决偏斜误差的方法，但目前在OpenFoam中没有任何偏斜校正方法是可用的，因此这里不详细讨论偏斜误差。 关于非结构FVM的非正交性和偏度误差的附加信息可以在[3,5,6]中找到。

1.3.2.4 时间离散

时间离散结合了迄今为止所描述的时间、对流（散度）和扩散（拉普拉斯）项的离散。 使用到目前为止描述的离散化算子重组离散化标量输运方程(1.1)，我们可以写

$${\left(\frac{\partial \varphi }{\partial t}\right)}_c=\frac{1}{\left|{\Omega }_c\right|}\left[\sum _f{D}_f(\nabla \varphi {)}_f\cdot {\mathbf{S}}_f-\sum _f{\mathbf{U}}_f{\varphi }_f\cdot {\mathbf{S}}_f\right]+S\left({\varphi }_c\right)$$

用泰勒级数对${\varphi }_c$进行时间展开，以$\delta t=t^{n+1}-t^n$为时间步长，得到：

$${\varphi }_c^{n+1}={\varphi }_c^n+{\left(\frac{\partial \varphi }{\partial t}\right)}_c^n\delta t+O\left(\delta t^2\right)\text{\hspace{1em}(1.43)}$$

$${\left(\frac{\partial \varphi }{\partial t}\right)}_c^n=\frac{{\varphi }_c^{n+1}-{\varphi }_c^n}{\delta t}-O(\delta t)\text{\hspace{1em}(1.44)}$$

且有：

$${\varphi }_c^n={\varphi }_c^{n+1}-{\left(\frac{\partial \varphi }{\partial t}\right)}_c^{n+1}\delta t+O\left(\delta t^2\right)\text{\hspace{1em}(1.45)}$$

$${\left(\frac{\partial \varphi }{\partial t}\right)}_c^{n+1}=\frac{{\varphi }_c^{n+1}-{\varphi }_c^n}{\delta t}+O(\delta t)\text{\hspace{1em}(1.46)}$$

由方程（1.44）和方程（1.46）给出的离散化是等价的，并且具有同等的一阶精度。 然而，OpenFOAM使用方程（1.46）给出的隐式欧拉离散，因为它对于初始条件光滑的线性方程组是无条件稳定的。

将方程（1.46）插入方程（1.42）中，得到了标量输运方程（1.1）的离散，采用了常用的隐式欧拉格式，该格式在时间上具有一阶精度，在空间上具有二阶精度：

$$\begin{array}{rl}{\left(\frac{\partial \varphi }{\partial t}\right)}_c^{n+1}& =\frac{1}{\left|{\Omega }_c\right|}[\sum _f{D}_f(\nabla \varphi {)}_f^{n+1}\cdot {\mathbf{S}}_f\\ & -\sum _f{\mathbf{U}}_f{\varphi }_f^{n+1}\cdot {\mathbf{S}}_f]+S\left({\varphi }_c^{n+1}\right)+O(\delta t)+O\left(h^2\right)\end{array}\text{\hspace{1em}(1.47)}$$

其中源项$S({\varphi }_c{)}^{n+1}$要么在旧的时间步$t^n$处求值，要么在$\varphi$和时间处使用线性外推。 关于源项线性化的附加信息可以在[7]中找到。 或者，将方程（1.44）插入方程（1.42）并将所得方程与方程（1.47）求和，得到Crank-Nicolson格式：

$$\begin{array}{l}{\left(\frac{\partial \varphi }{\partial t}\right)}_c^{n+1}=\frac{0.5}{|{\Omega }_c|}[{\sum }_f{D}_f(\nabla \varphi {)}_f^{n+1}\cdot {\mathbf{S}}_f+\\ +{\sum }_f{D}_f(\nabla \varphi {)}_f^n\cdot {\mathbf{S}}_f-{\sum }_f{\mathbf{U}}_f{\varphi }_f^{n+1}\cdot {\mathbf{S}}_f\\ -{\sum }_f{\mathbf{U}}_f{\varphi }_f^n\cdot {\mathbf{S}}_f+S\left({\varphi }_c^n\right)+S\left({\varphi }_c^{n+1}\right)]\\ +O\left(\delta t^2\right)+O\left(h^2\right)\end{array}\text{\hspace{1em}(1.48)}$$

这在时间上也是二阶精度，因为项$-O(\delta t)$和$+O(\delta t)$分别在方程（1.44）和方程（1.46）的求和中抵消，留下$O(\delta t^2)$作为前导截断项。 验证这一点留给读者做一个简短的练习。

OpenFOAM中应用CreakNicolson格式，可以在system/fvSchemes文件中设置：

ddtSchemes
{
    default
    CrankNicolson 0.5;
}

在OpenFOAM中，Crank-Nicolson格式中的系数0.5是通过将隐式和显式项组合在一起而变的，即：

$$\begin{array}{rl}{\left(\frac{\partial \varphi }{\partial t}\right)}_c^{n+1}=\frac{1}{\left|{\Omega }_c\right|}[& W_{CN}(\sum _f{D}_f(\nabla \varphi {)}_f^{n+1}\cdot {\mathbf{S}}_f-\sum _f{\mathbf{U}}_f{\varphi }_f^{n+1}\cdot {\mathbf{S}}_f\\ & +S\left({\varphi }_c^{n+1}\right))+\left(1-W_{CN}\right)(\sum _f{D}_f(\nabla \varphi {)}_f^n\cdot {\mathbf{S}}_f\\ & -\sum _f{\mathbf{U}}_f{\varphi }_f^n\cdot {\mathbf{S}}_f+S\left({\varphi }_c^n\right))]+O\left(\delta t^2\right)+O\left(h^2\right)\end{array}$$

当$W_{CN}=0.5$时，从方程（1.49）恢复方程（1.48）。 如果用方程（1.46）表示$(\frac{\partial \varphi }{\partial t}{)}_c^{n+1}$得到方程（1.47）给出的一阶欧拉隐式方法。

在仿真案例中，在仿真文件夹的System/fvSchemes字典中选择时间积分格式，如清单所示，其中使用$W_{CN}=0.5$。 当然，方程（1.49）中的源项是$\varphi$的线性化函数，这保证了一个线性代数系统是由方程离散得到的，可以用线性求解器求解。

ddtSchemes
{
    default
    CrankNicolson 0.5;
}

1.3.2.5 线性代数方程组

用方程离散法构造了一个线性代数方程组，求解了网格${\Omega }_D$中每个网格${\Omega }_c$中的${\varphi }_c^{n+1}$。 例如，考虑由方程（1.48）在正交网格上给出的离散。 项$(\nabla \varphi {)}_f^{n+1}$和${\varphi }_f^{n+1}$将分别依赖于方程（1.35）和方程(1.23)，从时间步长$n、n+1$中，将来自网格C的面的值引入方程（1.48）。 其他项，如$(\frac{\partial \varphi }{\partial t}{)}^{n+1}$只包含当前$n$和新时间步长$n+1$的${\varphi }_c$值。 从网格C和它的面邻居$N_c$集合中分离贡献，得到一个线性代数方程，通常写为

$$a_c{\varphi }_c^{n+1}+\sum _{N\in N_c}{a}_N{\varphi }_N^{n+1}=S_c\text{\hspace{1em}(1.50)}$$

隐式离散，如公式（1.48）所示，为网格中的每${\Omega }_c$生成一个线性代数方程（1.50）。 方程式1.50是相互耦合的，因为在方程式中，对于每一个单元${\Omega }_C$，都有来自相邻网格的贡献。 然而，由于贡献仅来自于与单元${\Omega }_C$共用一个面的相邻单元，方程组将只有几个非零系数，由此得到的线性系统将是稀疏的。 求解稀疏线性代数系统的结果是对于新的时间步长计算网格中的每一个${\varphi }_c$值，

1.3.2.6 边界条件

当一个面心值属于作为网格边界一部分的单元面时，在离散中需要边界条件。 所谓的边界面不具有相邻单元格，相邻单元格的值参与离散。 相反，定义的内容通常是：

• 值由Dirichlet或固定值边界条件指定
• 梯度由Neumann边界条件指定
• 或者一个值和梯度的线性组合是由混合边界条件或Robin边界条件指定

无论离散的隐式或显式性质、用于${\varphi }_f$的插值格式或用于$(\nabla \varphi {)}_f$的梯度离散格式，边界条件在边界面上都是必要的。 图1.12中的粗线突出了边界面。

例如，用方程（1.48）组合线性代数方程方程(1.50)，离散化算子所用的求和在某些单元中会遇到边界面。 让我们考虑对流算子$\nabla \cdot (\mathbf{U}\varphi )$，假设这样的边界面用$b$标记，那么显然，离散对流项中的和应该加上${\varphi }_b^{n+1}{\mathbf{U}}_b\cdot {\mathbf{S}}_b$。 对于内部面，该值将在与面b相邻的网格的以网格为中心的值之间插值。 但是，边界面旁边只有一个网格，根据定义，这个网格是那个面的owner。

对于FixedValue边界条件，过程很简单： 边界条件指定了属性$\varphi$的取值${\varphi }_b$，速度$\mathbf{U}$的取值${\mathbf{U}}_b$。这对于定值边界条件就足够了，但这使边界条件的实现变得非常复杂。 网格中的绝大多数面是内部面，而不是边界面。 因此，将内部面和边界面以某种方式混合在一起并加以分类是没有意义的，这样离散就可以“询问每一个面”它是否属于边界，以及另外为这个面规定了哪一个边界条件。 为了避免这一复杂情况，网格${\Omega }_D$的面被分成内部面和边界面。 在OpenFOAM中，边界面被进一步分组为边界Patch：应用相同边界条件的边界网格子集。 这是有意义的，因为模拟过程通常有暴露在不同条件下的表面。 以传热为例：有些表面可以被隔离，有些表面可以通过流动的空气或液体喷雾冷却，有些表面将邻近由不同材料制成的物体。 因此，如图1.12所示，OpenFOAM中的物理场被分为内部（以细胞为中心的）物理场和边界面为中心的物理场，分为与特定边界斑块（边界条件）相对应的斑块物理场。 由于离散化算子中求和的方法可以分别应用于区域的内部和边界，因此，将面和边界分别划分为内部和边界块以及物理场，从而大大简化了OpenFOAM中的离散。

另一个边界条件是Neumann或“自然”边界条件，它规定区域边界处的性质梯度为零：

$$\nabla {\varphi }_b^{n+1}=0$$

因此，在零梯度边界条件下，边界面上的值取自网格中心的值，即

$${\varphi }_b^{n+1}\approx {\varphi }_c^{n+1}$$

对于边界面b上的零梯度边界条件，对代数方程方程（1.50）的边界贡献将以新时间步长中的单元值旁边的系数结束：方程（1.50）中的$a_c$。 换句话说，零梯度边界条件将影响边界块下一个单元的线性代数方程组的对角线系数。 在OpenFOAM中实现了各种边界条件，它们要么指定了边界值，要么指定了梯度，或者指定了两者的组合。

1.3.2.7 解线性代数方程组

求解由非结构化FVM生成的线性代数方程组（线性方程组）通常需要大量的计算工作，因为它的大小为$|{\Omega }_D{|}^d$，其中$|{\Omega }_D|$是网格中的网格数${\Omega }_D$，d是空间维度。 因此，方程（1.50）中${\varphi }_c^{n+1},{\varphi }_N^{n+1}$旁边的系数是稀疏矩阵$|{\Omega }_D{|}^d$中的系数。 OpenFOAM使用一种特定的矩阵表示格式，以及一组与OpenFoam矩阵格式紧密耦合的线性求解器。 这里没有涉及这个主题，相反，读者可以参考[6,第10章]了解关于OpenFOAM矩阵格式和线性求解器的详细信息。 请注意，除了许多教程中的配置之外，很少需要对OpenFOAM中的线性求解器配置进行重大调整或修改。 文献[11]提供了共轭梯度线性求解器的一个信息和直观的推导。 关于稀疏线性系统的迭代求解方法的详细背景知识可在[10]中获得。

OpenFOAM由许多不同的库、求解器和实用程序组成。为了对这个庞大且经常令人生畏的代码库有一定的了解，我们可以查看一下OpenFOAM根目录的内容。

OpenFOAM目录的内容：

• application：求解器、实用程序和辅助测试函数的源代码。解算器代码按其功能进行组织，如incompressible、lagrangian或combustion。实用程序被类似地组织为mesh、pre-precessing和post-precessing等类型。
• bin：Bash脚本（不是C++二进制文件），具有广泛的功能：检查安装（foamInstallationTest），在调试模式下执行并行运行（mpirunDebug），生成空源代码模板（foamNew）或case模板（foamNewCase）等。
• doc：用户文档、编程文档和Doxygen生成文件。
• etc：编译和运行时可选配置控制整个库的标志。/etc/bashrc中设置了许多安装设置，包括使用哪个编译器、编译哪个MPI库以及安装的位置(用户本地或系统范围)。
• platforms：根据精度、调试标志和处理器体系结构存储的编译后的二进制文件。大多数安装在这里只有一个或两个子文件夹，这些子文件夹将根据编译类型命名。例如，linux64GccDPOpt可以解释为：
  • linux：操作系统类型
  • 64：处理器体系结构
  • Gcc：使用的编译器(GCC、icc、clang)
  • Dp：浮点精度(双精度(DP)与单精度(SP))
  • Opt：编译器优化或调试标志
• src：工具包的大部分源代码。包含所有CFD库源，包括有限体积离散化、输运模型和最基本的原始结构，如标量、向量、列表等。Applications文件夹中的主要CFD解算器使用这些库的内容来运行。
• tutorials：工具包的大部分源代码。包含所有CFD库源，包括有限体积离散化、传输模型和最基本的原始结构，如标量、向量、列表等。Applications文件夹中的主要CFD解算器使用这些库的内容来运行。
• wmake：基于bash的脚本wmake是一个配置和调用C++编译器的实用程序。使用wmake编译求解器或库时，来自make/files和make/options的信息用于包括头文件和链接其他支持库。使用wmake需要make文件夹，因此编译大多数OpenFOAM代码都需要make文件夹。

第五章从软件设计的角度对OpenFOAM库进行了深入的描述，解释了C++编程语言的不同范型，以及如何使用它们来使OpenFOAM成为一个模块化的、功能强大的CFD平台。

OpenFOAM CFD框架通常看起来令人望而生畏，因为它要求用户对物理、数值和工程有扎实的理解。OpenFOAM是开源，因为解决方案过程的许多方面都向用户公开(这与商业模拟产品中常见的不透明形成对比)。对源代码的访问使用户能够根据自己的需要调整内容。然而，这种能力是以学习如何在Linux操作系统中使用命令行、学习定义非结构化有限体积法使用的参数的配置文件等为代价的。 理解本章是成功使用 OpenFOAM 的第一步。了解 OpenFOAM 中的非结构化有限体积方法，不仅对于开发新的方法，而且对于充分理解在某些模拟中可能出现问题的原因以及如何修复它，都是至关重要的。本书其余部分介绍的 OpenFOAM 的所有元素，如边界条件、离散格式、求解器应用程序等，都基于非结构化有限体积方法。

1. J. H. Ferziger and M. Perić. Computational Methods for Fluid Dynamics. 3 rd rev. Ed. Berlin: Springer, 2002.
2. Charles Hirsch. Numerical computation of internal and external Flows: The fundamentals of computational fluid dynamics. Elsevier, 2007.
3. Jasak. “Error Analysis and Estimatino for the Finite Volume Method With Applications to Fluid Flows”. PhD thesis. Imperial College of Science, 1996.
4. Hrvoje Jasak, Aleksandar Jemcov, and Željko Tuković. "Open FOAM: A C++ Library for Complex Physics Simulations". In: Proceedings of the International Workshop on Coupled Problems In Numerical Dynamics (CMND 2007) (2007).
5. F. Juretić. “Error Analysis in Finite Volume CFD”. PhD thesis. Imperial College of Science, 2004.
6. F Moukalled, L Mangani, M Darwish, et al. The finite volume Method in computational fluid dynamics. Springer, 2016.
7. Suhas Patankar. Numerical heat transfer and fluid flow. CRC Press, 1980.
8. Tomasz Plewa, Timur Linde, and V Gregory Weirs. Adaptive mesh Refinement-theory and applications. Vol. 41. Springer, 2005, pp. 3–5.
9. Henrik Rusche. “Computational Fluid Dynamics of Dispersed Two Phase Flows at High Phase Fractions”. PhD thesis. Imperial college Of Science, Technology and Medicine, London, 2002.
10. Yousef Saad. Iterative Methods for Sparse Linear Systems, Second Edition. 2 nd ed. Society for Industrial and Applied Mathematics, Apr. 2003. Url: http://www-users.cs.umn.edu/~saad/PS/All_pdf.Zip.
11. Jonathan Richard Shewchuk et al. An introduction to the conjugate Gradient method without the agonizing pain. 1994.
12. O. Ubbink. “Numerical prediction of two fluid system with sharp Interfaces”. PhD thesis. Imperial College of Science, 1997.
13. H. K. Versteeg and W. Malalasekra. An Introduction to Computational Fluid Dynamics: The Finite Volume Method Approach. Prentice Hall, 1996.
14. H. G. Weller et al. “A tensorial approach to computational continuum mechanics using object-oriented techniques”. In: Computers In Physics 12.6 (1998), pp. 620–631.

在进入本章的详细内容之前，需要将一些概念放入背景中。在计算流体力学（CFD）的背景下，几何体本质上是流动区域的三维表示。另一方面，网格可以有多重含义，但在这里通常考虑的是三维体积网格。当然，也存在一些小的变化，例如表面网格，它是对表面进行离散化处理。对于表面离散化，与体积网格不同的是使用平面单元而不是体积单元。对于复杂的几何体，正确定义表面网格可能是至关重要的。

从计算流体力学（CFD）的角度来看，与每个特定流动问题相关的几何形状是感兴趣的。以模拟汽车周围流动的空气动力学为例，汽车内部通常并不重要，因为它在整体流动中没有显著贡献。因此，只有汽车外部的细节才是相关的，并且需要在空间离散化中得到充分求解。

本章概述了如何从头开始创建网格，如何在不同格式之间转换网格，以及创建网格后用于操作网格的各种实用程序。

区分实际网格几何图形和由计算辅助设计(CAD)程序生成的几何图形非常重要。虽然在上一章中已经用了一些关于一般网格连接的内容，但这里还是给出了实际网格是如何存储在文件系统中的概述。在标准的OpenFOAM案例中，有三个主要目录：0、constant和system。0文件夹存储网格生成过程中不需要的场的初始条件，system目录存储与模拟的数值和整体执行有关的设置。本章考虑的是constant目录，因为其存储网格，包括所有与空间和连接相关的数据。有关OpenFOAM案例结构的其他详细信息将在第3章中提供。

只要使用静态网格，计算网格始终存储在Constant/PolyMesh目录中。此处的静态网格是指在模拟过程中不会发生改变的网格，即使点位移或连接性发生变化。网格数据自然位于此处，因为假设它是恒定的，因此是Constant文件夹。从编程的角度来看，它被描述为PolyMesh，这是对OpenFOAM网格及其所有功能和限制的一般描述。对于给定的静态网格情况，网格数据将存储在常量/多边形网格中。这里找到的典型网格数据文件包括：points、faces、owner、neighbour及boundary。当然，包含的数据必须有效，才能正确定义网格。

在下面的讨论中，以potentialFoam求解器的PitzDaily教程为例，可以通过发出以下命令找到该示例

?>  tut
?>  cd basic/potentialFoam/pitzDaily

检查PolyMesh目录的内容后，很明显其还不包含所需的网格数据。本教程中仅提供了blockMeshDict，在case目录中执行blockMesh会生成网格和关联的连接数据：

?>  ls constant/polyMesh
blockMeshDict boundary
?>  blockMesh
?>  ls constant/polyMesh
blockMeshDict boundary faces neighbour owner points

以前使用过CFD代码的用户，特别是使用基于结构化网格的代码的用户，可能会错过每个网格单元的寻址。

OpenFOAM中的非结构化FVM方法不是基于每个单元构建网格，而是基于每个网格面构建网格。以下列表说明了Constant/PolyMesh中每个文件的用途。

1、points

定义向量场中网格的所有节点，其中它们在空间中的位置以米为单位。这些点并非细胞中心$\mathbf{C}$，而是网格的角点。要将网格在正x方向上平移1米，必须相应地更改每个节点。不需要为此更改polyMesh子目录中的任何其他结构，节2.4介绍了这一点。

通过使用文本编辑器打开相应的文件，可以更仔细地查看这些点。为了限制输出，忽略标题，只显示前几行：

?>  head -25 constant/polyMesh/points | tail -7
25012 // 节点数量
((-0.0206 0 -0.0005) // 点0的坐标
(-0.01901716308 0 -0.0005) // 点1的坐标
(-0.01749756573 0 -0.0005)
(-0.01603868134 0 -0.0005)
(-0.01463808421 0 -0.0005)

该文件包含一个包含25012个点的列表。此列表不需要以任何方式排序。此外，列表中的所有节点都是唯一的，这意味着相同的点坐标不能多次出现。访问和寻址这些点是通过vectorField中的列表位置从0开始执行的。该位置存储为label。

2、faces

根据点在点向量场中的位置合成网格面，并将其存储在labelListList中。这是一个嵌套列表，每个面包含一个元素。这些元素中的每一个都有自己的labelList，存储用于构造面的points的标签。图2.1显示了labelListList的结构。

每个网格面必须至少由三个节点组成，其大小后面紧跟一个点标签列表。在网格面上，每个节点都通过一条直边与其相邻的点相连[4]。使用定义网格面的节点，可以计算表面积矢量${\mathbf{S}}_f$，其方向由右手定则确定。

同样，仅显示faces文件的前几行以保持简短：

?>  head -25 constant/polyMesh/faces | tail -7
49180 // 网格面的数量
(
4(1 20 172 153) // 网格面0包含4个节点，其标签为(1 20 172 153)
4(19 171 172 20)
4(2 21 173 154)
4(20 172 173 21)
4(3 22 174 155)
...
)

从输出的第一行可以看出，网格由49180个网格面组成，上面只显示了其中的一个子集。与面列表的长度49180类似，每个labelList的长度在列表开始之前声明。因此，此处显示的所有面都是从4个点构建的，这些点由它们在点列表中的位置表示。

3、owner

Owner也是一个与存储面的列表具有相同维度的labelList。由于面已构建并存储在面列表中，因此必须定义它们与体网格的从属关系。根据定义，一个网格面只能在两个相邻的网格之间共享。owner列表存储哪个面属于哪个网格，而这是根据网格标签来决定的。具有较低单元格标签的网格拥有该面，其另一个面被视为neighbor。它指示代码第一个面(列表中的索引0)属于标签存储在该位置的网格。

查看下面的owner文件可知，网格面0、1、2、3分别归属于网格单元0、0、1、1。owner文件中的网格面数量与faces文件中的网格面数量一致。

?>  head -25 constant/polyMesh/owner | tail -7

49180
(0
0
1
1

同样，上一章解释了owner-neighbour寻址的工作原理。

4、neighbour

neighbour必须与owner列表结合起来一起考虑，其与owner列表相反。neighbour存储相邻的网格，而不是定义哪个网格拥有每个特定的面。将owner文件与neighbour文件进行比较，可以发现它们的主要区别：owner列表要短得多。这是因为边界面没有相邻的网格。

?>  head -25 constant/polyMesh/neighbour | tail -7

24170
(1
18
2
19

5、boundary

边界包含有关嵌套子字典列表中网格边界的所有信息。边界通常被称为patch或边界patch。与之前的网格组件类似，仅显示了一些相关行：

?>  head -25 constant/polyMesh/boundary | tail -8

5
(
    inlet 
    {
        type        patch;
        nFaces        30;
        startFace    24170;
    }

对于本节中使用的 pitzDaily 示例，边界文件包含 5 个patch的列表。每个patch由一个字典表示，由patch名开始。字典中包含的信息包括：patch类型、面数和起始面。由于面列表的排序，可以使用此约定快速轻松地处理属于某个面片的面。

边界面的寻址方法如图 2.2 所示。根据设计，所有没有neighbor的面都被存放在faces列表的末尾，根据他们的owner patch对其进行排序。所有作为边界面的面都必须被边界描述所覆盖。

从用户的角度来看，点和面、owner和neighbor都不需要手动接触或操作。如果手动更改它们，这肯定会破坏网格。但是，根据工作流程，可能需要针对某些设置更改边界文件。更改边界文件的最可能原因是更改patch名称或类型。在此处进行此更改可能比重新运行相应的网格生成器要容易得多。

现在解释了 OpenFOAM 网格的基本结构，接下来将回顾边界类型。有几种类型可以分配给边界，其中一些比其他更常见。区分边界（或patch）和边界条件（见图 2.3）很重要。

patch是计算域的外边界，其在边界文件中指定，因此是一种拓扑属性。边界和 CAD 几何之间的逻辑联系是两者的表面应尽可能相同。用网格拓扑表示，它是一组面，只有一个owner网格，没有neighbor网格。与patch相反，边界条件分别应用于每个物理场（U、p 等）的patch。patch类型有：

• Patch。大多数patch（边界）可以用patch类型来描述，因为它是最一般的描述。 Neumann、Dirichlet 或 Cauchy 边界条件都可以应用于这种类型的边界。
• wall。如果patch被定义为wall，这并不意味着没有流体通过过该边界。它仅使湍流模型能够正确地将壁面函数应用于该patch（参见第 7 章），仍然需要通过速度边界条件明确定义防止流体通过类型为wall的patch。
• symmetryPlane。将patch类型设置为 SymmetryPlane 声明它充当对称平面。除了 SymmetryPlane 之外，不能对其应用其他边界条件，并且必须将其应用于所有物理场。
• empty。在二维模拟的情况下，这种类型应该应用于“平面内”的patch。与 SymmetryPlane 类型类似，这些patch的边界条件也必须为所有物理场设置为empty。不会对这些patch应用其他边界条件。两个empty边界之间的所有网格边必须平行，否则无法进行精确的二维模拟。
• cyclic。如果一个几何结构由多个相同的部件组成（例如螺旋桨叶片或涡轮叶片），则只需将其中一个组件离散化并将其视为位于相同组件之间。对于四叶片螺旋桨，这意味着只有一个叶片是网格化的（90° 网格），并通过将cyclic边界类型分配给具有切线方向法线的patch。然后这些patch将充当物理耦合。
• wedge。这种边界类型类似于cyclic边界，只是专门为形成小的（例如≤5°）楔形的循环边界设计的。

从执行和兼容性的角度来看，polyMesh 结构的创建方式并不重要，只要网格数据本身是有效的即可。虽然 OpenFOAM 打包了各种网格生成工具，但只要可以进行有效的转换或输出，就可以使用外部第三方网格生成器。

除了上面提到的 OpenFOAM 网格的基本核心组件之外，还有各种可选的网格结构，它们只能用于特定的应用程序。由于它们是可选的，因此无论案例设置如何，它们都可以存在。 OpenFOAM 应用程序会根据需要读取它们，并在它们丢失时向用户报告。

6、Sets及Zones

作为用户，当用户担心时，很容易被 OpenFOAM 中的区域(zones)和集合(sets)以及两者都非常相似的事实所混淆：他们选择网格实体。对使用哪个问题的非常简短的回答是：使用区域，正如 Hrvoje Jasak 通过 Twitter 简要解释的那样（参见图 2.4）。

然而，这仅与求解器应用程序真正相关。如果应用程序以预处理或后处理为中心，则任何一个都可以。Set本质上是 labelHashSets，而Zone继承自 labelLists。两者都可以将任何网格实体（点、面或单元）存储在类似于列表的数据结构中。主要区别在于网格实体的内部处理，特别是在具有拓扑网格变化的并行模拟的情况下。在这种情况下，必须相应地更新列表中的地址，并且只有Zone提供这种方法。

选择通常由工具 setSet 或 topoSet 执行，它们都可以选择网格的子集并对其进行布尔运算。一般来说，这两个实用程序都可以将区域转换为集合，反之亦然。可以为任何网格实体（单元、点或面）创建集合或区域，但 cellSet 和 cellZone 是最常用的两个。区域作为普通字典存储在 constant/polyMesh 中，而集合存储在 constant/polyMesh 的 sets 子目录中。区域和集合以相同的方式存储在文件系统中：作为相应网格实体的标签的长列表。

我们已经发布了一些博客文章，其中包含有关区域和集合如何组装的一定程度的信息 [2, 1]。

2.1.1 CAD几何

导入外部 CAD 软件中生成的几何图形是CFD 工程师的常规任务。在 OpenFOAM 中，这通常使用 snappyHexMesh 执行，但是，稍后将解释这个网格生成器的用法。目前唯一重要的概念是本节仅处理Stereolithography(STL) 文件的导入。支持其他文件类型并以类似的方式工作。 STL 是一种文件格式，可以以三角面片方式存储几何图形的表面。二进制和 ASCII 编码文件都是可能的，但为了简单起见，我们将使用 ASCII 编码。

作为 STL 文件的示例，以下代码段显示了仅由一个三角形组成的 STL 曲面：

solid TRIANGLE
    facet normal -8.55322e-19 -0.950743 0.30998
        outer loop
            vertex -0.439394 1.29391e-18 -0.0625
            vertex -0.442762 0.00226415 -0.0555556
            vertex -0.442762 1.29694e-18 -0.0625
        endloop
    endfacet
endsolid TRIANGLE

在此示例中，仅定义了一个名为 TRIANGLE 的实体。一个 STL 文件可能包含多个实体，这些实体一个接一个地定义。组成表面的每个三角形都有一个法线向量和三个点。

使用 ASCII STL 文件的缺点是它们的文件大小会随着表面分辨率的增加而迅速增长。边没有明确包含，因为文件中只存储了三角形。因此，从 STL 中识别和提取特征边缘有时是一项具有挑战性的任务。

使用 STL 作为文件格式的一个优点是可以获得三角形表面网格，根据定义，它总是具有平面表面组件（三角形）。

有很多专门为 OpenFOAM 设计的开源网格生成器，分布在两个主要开发分支（vanilla OpenFOAM 和foam-extend）中。这包括 blockMesh、snappyHexMesh、foamyHexMesh、foamyQuadMesh 和 cfMesh。还有一些其他的工具，如extrudeMesh 和 extrude2DMesh，但在本节中没有讨论，因为多数 OpenFOAM 用户不使用它们。此外，它们主要属于网格实用程序，而不是本章讨论的核心网格生成器。本节将简要介绍 blockMesh 和 snappyHexMesh，并回顾它们的用法和工作原理。一般来说，网格生成器的目的是以用户友好的方式生成上一节中描述的 polyMesh 数据结构。两个网格生成器具有相似的输入和输出，因为它们读取字典文件并将最终网格写入constant/polyMesh文件中。

2.2.1 blockMesh

当调用可执行的 blockMesh 时，会自动从 constant/polyMesh 目录中读取blockMeshDict ，因此该目录文件必须存在。

blockMesh 生成块结构的六面体网格，然后将其转换为 OpenFOAM 所需的任意非结构化格式。使用 blockMesh 为复杂的几何图形生成网格通常是一项非常乏味和困难的任务，有时甚至是不可能的。对于复杂的几何图形，用户生成 blockMeshDict 所花费的精力会大大增加。因此，通常只使用 blockMesh 生成简单的网格，然后将实际几何的离散化转移到 snappyHexMesh。这使得 blockMesh 成为生成网格的好工具，这些网格要么由相当简单的几何体组成，要么为 snappyHexMesh 生成背景网格。

blockMesh 用于构建网格的示例块如图 2.5 所示。每个块由称为顶点的 8 个角组成。六面体块是从这些角点及边线构建的。

如图 2.5 所示，将顶点相互连接。最后，块的表面由patch定义，尽管这些patch仅必须为没有相邻块的块边界显式指定。两个块之间的边界不得在patch定义中列出，因为根据定义，它们不是patch。特定边上的节点的长度和数量必须匹配，以保持拓扑一致。实际模拟的边界条件稍后将应用于这些patch。

请注意，可以生成少于 8 个顶点的块并且在patch上具有不匹配的节点（请参阅 [4]），但是，本指南不涵盖这一点。块的边缘默认为直线，但可以替换为不同的线类型，例如圆弧、多段线或样条线。选择例如圆弧确实会影响块边缘拓扑的形状，但该边缘上最终网格点之间的连接保持直线（见图 2.6）。

1、坐标系统

最终网格在全局（右手）坐标系中构建，该坐标系是笛卡尔坐标系并与主要坐标轴对齐：x、y 和 z。当块必须在空间中任意对齐和定位时，这会导致问题。为了避免这个问题，每个块都分配有自己的右手坐标系，根据定义，这并不要求三个轴正交。这三个轴标记为 x1、x2、x3（参见 [4] 和图 2.5）。根据图 2.5 所示的符号定义局部坐标系：顶点 0 定义原点，顶点对 (0, 1) 表示 x1，而 x2 和 x3 由顶点对 (0, 3) 和 (0, 4)分别定义。

2、节点分布

在网格划分过程中，每个块被细分为单元。单元由块坐标系的三个坐标轴中每个坐标轴的边缘上的节点定义，并遵循下式给出的关系： $$n_{cells}=n_{nodes}-1\text{\hspace{1em}(2.1)}$$ 用户可以在 blockMeshDict 中定义在某条边上将出现多少个单元格。边上的单元可以均匀分布，也可以基于分级分布在非均匀分布上。存在两种类型的分级：simpleGrading 和 edgeGrading based。 simpleGrading 描述了基于特定边缘上最后一个网格与第一个网格的大小比对边缘的分级（见图 2.7）：

$$e_r=\frac{{\delta }_e}{{\delta }_s}\text{\hspace{1em}(2.2)}$$ 如果 $e_r=1$，所有节点都在该指定边上均匀分布，则不存在分级。当膨胀比 $e_r>1$ 时，节点间距从边的开始到结束增加。从 blockMesh 的 C++ 源代码中可以发现，由用户 (er) 定义的扩展比由以下关系缩放： $$r=e_r^{\frac{1}{1-n}}\text{\hspace{1em}(2.3)}$$

其中$n$表示该指定边上的节点数。通第$i$个节点在边上的相对位置可以通过下式计算：

{\lambda(r, i)=\frac{1-r^{i}}{1-r^{n}} \quad \text { with } \lambda \in[0,1]}\tag{2.4}\label{2.4}

尽管对于BlockMeshdict中的所有块来说，这看起来太费力了，但当需要在相邻的两个块之间平滑地转换网格大小时，这就派上了用场。 在许多情况下，简单的试错通常就足够了。

3、为最小示例定义字典

作为如何正确设置BlockMeshdict的一个小例子，对一个体积为1m3的立方体进行了离散化。 可以在示例案例存储库的chapter2/blockMesh目录中找到准备好的案例。字典由一个关键字和四个子字典组成。 第一个关键字是convertTometers，通常为1。 所有的点位置都是按这个因子缩放的，如果几何尺寸很大或很小，这就派上了用场。 在任何这种情况下，我们最终都会键入大量的前导零或后导零，这是一项乏味的任务。 通过相应地设置convertTometers，我们可以节省一些键入。blockmeshdict的第一个相关行是：

convertToMeters 1;

其次，必须定义顶点。 重要的是要记住BlockMesh中的顶点不同于所创建的PolyMesh中的点，尽管它们的定义相当相似。 对于单位立方体示例，顶点定义为

vertices
(
    (0 0 0)
    (1 0 0)
    (1 1 0)
    (0 1 0)
    (0 0 1)
    (1 0 1)
    (1 1 1)
    (0 1 1)
);

通过查看上面的定义，可以清楚地看出该语法是一个列表，类似于polymesh定义中的点列表。 这是由于在OpenFoam中圆括号指示列表，而花括号则定义字典。 前四行定义了$x_3=0$平面中的所有四个顶点，下面的行定义了$x_3=1$平面中的所有四个顶点。 与polymesh中的点类似，每个元素都是通过其在列表中的位置而不是坐标来访问的。 请注意，每个顶点必须是唯一的，因此在列表中只出现一次。

下一步，必须定义这些块。 图2.5可以作为参考。 单位多维数据集的示例块定义如下所示：

blocks
(
    hex (0 1 2 3 4 5 6 7) (10 10 10) simpleGrading (1 1 1);
);

同样，由于圆括号的存在，这是一个包含块的列表，而不是字典。 这个定义乍一看可能有点奇怪，但实际上很直截了当。 第一个单词hex和第一组包含八个数字的圆括号告诉blockMesh在顶点0到7中生成一个六面体。 这些顶点正是上面顶点部分中指定的顶点，并通过它们的标签进行访问。 它们的顺序不是任意的，而是由如下所示的局部块坐标系定义的：

1. 对于局部$x_3=0$平面，列出所有从原点开始并按照右手坐标系移动的四个顶点标记。
2. 对局部$x_3\ne 0$平面执行相同操作

通过打乱特定块定义中顶点列表的顺序可以获得有效的块定义。 生成的块将看起来扭曲或使用不正确的全局坐标方向。 一旦执行blockMesh和checkMesh并在后处理器（例如Paraview）中分析网格，就检测到这一点。

注：checkMesh是一种原生的OpenFOAM工具，用于根据各种标准检查网格的完整性和质量。如果checkMesh的输出指出网格不正常，则必须对其进行改进。

第二组圆形括号定义了在块的每个特定方向上分布多少单元格。 在这个示例中，块每个方向有10个网格格。 若修改$x_1$方向为2个网格，$x_2$中为20个网格，$x_3$中为1337个单元格，则块定义可以写为：

hex (0 1 2 3 4 5 6 7) (2 20 1337) simpleGrading (1 1 1);

最后剩下的是SimpleGrading部分，与圆括号中的最后一组数字结合在一起。 如前所述，这是定义分级（或扩展比率）的最简单的方法。 在本例中，关键字simpleGrading定义了三个局部坐标系轴方向上所有四条边的等级分为相同。 因此simpleGrading后面括号中的三个数字分别定义了$x_1$、$x_2$和$x_3$方向的分级。 不过，有时这还不够通用。 这就是可以使用edgeGrading的地方。 这种更先进的分级方法本质上与simpleGrading相同，但可以明确指定六面体上12条边中每条边的分级，此时最后一组括号不会列出3个数字，而是3乘以4。 现在，每个边都可以单独设置。

保存blockMeshdict文件并在之后执行blockMesh，将得到一个与blockMeshdict中定义的类似的有效网格。 但是blockMesh会警告未定义的patches，这些patch默认被放入DefaultFaces中。

手动指定patch是通过在名为patches的列表中定义它们来完成的，对于示例patch 0：

patches
(
    XMIN
    {
        type patch;
        faces
        (
            (4 7 3 0)
        );
    }
);

这将指示blockMesh根据从顶点4、7、3和0构造的面生成名为Xmin的patch类型的patch。在内部，patch名称定义为word，并且此数据类型会定期显示在错误消息中。不过，顶点的排序方式不是任意的。它们需要从块内部看以顺时针方向指定。图2.8显示了示例的单位立方体的图像，该立方体由1000个小立方体组成，带有突出显示的XMIN，YMIN和ZMAX的patch。用于生成此网格的文件可以在chapter2/blockMesh下的示例存储库中找到。

如前所述，缺省情况下，块的边是行，因此包含边定义的列表是可选的。 与上面定义的块和patch非常类似，通过弧线而不是默认线连接两个顶点将如下所示：

edges
(
    arc 0 1 (0.5 -0.5 0)
);

包含边定义的列表中的每个项都以指示边类型的关键字开始，然后是开始和结束顶点的标签。 在本例中，这条线被构造弧所需的第三个点闭合。 对于任何其他边形状（例如折线或样条），该点将被一系列支撑点所取代。

插入上面列出的代码如何改变单位立方体的形状（参见图2.8)，在图2.9中给出了一个示例。

要继续SnappyHexMesh部分，需要生成一个由每个方向的50个单元格组成的单位立方体。

2.2.2 snappyHexMesh

与BlockMesh相比，snappyHexMesh可能不需要太多繁琐的工作，比如添加和连接块。 另一方面，对最终网格的控制较少。 使用snappyHexMesh，可以轻松生成六面体占优网格，只需要两件事：一个六面体背景网格，第二个是一个或多个兼容曲面格式的几何图形。 snappyHexMesh支持由各种体积形状定义的局部网格细化（见表2.1)，边界层单元（棱镜和多面体）的应用以及并行执行。

snappyHexmesh是一个复杂的程序，由大量的控制参数控制。 详细描述这些超出了本书的范围。 请结合本书阅读[4]，以获得对snappyHexMesh更深入的讨论。 其他信息可以在这里[5]找到。

snappyHexMesh的执行流程可以分为三个主要步骤，然后依次执行。可以通过在snappyHexMeshDict的开头将相应的关键字设置为false来禁用这些步骤中的每个步骤。这三个步骤可以概括如下：

1. castellateMesh。这是第一阶段，执行两个主要操作。首先，它将几何形状添加到网格中，并删除不在流动域内的单元。其次，根据用户的指定对现有单元格进行拆分和细化。结果是一个仅由或多或少类似于几何形状的六面体组成的网格。但是，应该放置在几何图形表面上的大多数网格点都不与它对齐。在图2.10中示出了在网格划分过程的这一阶段的后面示例的屏幕截图。
2. snap。通过执行捕捉步骤，表面附近的网格点被移动到表面上。 这一点在图2.11中可见一斑。 在此过程中，这些网格的拓扑结构可能由六面体变为多面体。 靠近表面的网格可能被删除或合并在一起。

3. addLayers。最后，在几何表面上引入额外的单元，通常用于细化近壁流动（见图2.12）。 将预先存在的单元从几何图形中移开，以便为额外的单元创造空间。 那些网格很可能是棱柱层。

上述所有设置以及更多设置都在system/snappyHexMeshDict中定义，该字典包含snappyHexMesh所需的所有参数。在网格化/snappyHexMesh下的OpenFOAM教程目录中可以找到一些有用的教程。与其他OpenFOAM字典相比，snappyHexMeshDict非常长，由许多层次结构级别组成，这些层次结构级别由嵌套的子词典表示。对于上述每个步骤 (假设您具有标准配置)，将一个时间步骤写入case目录。这三个步骤中的每一个将在下一节中单独讨论。

2.2.2.1 Cell levels

网格级别用于描述背景网格单元格的细化状态。 启动snappyHexMesh时，将读取背景网格，并将所有单元格分配为单元格级别0（图2.12中为蓝色单元格）。

如果一个网格被细化一个级别，每个边缘被减半，从前一个“父”单元格生成八个单元格。 这种细化方法基于八叉树，只适用于六面体，这就是为什么snappyHexmesh需要六面体背景网格。 使用snappyHexmesh不可能只在一个方向上细化单元格，因为八叉树无法覆盖这一点。 因此根据定义，它们在所有三个空间方向上都是一致的。

2.2.2.2 定义几何图形

在snappyHexMeshdict中定义任何内容，constant/polymesh中的现有网格将自动读取并用作背景网格。 如果没有这样的网格可用，或者如果它不是纯粹基于六面体，snappyHexmesh将无法运行。 对于外部流动模拟，由背景网格定义的外部边界不如内部流动重要。 因此，它们可以保持由背景网格定义，而不需要在它们上花费更多的工作。 另一方面，在内部流动模拟中，背景网格的外部形状是由实际几何形状定义的，因此它是不受关注的。

STL几何图形可以使用几乎任何CAD程序生成。 Paraview可用于生成基本形状的STL表示，如圆柱、球体或锥体。 在sources菜单下，可以使用file菜单下的save data项导出各种形状。

对于现实世界的几何学，当然有各种方法来生成曲面网格并将其存储为STL。 然而，请记住，表面网格的质量对获得良好的体积网格至关重要。

作为一个简单的例子，在前一节中准备的单位立方体网格被重用，并在其中插入一个球体。 球体是使用STL文件生成的，而不是表2.1中列出的形状。 加载STL几何图形可以以直接的方式完成，只需将几何图形复制到case的constant/triSurface，并在snappyHexMeshdict中添加以下几何图形子字典。 一个这样的示例如下所示：

geometry
{
    sphere.stl // Name of the STL file
    {
        type triSurfaceMesh; // Type that deals with STL import
        name SPHERE; // Name access the geometry from now on
    }
}

上面的行告诉snappyHexMesh从constant/triSurface中读取spher.stl作为triSurfaceMesh，并将该STL中包含的几何体引用为sphere。一些简单的几何对象可以在不需要打开任何CAD程序的情况下直接在snappyHexMesh中构建。表2.1列出了这些几何形状。

表2.1中列出的任何形状都可以在geometry子字典中构造，方法是简单地追加到现有的子字典中。例如，要将一个box添加到几何子字典中，该子字典由最小点和最大点构成。使用此方法时，不可能直接旋转长方体，它将始终与坐标轴对齐。

smallerBox
{
    type searchableBox;
    min (0.2 0.2 0.2);
    max (0.8 0.8 0.8);
}

与STL定义类似，定义searchableBox的子字典的前导字符串是用于稍后访问该几何体的名称。有时候，我们希望用表2-1中列出的形状组合一个几何体，但要将其视为一个几何体而不是多个几何体。这是可以使用searchableSurfaceCollection的地方。通过在已经存在的几何体组件上使用此方法，可以组合、旋转、平移和缩放曲面。在任何情况下，将SPHERE和smallerBox合并为一个，并将fancybox放大2倍，将如下所示：

geometry
{
    ...
    fancyBox
    {
        type searchableSurfaceCollection;
        mergeSubRegions true;
        SPHERE2
        {
            surface SPHERE
            scale (1 1 1);
        }
        smallerBox2
        {
            surface smallerBox;
            scale (2 2 2);
        }
    }
}

2.2.2.3 设置castellatedMesh步骤

这是执行snappyHexMesh期间三个步骤中的第一个步骤。它包括以下两个主要步骤：根据用户规范分割单元并删除网格区域之外的单元。图2.13给出了这一过程的示意图。

现有的背景网格（图2.13中的黑色）从constant/polyMesh中读取。根据snappyHexMeshDict的castellatedMeshControls子字典中的参数细化网格。区分由几何曲面定义的细化和体积细化是很重要的。曲面细分可确保表示几何图形的边界面细分到定义的层级。需要注意的是，这不仅会影响拥有特定单元的单元，还会影响邻近的单元。因此，表面细化可能看起来有点类似于体积细化，但是，它是明显不同的。对SPHERE应用这样的曲面细化由castellatedMeshControls中的条目控制，如下所示：

castellatedMeshControls
{
    ...
    refinementSurfaces
    {
        SPHERE // Name of the surface
        {
            level (1 1); // Min and max refinement level
        }
    }
    ...
}

这会将球体的表面细化到1级。圆括号之间的两个数字定义了此曲面的最小和最大细化级别。snappyHexMesh根据曲面曲率在两者之间进行选择：高度弯曲的表面区域被细化到较高水平，较小弯曲的表面区域被细化到较低水平。

snappyHexMesh中的细化不限于曲面定义。在几何子字典中定义的任何几何也可以用作体积细化的定义形状。这些体积细化称为refinementRegions，并在castellatedMesh控件中的同名子字典中定义。与refinementSurfaces相比，refinementRegions提供了更高级别的多功能性，因此需要定义更多选项。

模式有三种选择：inside、outside和distance。 顾名思义，inside只影响选定几何体内的单元格，而outside则恰恰相反。 第三个选项，距离，是两者的组合，并在表面的向外和向内法线方向上计算。 除了modes，还有一个levels选项，它比refinementSurface更复杂。 从名称上已经可以猜到，它确实支持任意数量的级别。 每个级别必须结合distance来定义。 随着在列表中位置的增加，级别必须减少，distance必须增加。 将smallerBox中的任何内容细化到级别1可以通过向castellatedMeshControls添加以下行来完成：

castellatedMeshControls
{
    ...
    refinementRegions
    {
        smallerBox // Geometry name
        {
            mode inside; // inside, outside, distance
            levels ((1E15 1)); // distance and level
        }
    }
    ...
}

上述代码使用$1\times 10^{15}$ m的距离，以便安全地选择几何体内部的所有网格单元。

如果不指定位于最终网格体积内的点，snappyHexMesh就无法确定用户要离散化球体的哪一侧。这就是为什么locationInMesh关键字也必须在castellatedMeshControls子字典中定义的原因。此点不得放置在背景网格的面上。对于单位立方体示例，该点定义为：

locationInMesh (0.987654 0.987654 0.987654);

下一步是调整snappyHexMeshDict中snap子字典的参数。

2.2.2.4 设置捕捉步骤

与snappyHexMesh的其他两个步骤相比，这不需要大量的用户输入。此步骤负责通过将新点引入网格并替换它们来将纯六面体网格面与几何体对齐（请参见图2.11）。这是一个高度迭代的过程，这就是为什么不需要太多用户交互的原因。snappyHexMeshDict的示例snapControls子字典如下：

snapControls
{
    nSmoothPatch 3;
    tolerance 2.0;
    nSolveIter 30;
    nRelaxIter 5;
    // Feature snapping
    nFeatureSnapIter 10;
    implicitFeatureSnap false;
    explicitFeatureSnap true;
    multiRegionFeatureSnap false;
}

只定义了迭代计数器、容差和标志。 一半的参数处理到几何体边缘的对齐，这不是本描述的一部分。 然而，对此的描述可以在[3]中找到。 根据具体情况，增加迭代计数器通常会得到更高质量的网格，但也会显著增加网格划分时间。

所有参数在OpenFOAM提供的snappyHexMeshDicts中有更详细的解释。

2.2.2.5 设置addLayers步骤

addLayers步骤的所有设置都在snappyHexMeshDict的addLayersControls子字典中定义。可以使用任何表面从挤出棱柱层，而不论其类型为何。首先，需要通过layers子字典指定每个边界要拉伸的单元层数。示例条目如下所示：

addLayersControls
{
    ...
    layers
    {
        "SPHERE_.*" // Patch name with regular expressions
        {
            nSurfaceLayers 3; // Number of cell layers
        }
    }
    ...
}

每个patch名称后跟一个包含nSurfaceLayers关键字的子字典。此关键字定义要拉伸的网格层数，因此后面跟一个整数，表示要拉伸的网格层数。在上面的示例中，使用正则表达式来匹配以SPHERE_开头的任何patch程序名称。在这种情况下，它只是球体本身，但是以这种方式使用通配符可以大大减少设置时间。最终网格的横截面如图2.12所示。

snappyHexMesh的各种参数（尤其是与层挤出相关的参数）需要进行调整，以获得符合要求的网格。下面简要解释其中的几个。

• relativeSizes可以针对以下值从绝对标注切换到相对标注。默认情况下为true。
• expansionRatio定义从一个网格层到下一个网格层的扩展因子。
• finalLayerThickness是最后一个网格层（距离壁面最远）相对于网格的下一个网格的厚度，或以绝对米为单位，具体取决于对relativeSizes参数的选择。
• minThickness如果层的厚度不能大于minThickness，则不会拉伸该层。

在该示例中，采用了以下所示的设置。

relativeSizes            true;
expansionRatio            1.0;
finalLayerThickness      0.5;
minThickness            0.25;

最后，必须在case目录中执行snappyHexMesh以开始网格化过程。每个步骤都会生成一个新的时间点目录，其中包含该特定阶段的网格。如果选择通过调整snappyHexMeshDict中的参数来更改网格，请记住在重新运行snappyHexMesh之前删除旧的时间点。

1. 另一个用于OpenFOAM的高质量网格生成器是enGrid，可以从www.example.com免费获得http://engits.eu/en/engrid。
2. cfMesh应用程序是一个分布式内存并行OpenFOAM网格工具，它与snappyHexMesh一样，将STL曲面作为输入。此软件包由Creative Fields Ltd.开发，可在www.c-fields.com上与文档一起下载。

2.2.3 cfMesh

cfMesh库是一个跨平台库，用于自动生成网格，它构建在OpenFOAM之上。它与OpenFOAM和foam-extend的所有最新版本兼容，并根据通用公共许可证（GPL）授权。该库由Franjo Juretić博士开发，并由Creative Fields Ltd.发布，可从http://www.c-fields.com 下载。

本节简要概述了库、控制网格生成过程的选项以及cfMesh附带的一些实用程序。它绝不是项目的完整文档。有关更多信息，请访问上述项目网页。

cfMesh库支持使用主库中的组件构建的各种3D和2D工作流，这些组件是可扩展的，可以组合到各种网格化工作流中。核心库基于网格修改器的概念，可通过消息传递接口（MPI）使用对称多处理器（SMP）与分布式内存并行（DMP）实现高效的并行化。此外，还特别注意内存使用，通过实现数据容器（列表、图形等）来保持较低的内存使用。其在网格化处理期间不需要许多动态存储器分配操作。

cfMesh中的网格化过程是自动进行的，需要输入三角剖分和包含各种网格化参数（设置）的字典。给定曲面网格和设置后，网格化过程将从控制台启动，并且自动运行，无需任何用户干预。该库经过优化，网格化工作流需要的设置较少，并且语法简单。目前，cfMesh可以在内部创建体积网格，请参见图2.14，它不需要几何封闭。

2.2.3.1 可用网格化工作流

所有工作流都针对共享内存计算机进行并行化处理，并在运行时使用所有可用的CPU内核。使用的核心数可由OMP_NUM_THREADS环境变量控制，该变量可设置为所需的核心数。

可用的网格划分工作流通过从输入几何图形和用户指定的设置创建所谓的网格模板来启动网格划分过程。 模板随后被调整以匹配输入几何图形。 将模板拟合到输入几何形状的过程被设计为能够容忍质量差的输入数据，这不需要封闭几何。 可用的工作流因模板中生成的网格类型而不同。

笛卡尔工作流生成的三维网格主要由六面体单元组成，不同尺寸单元之间的过渡区域为多面体。 它是通过在shell窗口中键入CartesianMesh开始的。 默认情况下，它生成一个边界层，可以根据用户的要求进一步细化。 此外，该工作流可以使用MPI并行化来运行，该并行化用于生成不适合于单个可用计算机内存的大型网格。

该工作流将生成2D笛卡尔网格。在控制台中键入cartesian 2DMesh可启动网格生成器。依预设，它会产生一个边界层，可进一步细分。此网格化工作流需要带状几何体，如图2.16a所示，该几何体在x-y平面中延伸并在z方向上拉伸。

四面体工作流生成由四面体单元组成的网格，如图2.17所示，并通过在控制台中键入tetMesh启动。默认情况下，它不生成任何边界图层，并且可以根据用户请求添加和优化边界层。

2.2.3.2 输入几何

cfMesh使用的几何尺寸需要以曲面三角剖分的形式定义。对于2D情况，几何图形以三角形带的形式给出，其边界边位于x-y平面中（不支持其他方向）。几何由下列对象组成：

• List of points。包含曲面三角剖分中的所有点。
• List of triangles。包含曲面网格中的所有三角形。
• Patches。在网格划分过程中转移到体网格上的实体。 曲面中的每一个三角形都分配给单个的patch，不能分配给多个patch。 每个patch由其名称和类型标识。 默认情况下，所有patch名称和类型都被转移到体积网格中，并且很容易用于定义模拟的边界条件。
• Facet subsets。在网格划分过程中没有转移到体网格上的实体。 它们用于定义网格划分设置。 每个面子集包含曲面网格中三角形的索引。 请注意，曲面网格中的三角形可以包含在多个子集中。 刻面子集可以由cfSuite生成，cfSuite是由Creative物理场有限公司开发的商业应用程序。
• Feature edges。特征边在网格化过程中被视为约束。三条或更多特征边相交的曲面点被视为角点。特征边可以通过surfaceFeatureEdges工具或cfSuite生成。

图2.18显示了一个带有高亮显示的patch的曲面网格。

由cfMesh转换的所有尖锐特征必须在网格化过程之前由用户定义。在网格划分过程中，两个面片之间的边界边缘（图2.19a）和特征边缘作为尖锐特征进行处理（见图2.19b）。三角剖分中的其他边不受约束。

网格的文件格式建议：fms、ftr和stl。此外几何可以在所有支持的格式导入带有OpenFOAM surfaceConvert效用。然而这三个建议格式支持定义的patch转移到默认体网格。其他格式还可以用于网格划分但是他们不支持定义输入几何和patch的面上产生的体积网格边界的结束在一个patch。

cfMesh的首选格式为fms，用于保存设置网格作业的所有相关信息。它将补片、子集和特征边存储在单个文件中。此外，它是唯一一种可以将所有几何实体存储到单个文件中的格式，强烈建议用户使用它。

2.2.3.3 字典和可用的设置

网格划分过程由位于案例的系统目录中的meshDict字典中提供的设置来指导。 对于使用MPI的并行网格划分，需要一个位于case的系统目录中的decomposePardict，并且用于并行运行的节点数必须与decomposepardict中的numberofsubdomains条目匹配。 decomposePardict中的其他条目不是必需的。 生成的卷网格写在constant/polyMesh目录中。 meshDict中可用的设置将在本节的其余部分中更详细地解释。

cfMesh库只需要两个强制设置就可以启动网格划分过程：

• surfaceFile。指向几何文件。 几何文件的路径相对于case目录的路径。
• maxCellSize。表示用于网格划分的默认网格大小。 它是域中生成的最大单元格大小。

2.2.3.4 细化设置

当统一的单元大小不令人满意时，cfMesh中有许多用于局部细化源的选项。

• boundaryCellSize选项用于优化边界处的单元格。这是一个全局选项，所请求的网格大小将应用于边界的所有位置。
• minCellSize是一个全局选项，可激活网格模板的自动细化。此选项在像元大于估计要素大小的区域中执行优化。此设置提供的标量值指定此过程可生成的最小像元大小。此选项对于快速模拟非常有用，因为它可以在复杂几何图形中生成网格，用户只需很少的工作。但是，如果需要高网格质量，它会在需要进行网格细化的位置提供提示。
• localRefinement允许在边界处进行局部细化区域。它是一个字典的字典，并且localRefinement主字典中的每个字典都由用于细化的几何图形中的面片或面子集命名。实体的请求像元大小由cellSize关键字和标量值控制，或通过指定additionalRefinementLevels关键字和相对于maxCellSize的所需细化数来控制。
• objectRefinement用于指定体积块内的细化区域。可用于优化的支持对象包括：直线、球体、长方体和截锥。它被指定为字典的字典，其中objectRefinement字典中的每个字典都表示用于细化的对象的名称。

cfMesh中实现的网格化工作流基于由内向外网格化，网格化过程从基于用户指定的单元尺寸生成所谓的网格模板开始。但是，如果网格尺寸局部大于几何要素尺寸，则可能导致该几何形状由网格填充。相反，如果指定的网格大小大于局部特征大小，则几何体中较薄部分的网格可能会丢失。

keepCellsIntersectingBoundary选项是一个全局选项，可确保模板中与边界相交的所有单元仍是模板的一部分。默认情况下，所有网格化工作流仅保留模板中完全位于几何体内部的单元。keepCellsIntersectingBoundary关键字后面必须跟1（活动）或0（非活动）。激活此选项可能会导致局部连接的网格超过间隙，该问题可以通过checkForGluedMesh选项解决，该选项后面还必须跟有1（活动）或0（非活动）。

keepCellsIntersectingPatches选项是在用户指定的区域中保留模板中的单元的选项。它是字典的字典，主字典中的每个字典都以patch或面子集命名。启用keepCellsIntersectingBoundary选项时，此选项处于非活动状态。

removeCellsIntersectingPatches选项是从模板中删除用户指定区域中的单元格的选项。 它是一个字典的字典，主字典中的每个字典都由一个patch或一个方面子集命名。 当KeepCellsIntersectingBoundary选项打开时，该选项处于活动状态。

cfMesh中的边界层从体积网格的边界面向内部拉伸，在网格划分过程之前不能拉伸。 此外，它们的厚度由边界处指定的单元尺寸控制，网格生成器倾向于生成与单元尺寸相似厚度的层。 CFMesh中的层可以跨越多个补丁，如果它们共享凹边或角，且价大于3。 此外，CFMESH从不对边界层的拓扑进行Breask,其最终几何形状依赖于光滑过程。 所有边界层设置都在BoundaryLayers字典中提供。 选项有：

• nlayers指定将在网格中生成的层数。 它不是强制性的。 如果没有指定，网格划分工作流将生成默认的层数，该层数要么为一，要么为零。
• thicknessRatio是连续两层厚度之间的比值。 它不是强制性的。 该比率必须大于或等于1。
• maxFirstLayerThickness确保第一边界层的厚度永远不会超过指定值。 它不是强制性的。
• patchBoundaryLayers设置是一个字典，用于指定单个斑块边界层的局部属性。

可以在名称与patch名称相同的字典中分别为每个面片指定nLayers，thicknessRatio与maxFirstLayerThickness选项。默认情况下，在面片处生成的层数由全局层数控制，或由在与现有面片一起形成连续层的任何面片处指定的最大层数控制。allowDiscontinuity选项可确保patch所需的层数不会扩展到同一层中的其他patch。

本节中提供的设置用于在网格生成过程中更改patch名称和类型。这些设置在renameBoundary字典中提供，其中包含以下选项：

• newPatchNames是renameBoundary字典内的字典。它包含带有应重命名的patch名称的词典。对于每个patch，可以使用以下设置指定新名称或新patch类型：
  • newName关键字后跟给定patch的新名称。该设置不是必需的。
  • type关键字后跟给定patch的新类型。该设置不是必需的。
• defaultName是除newPatchNames字典中指定的patch之外的所有patch的新名称。该设置不是必需的。
• defaultType为所有patch（newPatchNames目录中指定的patch除外）设置新类型。该设置不是必需的。

2.2.3.5 cfMesh中的各种实用程序

目前，cfMesh项目提供了以下实用程序：

• FLMAToSurface将几何体从AVL的flma格式转换为cfMesh可读的格式。输入文件中定义的单元选择将作为面子集进行传输。
• FPMAToMesh是用于从AVL的fpma格式导入体积网格的实用程序。在输入网格上定义的选择将作为子集传递。
• copySurfaceParts将指定多面子集中的曲面多面复制到新曲面网格中。
• extrudeEdgesInto2DSurface将几何图形中作为特征边写入的边拉伸到生成2D网格所需的三角形带中。生成的三角形存储在单个面片中。
• meshToFPMA将网格转换为AVL的fpma格式。
• patchesToSubsets将几何中的曲面片转换为多面子集。
• preparePar创建MPI并行化所需的处理器目录。处理器目录的数目取决于在decomposeParDict中指定的numberOfSubdmains。
• removeSurfaceFacets是用于从曲面网格中移除镶嵌面的工具。应移除的面由面片名称或面子集给出。
• subsetToPatch在曲面网格中创建由给定小平面子集中的小平面组成的patch。
• SurfaceFeatureEdges用于生成几何中的特征边。 如果输出是FMS文件，则生成的边缘存储为特征边缘。 否则，它将生成以所选特征边缘为界的patch。
• SurfaceGenerateBoundingBox在几何图形周围生成一个box。 它不会解决自交集，以防盒子与几何的其余部分相交。

虽然blockMesh和snappyHexMesh是强大的网格生成工具，但用户可能经常使用第三方网格划分工具来定义和离散更复杂的流域。

2.3.1 从第三方网格软件包转换

许多高级外部网格化实用程序在网格生成过程中为用户提供了额外的控制级别。这包括可选元素类型、拟合边界层网格和长度比例控制等。某些网格生成器可以直接导出为功能OpenFOAM网格格式。下面列出了OpenFOAM-3.0中支持转换的网格格式的汇编：

• Ansys
• CFX
• Fluent
• GMSH
• Gambit
• Ideas
• Kiva
• Netgen
• Plot3D
• Star-CD
• tetgen
• KIVA

导入实用程序的功能以及它们使用的命名法差别很大。Fluent导入工具将内部边界转换为faceSet，而其他工具完全忽略此类特征。

由于许可证问题，用于从Star CCM+导入网格的网格转换工具及其相关库需要手动下载和编译，而不是通过Allrun脚本。

如果您的特定网格化软件未在上述列表中提及，则它很可能能够将网格导出为支持的中间格式。

上面提到的所有转换实用程序的源代码可以在这里找到:$W{M}_{P}ROJEC{T}_{D}IR/applications/utilities/mesh/conversion/。用户还可以选择使用foamMeshToFluent和foamToStarMesh实用工具将OpenFOAM网格转换为Fluent或Star-CD网格格式。这对于导出从前面提到的snappyHexMesh实用程序生成的网格特别有用。网格转换过程通常非常简单，不同转换实用程序之间的语法变化很小。出于这个原因，将只给出一个例子，并将基于fluentMeshToFoam转换实用程序。要开始这个过程，请将教程案例复制到选择的目录中，本教程基于icoFoam求解器现有的网格转换教程。‘‘‘bash?>cp-r$FOAM_TUTORIALS/incompressible/icoFoam/elbow/ meshConversionTest ?> cd meshConversionTest

转换网格非常简单，只需运行转换实用程序并将网格文件作为参数传递，该参数必须存在于目录中。 在转换过程中，Utility将向控制台输出补丁名称和网格统计信息。 Polymesh目录中包含的文件将相应更新。 

```bash
?>  fluentMeshToFoam elbow.msh

重要的是要记住，导入的网格只和导出的网格一样好。 对于Fluent网格，由于OpenFoam只支持三维网格，所以不可能导入2D网格。 导入完成后，需要更新case，以反映初始和边界条件文件中的新修补程序名称。 所有现有的补丁可以从导入工具的输出中收集，也可以通过编辑器打开constant/polymesh/boundary手动查找。 对于本教程，U和P物理场是为这个特定的网格预先配置的。

在导入过程中缩放网格就像在命令中添加选项和缩放因子一样简单。 为了本教程，网格应该缩小一个数量级。

?>  fluentMeshToFoam -scale 0.1 elbow.msh

在许多第三方网格划分工具中构造网格时，用户通常可以为每个特定的贴片分配边界条件类型，如入口、出口、壁面等。 转换过程将尝试将某些边界条件格式匹配到相应的OpenFoam格式，但不能保证边界条件转换的成功或准确性。 检查转换是否正确解析了流信息是至关重要的。 要检查这一点，请检查constant/polymesh/boundary，并在新转换的网格上运行checkmesh。

2.3.2 从二维网格到轴对称网格的转换

为了将网格转换为轴对称网格，必须满足以下要求。 网格必须已经是一个有效的OpenFOAM网格，它必须只有一个网格“厚”。 后一个要求对OpenFOAM中的所有二维网格都有效。 由于icoFoam的腔体示例满足所有这些要求，因此本教程使用它。 它位于$foam_tutorials/incompressible/icofoam/cavity。

由于网格形状为矩形并且只有一个单元，因此可以创建非常基本的几何体。在OpenFOAM中，轴对称网格具有以下属性：网格为一个单元“厚”，并绕对称轴旋转以形成5°楔形。楔形体的两个成角度边界被视为楔形体类型的单独patch。

makeAxialMesh的来源可在OpenFOAM wiki上找到：http://openfoamwiki.net/index.php/Contrib_MakeAxialMesh。按照其中的说明下载并编译该实用程序。

接下来的步骤是在您选择的工作目录中创建案例文件夹的副本，重命名该目录以避免将来发生任何混淆，然后创建2D基础网格。

?>  cp -r <span class="katex"><span class="katex-html" aria-hidden="true"><span class="base"><span class="strut" style="height:1em;vertical-align:-0.25em;"></span><span class="mord mathnormal" style="margin-right:0.02778em;">FO</span><span class="mord mathnormal">A</span><span class="mord"><span class="mord mathnormal" style="margin-right:0.10903em;">M</span><span class="msupsub"><span class="vlist-t vlist-t2"><span class="vlist-r"><span class="vlist" style="height:0.3283em;"><span style="top:-2.55em;margin-left:-0.109em;margin-right:0.05em;"><span class="pstrut" style="height:2.7em;"></span><span class="sizing reset-size6 size3 mtight"><span class="mord mathnormal mtight" style="margin-right:0.13889em;">T</span></span></span></span><span class="vlist-s">​</span></span><span class="vlist-r"><span class="vlist" style="height:0.15em;"><span></span></span></span></span></span></span><span class="mord mathnormal" style="margin-right:0.10903em;">U</span><span class="mord mathnormal" style="margin-right:0.00773em;">TOR</span><span class="mord mathnormal" style="margin-right:0.07847em;">I</span><span class="mord mathnormal">A</span><span class="mord mathnormal">L</span><span class="mord mathnormal" style="margin-right:0.05764em;">S</span><span class="mord">/</span><span class="mord mathnormal">u</span><span class="mord mathnormal">t</span><span class="mord mathnormal">i</span><span class="mord mathnormal" style="margin-right:0.01968em;">l</span><span class="mord mathnormal">i</span><span class="mord mathnormal">t</span><span class="mord mathnormal">i</span><span class="mord mathnormal">es</span><span class="mord">/</span><span class="mord mathnormal">in</span><span class="mord mathnormal">co</span><span class="mord mathnormal">m</span><span class="mord mathnormal">p</span><span class="mord mathnormal">ress</span><span class="mord mathnormal">ib</span><span class="mord mathnormal" style="margin-right:0.01968em;">l</span><span class="mord mathnormal">e</span><span class="mord">/</span><span class="mord mathnormal">i</span><span class="mord mathnormal">co</span><span class="mord mathnormal" style="margin-right:0.13889em;">F</span><span class="mord mathnormal">o</span><span class="mord mathnormal">am</span><span class="mord">/</span><span class="mord mathnormal">c</span><span class="mord mathnormal">a</span><span class="mord mathnormal" style="margin-right:0.03588em;">v</span><span class="mord mathnormal">i</span><span class="mord mathnormal">t</span><span class="mord mathnormal" style="margin-right:0.03588em;">y</span><span class="mord">.</span><span class="mclose">?</span><span class="mspace" style="margin-right:0.2778em;"></span><span class="mrel">&gt;</span><span class="mspace" style="margin-right:0.2778em;"></span></span><span class="base"><span class="strut" style="height:0.8889em;vertical-align:-0.1944em;"></span><span class="mord mathnormal">m</span><span class="mord mathnormal" style="margin-right:0.03588em;">v</span><span class="mord mathnormal">c</span><span class="mord mathnormal">a</span><span class="mord mathnormal" style="margin-right:0.03588em;">v</span><span class="mord mathnormal">i</span><span class="mord mathnormal">t</span><span class="mord mathnormal" style="margin-right:0.03588em;">y</span><span class="mord mathnormal">a</span><span class="mord mathnormal">x</span><span class="mord mathnormal">i</span><span class="mord mathnormal" style="margin-right:0.05764em;">S</span><span class="mord mathnormal" style="margin-right:0.03588em;">y</span><span class="mord mathnormal">m</span><span class="mord mathnormal" style="margin-right:0.07153em;">C</span><span class="mord mathnormal">a</span><span class="mord mathnormal" style="margin-right:0.03588em;">v</span><span class="mord mathnormal">i</span><span class="mord mathnormal">t</span><span class="mord mathnormal" style="margin-right:0.03588em;">y</span><span class="mclose">?</span><span class="mspace" style="margin-right:0.2778em;"></span><span class="mrel">&gt;</span><span class="mspace" style="margin-right:0.2778em;"></span></span><span class="base"><span class="strut" style="height:0.8889em;vertical-align:-0.1944em;"></span><span class="mord mathnormal">c</span><span class="mord mathnormal">d</span><span class="mord mathnormal">a</span><span class="mord mathnormal">x</span><span class="mord mathnormal">i</span><span class="mord mathnormal" style="margin-right:0.05764em;">S</span><span class="mord mathnormal" style="margin-right:0.03588em;">y</span><span class="mord mathnormal">m</span><span class="mord mathnormal" style="margin-right:0.07153em;">C</span><span class="mord mathnormal">a</span><span class="mord mathnormal" style="margin-right:0.03588em;">v</span><span class="mord mathnormal">i</span><span class="mord mathnormal">t</span><span class="mord mathnormal" style="margin-right:0.03588em;">y</span><span class="mclose">?</span><span class="mspace" style="margin-right:0.2778em;"></span><span class="mrel">&gt;</span><span class="mspace" style="margin-right:0.2778em;"></span></span><span class="base"><span class="strut" style="height:0.9805em;vertical-align:-0.2861em;"></span><span class="mord mathnormal">b</span><span class="mord mathnormal" style="margin-right:0.01968em;">l</span><span class="mord mathnormal">oc</span><span class="mord mathnormal" style="margin-right:0.03148em;">k</span><span class="mord mathnormal" style="margin-right:0.10903em;">M</span><span class="mord mathnormal">es</span><span class="mord mathnormal">h</span><span class="mord">‘‘‘</span><span class="mord cjk_fallback">对于轴对称网格，</span><span class="mord mathnormal">m</span><span class="mord mathnormal">o</span><span class="mord mathnormal" style="margin-right:0.03588em;">v</span><span class="mord mathnormal">in</span><span class="mord mathnormal" style="margin-right:0.03588em;">g</span><span class="mord mathnormal">Wa</span><span class="mord mathnormal" style="margin-right:0.01968em;">ll</span><span class="mord cjk_fallback">面片被用作对称轴（参见图</span><span class="mord">2.24</span><span class="mord cjk_fallback">）。此外，单个</span><span class="mord mathnormal" style="margin-right:0.10764em;">f</span><span class="mord mathnormal">ro</span><span class="mord mathnormal">n</span><span class="mord mathnormal">t</span><span class="mord mathnormal">A</span><span class="mord mathnormal">n</span><span class="mord mathnormal">d</span><span class="mord mathnormal" style="margin-right:0.05017em;">B</span><span class="mord mathnormal">a</span><span class="mord mathnormal">c</span><span class="mord mathnormal" style="margin-right:0.03148em;">k</span><span class="mord cjk_fallback">面片将被分割，并充当楔形体的两个边界（</span><span class="mord mathnormal" style="margin-right:0.10764em;">f</span><span class="mord mathnormal">ro</span><span class="mord mathnormal">n</span><span class="mord mathnormal">t</span><span class="mord mathnormal">A</span><span class="mord mathnormal">n</span><span class="mord mathnormal">d</span><span class="mord mathnormal" style="margin-right:0.05017em;">B</span><span class="mord mathnormal">a</span><span class="mord mathnormal">c</span><span class="mord"><span class="mord mathnormal" style="margin-right:0.03148em;">k</span><span class="msupsub"><span class="vlist-t vlist-t2"><span class="vlist-r"><span class="vlist" style="height:0.1514em;"><span style="top:-2.55em;margin-left:-0.0315em;margin-right:0.05em;"><span class="pstrut" style="height:2.7em;"></span><span class="sizing reset-size6 size3 mtight"><span class="mord mathnormal mtight">n</span></span></span></span><span class="vlist-s">​</span></span><span class="vlist-r"><span class="vlist" style="height:0.15em;"><span></span></span></span></span></span></span><span class="mord mathnormal">e</span><span class="mord mathnormal" style="margin-right:0.03588em;">g</span><span class="mord cjk_fallback">与</span><span class="mord mathnormal" style="margin-right:0.10764em;">f</span><span class="mord mathnormal">ro</span><span class="mord mathnormal">n</span><span class="mord mathnormal">t</span><span class="mord mathnormal">A</span><span class="mord mathnormal">n</span><span class="mord mathnormal">d</span><span class="mord mathnormal" style="margin-right:0.05017em;">B</span><span class="mord mathnormal">a</span><span class="mord mathnormal">c</span><span class="mord"><span class="mord mathnormal" style="margin-right:0.03148em;">k</span><span class="msupsub"><span class="vlist-t vlist-t2"><span class="vlist-r"><span class="vlist" style="height:0.1514em;"><span style="top:-2.55em;margin-left:-0.0315em;margin-right:0.05em;"><span class="pstrut" style="height:2.7em;"></span><span class="sizing reset-size6 size3 mtight"><span class="mord mathnormal mtight">p</span></span></span></span><span class="vlist-s">​</span></span><span class="vlist-r"><span class="vlist" style="height:0.2861em;"><span></span></span></span></span></span></span><span class="mord mathnormal">os</span><span class="mord cjk_fallback">）。在命令行中输入的参数反映了这一点：</span><span class="mord">‘‘‘</span><span class="mord mathnormal">ba</span><span class="mord mathnormal">s</span><span class="mord mathnormal">h</span><span class="mclose">?</span><span class="mspace" style="margin-right:0.2778em;"></span><span class="mrel">&gt;</span><span class="mspace" style="margin-right:0.2778em;"></span></span><span class="base"><span class="strut" style="height:0.7778em;vertical-align:-0.0833em;"></span><span class="mord mathnormal" style="margin-right:0.03148em;">mak</span><span class="mord mathnormal">e</span><span class="mord mathnormal">A</span><span class="mord mathnormal">x</span><span class="mord mathnormal">ia</span><span class="mord mathnormal" style="margin-right:0.10903em;">lM</span><span class="mord mathnormal">es</span><span class="mord mathnormal">h</span><span class="mspace" style="margin-right:0.2222em;"></span><span class="mbin">−</span><span class="mspace" style="margin-right:0.2222em;"></span></span><span class="base"><span class="strut" style="height:0.8889em;vertical-align:-0.1944em;"></span><span class="mord mathnormal">a</span><span class="mord mathnormal">x</span><span class="mord mathnormal">i</span><span class="mord mathnormal">s</span><span class="mord mathnormal">m</span><span class="mord mathnormal">o</span><span class="mord mathnormal" style="margin-right:0.03588em;">v</span><span class="mord mathnormal">in</span><span class="mord mathnormal" style="margin-right:0.03588em;">g</span><span class="mord mathnormal">Wa</span><span class="mord mathnormal" style="margin-right:0.01968em;">ll</span><span class="mspace" style="margin-right:0.2222em;"></span><span class="mbin">−</span><span class="mspace" style="margin-right:0.2222em;"></span></span><span class="base"><span class="strut" style="height:0.8889em;vertical-align:-0.1944em;"></span><span class="mord mathnormal" style="margin-right:0.02691em;">w</span><span class="mord mathnormal">e</span><span class="mord mathnormal">d</span><span class="mord mathnormal" style="margin-right:0.03588em;">g</span><span class="mord mathnormal">e</span><span class="mord mathnormal" style="margin-right:0.10764em;">f</span><span class="mord mathnormal">ro</span><span class="mord mathnormal">n</span><span class="mord mathnormal">t</span><span class="mord mathnormal">A</span><span class="mord mathnormal">n</span><span class="mord mathnormal">d</span><span class="mord mathnormal" style="margin-right:0.05017em;">B</span><span class="mord mathnormal">a</span><span class="mord mathnormal">c</span><span class="mord mathnormal" style="margin-right:0.03148em;">k</span><span class="mord">‘‘‘</span><span class="mord cjk_fallback">该工具将创建一个新的时间目录（在本例中为</span><span class="mord">0.005</span><span class="mord cjk_fallback">）来存储转换后的网格。如果创建未按预期工作，则只需删除此目录，并再次恢复基本网格。案例目录现在应包含如下所示的文件夹：</span><span class="mord">‘‘‘</span><span class="mord mathnormal">ba</span><span class="mord mathnormal">s</span><span class="mord mathnormal">h</span><span class="mclose">?</span><span class="mspace" style="margin-right:0.2778em;"></span><span class="mrel">&gt;</span><span class="mspace" style="margin-right:0.2778em;"></span></span><span class="base"><span class="strut" style="height:0.8889em;vertical-align:-0.1944em;"></span><span class="mord mathnormal" style="margin-right:0.01968em;">l</span><span class="mord mathnormal">s</span><span class="mord">00.005</span><span class="mord mathnormal">co</span><span class="mord mathnormal">n</span><span class="mord mathnormal">s</span><span class="mord mathnormal">t</span><span class="mord mathnormal">an</span><span class="mord mathnormal">t</span><span class="mord mathnormal">sys</span><span class="mord mathnormal">t</span><span class="mord mathnormal">e</span><span class="mord mathnormal">m</span><span class="mord">‘‘‘</span><span class="mord cjk_fallback">此时，网格已弯曲成</span><span class="mord">5°</span><span class="mord cjk_fallback">楔形，如图</span><span class="mord">2.25</span><span class="mord cjk_fallback">所示。但是，</span><span class="mord mathnormal">m</span><span class="mord mathnormal">o</span><span class="mord mathnormal" style="margin-right:0.03588em;">v</span><span class="mord mathnormal">in</span><span class="mord mathnormal" style="margin-right:0.03588em;">g</span><span class="mord mathnormal">Wa</span><span class="mord mathnormal" style="margin-right:0.01968em;">ll</span><span class="mord cjk_fallback">面片中的面仍然存在，但它们现在被挤压为面面积接近于零的面。</span><span class="mord mathnormal" style="margin-right:0.03148em;">mak</span><span class="mord mathnormal">e</span><span class="mord mathnormal">A</span><span class="mord mathnormal">x</span><span class="mord mathnormal">ia</span><span class="mord mathnormal" style="margin-right:0.10903em;">lM</span><span class="mord mathnormal">es</span><span class="mord mathnormal">h</span><span class="mord cjk_fallback">变换点位置，但不改变网格连接。因此，对称面片没有指定面（</span><span class="mord mathnormal">n</span><span class="mord mathnormal" style="margin-right:0.13889em;">F</span><span class="mord mathnormal">a</span><span class="mord mathnormal">ces</span><span class="mspace" style="margin-right:0.2778em;"></span><span class="mrel">=</span><span class="mspace" style="margin-right:0.2778em;"></span></span><span class="base"><span class="strut" style="height:1em;vertical-align:-0.25em;"></span><span class="mord">0</span><span class="mord cjk_fallback">），必须移除。</span><span class="mclose">!</span><span class="mopen">[</span><span class="mord cjk_fallback">图</span><span class="mord">2.25</span><span class="mord cjk_fallback">：</span><span class="mord mathnormal" style="margin-right:0.03148em;">mak</span><span class="mord mathnormal">e</span><span class="mord mathnormal">A</span><span class="mord mathnormal">x</span><span class="mord mathnormal">ia</span><span class="mord mathnormal" style="margin-right:0.10903em;">lM</span><span class="mord mathnormal">es</span><span class="mord mathnormal">h</span><span class="mord cjk_fallback">楔形块变换前后的</span><span class="mord">2</span><span class="mord mathnormal" style="margin-right:0.02778em;">D</span><span class="mord cjk_fallback">型腔网格</span><span class="mclose">]</span><span class="mopen">(</span><span class="mord mathnormal">h</span><span class="mord mathnormal">ttp</span><span class="mord mathnormal">s</span><span class="mspace" style="margin-right:0.2778em;"></span><span class="mrel">:</span><span class="mspace" style="margin-right:0.2778em;"></span></span><span class="base"><span class="strut" style="height:1em;vertical-align:-0.25em;"></span><span class="mord">//</span><span class="mord mathnormal" style="margin-right:0.10764em;">f</span><span class="mord mathnormal">a</span><span class="mord mathnormal">ee</span><span class="mord">0</span><span class="mspace" style="margin-right:0.2222em;"></span><span class="mbin">−</span><span class="mspace" style="margin-right:0.2222em;"></span></span><span class="base"><span class="strut" style="height:0.8389em;vertical-align:-0.1944em;"></span><span class="mord">1253397841.</span><span class="mord mathnormal">cos</span><span class="mord">.</span><span class="mord mathnormal">a</span><span class="mord mathnormal">p</span><span class="mspace" style="margin-right:0.2222em;"></span><span class="mbin">−</span><span class="mspace" style="margin-right:0.2222em;"></span></span><span class="base"><span class="strut" style="height:1em;vertical-align:-0.25em;"></span><span class="mord mathnormal">c</span><span class="mord mathnormal">h</span><span class="mord mathnormal">e</span><span class="mord mathnormal">n</span><span class="mord mathnormal" style="margin-right:0.03588em;">g</span><span class="mord mathnormal">d</span><span class="mord mathnormal">u</span><span class="mord">.</span><span class="mord mathnormal">m</span><span class="mord mathnormal" style="margin-right:0.03588em;">y</span><span class="mord mathnormal" style="margin-right:0.03588em;">q</span><span class="mord mathnormal">c</span><span class="mord mathnormal" style="margin-right:0.01968em;">l</span><span class="mord mathnormal">o</span><span class="mord mathnormal">u</span><span class="mord mathnormal">d</span><span class="mord">.</span><span class="mord mathnormal">co</span><span class="mord mathnormal">m</span><span class="mord">/</span><span class="mord mathnormal" style="margin-right:0.02691em;">w</span><span class="mord mathnormal">i</span><span class="mord mathnormal" style="margin-right:0.04398em;">z</span><span class="mord mathnormal">im</span><span class="mord mathnormal" style="margin-right:0.03588em;">g</span><span class="mord">/202212030308697.</span><span class="mord mathnormal">p</span><span class="mord mathnormal">n</span><span class="mord mathnormal" style="margin-right:0.03588em;">g</span><span class="mclose">)</span><span class="mord cjk_fallback">在这种情况下，建议使用</span><span class="mord mathnormal">co</span><span class="mord mathnormal" style="margin-right:0.01968em;">ll</span><span class="mord mathnormal">a</span><span class="mord mathnormal">p</span><span class="mord mathnormal">se</span><span class="mord mathnormal" style="margin-right:0.05764em;">E</span><span class="mord mathnormal">d</span><span class="mord mathnormal" style="margin-right:0.03588em;">g</span><span class="mord mathnormal">es</span><span class="mord cjk_fallback">工具。它需要两个必需的命令行参数：边长和合并角：</span><span class="mord">‘‘‘</span><span class="mord mathnormal">ba</span><span class="mord mathnormal">s</span><span class="mord mathnormal">h</span><span class="mclose">?</span><span class="mspace" style="margin-right:0.2778em;"></span><span class="mrel">&gt;</span><span class="mspace" style="margin-right:0.2778em;"></span></span><span class="base"><span class="strut" style="height:0.8889em;vertical-align:-0.1944em;"></span><span class="mord mathnormal">co</span><span class="mord mathnormal" style="margin-right:0.01968em;">ll</span><span class="mord mathnormal">a</span><span class="mord mathnormal">p</span><span class="mord mathnormal">se</span><span class="mord mathnormal" style="margin-right:0.05764em;">E</span><span class="mord mathnormal">d</span><span class="mord mathnormal" style="margin-right:0.03588em;">g</span><span class="mord mathnormal">es</span><span class="mspace" style="margin-right:0.2778em;"></span><span class="mrel">&lt;</span><span class="mspace" style="margin-right:0.2778em;"></span></span><span class="base"><span class="strut" style="height:1em;vertical-align:-0.25em;"></span><span class="mord mathnormal">e</span><span class="mord mathnormal">d</span><span class="mord mathnormal" style="margin-right:0.03588em;">g</span><span class="mord mathnormal">e</span><span class="mord mathnormal" style="margin-right:0.01968em;">l</span><span class="mord mathnormal">e</span><span class="mord mathnormal">n</span><span class="mord mathnormal" style="margin-right:0.03588em;">g</span><span class="mord mathnormal">t</span><span class="mord mathnormal">h</span><span class="mopen">[</span><span class="mord mathnormal">m</span><span class="mclose">]</span><span class="mspace" style="margin-right:0.2778em;"></span><span class="mrel">&gt;&lt;</span><span class="mspace" style="margin-right:0.2778em;"></span></span><span class="base"><span class="strut" style="height:1em;vertical-align:-0.25em;"></span><span class="mord mathnormal">m</span><span class="mord mathnormal" style="margin-right:0.02778em;">er</span><span class="mord mathnormal" style="margin-right:0.03588em;">g</span><span class="mord mathnormal">e</span><span class="mord mathnormal">an</span><span class="mord mathnormal" style="margin-right:0.03588em;">g</span><span class="mord mathnormal" style="margin-right:0.01968em;">l</span><span class="mord mathnormal">e</span><span class="mopen">[</span><span class="mord mathnormal">d</span><span class="mord mathnormal">e</span><span class="mord mathnormal" style="margin-right:0.03588em;">g</span><span class="mord mathnormal">rees</span><span class="mclose">]</span><span class="mspace" style="margin-right:0.2778em;"></span><span class="mrel">&gt;</span><span class="mspace" style="margin-right:0.2778em;"></span></span><span class="base"><span class="strut" style="height:0.6944em;"></span><span class="mord">‘‘‘</span><span class="mord cjk_fallback">在许多应用中，边长为</span></span></span></span>1 × 10^{−8}$米，合并角为179°时，可以正确识别和移除最近塌陷的面。在网格边长比例极小的某些情况下，可能需要较小的边长以避免误报和无意中删除有效边。对于本示例，使用所示参数执行collapseEdges不会出现问题。

```bash
?>  collapseEdges -latestTime 1e-8 179

对于一些最终的内务处理，建议从边界列表中删除现在为空的面片。打开constant/polyMesh/boundary并删除movingWall与frontAndBack条目。请注意，它们被列为包含零个面：nFaces 0;。将边界列表大小变更为3，以反映这两项删除。边界文件现在应类似于：

3
(
    fixedWalls
    {
        type            wall;
        nFaces            60;
        startFace        760;
    }
    frontAndBack_pos
    {
        type            wedge;
        nFaces           400;
        startFace        820;
    }
    frontAndBack_neg
    {
        type            wedge;
        nFaces          400;
        startFace       1220;
    }
)

此时，可以使用autoPatch工具将fixedWalls面片分割为3个单独的面片。这将查看连续的面片，并尝试根据给定的特征角度确定适当的位置来分割它。

在这种情况下，任何形成大于30°角的patch边都可以被分割以形成新的patch。在分配边界条件时，这提供了更大的灵活性。

?>  autoPatch -latestTime 30

分割后将重命名面片。-latestTime标志将仅读取可用的最新时间步。分割网格面会储存在另一个时间步长目录中，而不是覆写时间步长。最后，应使用checkMesh工具检查网格是否存在错误，这应被视为最佳实践的一般规则：更改网格时始终运行checkMesh。

可以在目录WM_PROJECT_DIR/applications/utilities/mesh中找到处理网格操作的实用程序应用程序（或只是简短的实用程序）。网格实用程序分为以下类别：生成、操纵、推进和转换。这种分类在最新版本中完全没有改变。生成网格并将其从不同格式转换为OpenFOAM格式已在第2.2节和第2.3节中描述。本节介绍如何在生成基础网格后操纵网格以及诸如网格细化等高级操作。

2.4.1 按指定标准细化网格

在此示例中，网格细化应用程序refineHexMesh用于细化interFoam解算器的damBreak教程的网格。这样做的目的是细化初始自由表面周围的区域，其中两相标志场（\alpha_{water}$）的梯度小于零（$\nabla(\alpha_{water}$<0）。首先，必须在您选择的工作目录中生成damBreak教程的本地副本。必须执行由Allrun脚本执行的所有实用程序，只有求解器未启动。这将生成网格并初始化$\alpha_{water}$场。‘‘‘bash?>cp-r$FOAM_TUTORIALS/multiphase/interFoam/laminar/damBreak . ?> cd damBreak ?> blockMesh ?> setFields

现在用blockMesh生成网格，并使用setFields前处理实用程序设置$\alpha_{water}$场。 setfields实用工具在3.2节中进行了描述。 基本的计算工具foamCalc可以用来计算和存储$\alpha_{water}$的梯度。

```bash
?>  foamCalc magGrad alpha.water

这将把梯度幅度的以单元为中心的标量场存储在名为magGradalphaWater.的初始时间目录0中。 要使用refineMesh应用程序根据梯度大小细化网格，必须将该实用程序的配置字典文件复制到dambreake案例的系统目录中。

?>  cp <span class="katex"><span class="katex-html" aria-hidden="true"><span class="base"><span class="strut" style="height:1em;vertical-align:-0.25em;"></span><span class="mord mathnormal" style="margin-right:0.02778em;">FO</span><span class="mord mathnormal">A</span><span class="mord"><span class="mord mathnormal" style="margin-right:0.10903em;">M</span><span class="msupsub"><span class="vlist-t vlist-t2"><span class="vlist-r"><span class="vlist" style="height:0.3283em;"><span style="top:-2.55em;margin-left:-0.109em;margin-right:0.05em;"><span class="pstrut" style="height:2.7em;"></span><span class="sizing reset-size6 size3 mtight"><span class="mord mathnormal mtight">A</span></span></span></span><span class="vlist-s">​</span></span><span class="vlist-r"><span class="vlist" style="height:0.15em;"><span></span></span></span></span></span></span><span class="mord mathnormal" style="margin-right:0.13889em;">PP</span><span class="mord">/</span><span class="mord mathnormal">u</span><span class="mord mathnormal">t</span><span class="mord mathnormal">i</span><span class="mord mathnormal" style="margin-right:0.01968em;">l</span><span class="mord mathnormal">i</span><span class="mord mathnormal">t</span><span class="mord mathnormal">i</span><span class="mord mathnormal">es</span><span class="mord">/</span><span class="mord mathnormal">m</span><span class="mord mathnormal">es</span><span class="mord mathnormal">h</span><span class="mord">/</span><span class="mord mathnormal">mani</span><span class="mord mathnormal">p</span><span class="mord mathnormal">u</span><span class="mord mathnormal" style="margin-right:0.01968em;">l</span><span class="mord mathnormal">a</span><span class="mord mathnormal">t</span><span class="mord mathnormal">i</span><span class="mord mathnormal">o</span><span class="mord mathnormal">n</span><span class="mord">/</span><span class="mord mathnormal">re</span><span class="mord mathnormal" style="margin-right:0.10764em;">f</span><span class="mord mathnormal">in</span><span class="mord mathnormal">e</span><span class="mord mathnormal" style="margin-right:0.10903em;">M</span><span class="mord mathnormal">es</span><span class="mord mathnormal">h</span><span class="mord">/</span><span class="mord mathnormal">re</span><span class="mord mathnormal" style="margin-right:0.10764em;">f</span><span class="mord mathnormal">in</span><span class="mord mathnormal">e</span><span class="mord mathnormal" style="margin-right:0.10903em;">M</span><span class="mord mathnormal">es</span><span class="mord mathnormal">h</span><span class="mord mathnormal" style="margin-right:0.02778em;">D</span><span class="mord mathnormal">i</span><span class="mord mathnormal">c</span><span class="mord mathnormal">t</span><span class="mord mathnormal">sys</span><span class="mord mathnormal">t</span><span class="mord mathnormal">e</span><span class="mord mathnormal">m</span><span class="mord">/</span><span class="mclose">?</span><span class="mspace" style="margin-right:0.2778em;"></span><span class="mrel">&gt;</span><span class="mspace" style="margin-right:0.2778em;"></span></span><span class="base"><span class="strut" style="height:1em;vertical-align:-0.25em;"></span><span class="mord mathnormal" style="margin-right:0.01968em;">l</span><span class="mord mathnormal">ssys</span><span class="mord mathnormal">t</span><span class="mord mathnormal">e</span><span class="mord mathnormal">m</span><span class="mord">/</span><span class="mord mathnormal">co</span><span class="mord mathnormal">n</span><span class="mord mathnormal">t</span><span class="mord mathnormal">ro</span><span class="mord mathnormal" style="margin-right:0.01968em;">l</span><span class="mord mathnormal" style="margin-right:0.02778em;">D</span><span class="mord mathnormal">i</span><span class="mord mathnormal">c</span><span class="mord mathnormal">t</span><span class="mord mathnormal" style="margin-right:0.10764em;">f</span><span class="mord mathnormal" style="margin-right:0.03588em;">v</span><span class="mord mathnormal" style="margin-right:0.05764em;">S</span><span class="mord mathnormal">c</span><span class="mord mathnormal">h</span><span class="mord mathnormal">e</span><span class="mord mathnormal">m</span><span class="mord mathnormal">esre</span><span class="mord mathnormal" style="margin-right:0.10764em;">f</span><span class="mord mathnormal">in</span><span class="mord mathnormal">e</span><span class="mord mathnormal" style="margin-right:0.10903em;">M</span><span class="mord mathnormal">es</span><span class="mord mathnormal">h</span><span class="mord mathnormal" style="margin-right:0.02778em;">D</span><span class="mord mathnormal">i</span><span class="mord mathnormal">c</span><span class="mord mathnormal">t</span><span class="mord mathnormal">d</span><span class="mord mathnormal">eco</span><span class="mord mathnormal">m</span><span class="mord mathnormal">p</span><span class="mord mathnormal">ose</span><span class="mord mathnormal" style="margin-right:0.13889em;">P</span><span class="mord mathnormal">a</span><span class="mord mathnormal" style="margin-right:0.02778em;">rD</span><span class="mord mathnormal">i</span><span class="mord mathnormal">c</span><span class="mord mathnormal">t</span><span class="mord mathnormal" style="margin-right:0.10764em;">f</span><span class="mord mathnormal" style="margin-right:0.03588em;">v</span><span class="mord mathnormal" style="margin-right:0.05764em;">S</span><span class="mord mathnormal">o</span><span class="mord mathnormal" style="margin-right:0.01968em;">l</span><span class="mord mathnormal">u</span><span class="mord mathnormal">t</span><span class="mord mathnormal">i</span><span class="mord mathnormal">o</span><span class="mord mathnormal">n</span><span class="mord mathnormal">se</span><span class="mord mathnormal" style="margin-right:0.13889em;">tF</span><span class="mord mathnormal">i</span><span class="mord mathnormal">e</span><span class="mord mathnormal" style="margin-right:0.01968em;">l</span><span class="mord mathnormal">d</span><span class="mord mathnormal" style="margin-right:0.02778em;">sD</span><span class="mord mathnormal">i</span><span class="mord mathnormal">c</span><span class="mord mathnormal">t</span><span class="mord">‘‘‘</span><span class="mord cjk_fallback">在</span><span class="mord mathnormal">sys</span><span class="mord mathnormal">t</span><span class="mord mathnormal">e</span><span class="mord mathnormal">m</span><span class="mord cjk_fallback">目录中的配置中，</span><span class="mord mathnormal">re</span><span class="mord mathnormal" style="margin-right:0.10764em;">f</span><span class="mord mathnormal">in</span><span class="mord mathnormal">eHe</span><span class="mord mathnormal">x</span><span class="mord mathnormal" style="margin-right:0.10903em;">M</span><span class="mord mathnormal">es</span><span class="mord mathnormal">h</span><span class="mord cjk_fallback">对某个网格集中的所有网格进行细化。</span><span class="mord">‘‘‘</span><span class="mord mathnormal">c</span><span class="mspace" style="margin-right:0.2222em;"></span><span class="mbin">+</span><span class="mspace" style="margin-right:0.2222em;"></span></span><span class="base"><span class="strut" style="height:1em;vertical-align:-0.25em;"></span><span class="mord">+</span><span class="mord">//</span><span class="mord mathnormal" style="margin-right:0.07153em;">C</span><span class="mord mathnormal">e</span><span class="mord mathnormal" style="margin-right:0.01968em;">ll</span><span class="mord mathnormal">s</span><span class="mord mathnormal">t</span><span class="mord mathnormal">ore</span><span class="mord mathnormal" style="margin-right:0.10764em;">f</span><span class="mord mathnormal">in</span><span class="mord mathnormal">e</span><span class="mpunct">;</span><span class="mspace" style="margin-right:0.1667em;"></span><span class="mord mathnormal">nam</span><span class="mord mathnormal">eo</span><span class="mord mathnormal" style="margin-right:0.10764em;">f</span><span class="mord mathnormal">ce</span><span class="mord mathnormal" style="margin-right:0.01968em;">ll</span><span class="mord mathnormal">se</span><span class="mord mathnormal">t</span><span class="mord mathnormal">se</span><span class="mord mathnormal">t</span><span class="mord mathnormal">c</span><span class="mord">0</span><span class="mpunct">;</span><span class="mspace" style="margin-right:0.1667em;"></span><span class="mord">‘‘‘</span><span class="mord cjk_fallback">创建此</span><span class="mord mathnormal" style="margin-right:0.07153em;">C</span><span class="mord mathnormal">e</span><span class="mord mathnormal" style="margin-right:0.05764em;">llS</span><span class="mord mathnormal">e</span><span class="mord mathnormal">t</span><span class="mord cjk_fallback">时，将存储在</span><span class="mord mathnormal">co</span><span class="mord mathnormal">n</span><span class="mord mathnormal">s</span><span class="mord mathnormal">t</span><span class="mord mathnormal">an</span><span class="mord mathnormal">t</span><span class="mord">/</span><span class="mord mathnormal">p</span><span class="mord mathnormal">o</span><span class="mord mathnormal" style="margin-right:0.01968em;">l</span><span class="mord mathnormal" style="margin-right:0.03588em;">y</span><span class="mord mathnormal">m</span><span class="mord mathnormal">es</span><span class="mord mathnormal">h</span><span class="mord cjk_fallback">中，</span><span class="mord mathnormal" style="margin-right:0.00773em;">R</span><span class="mord mathnormal">e</span><span class="mord mathnormal" style="margin-right:0.10764em;">f</span><span class="mord mathnormal">in</span><span class="mord mathnormal">eHe</span><span class="mord mathnormal">x</span><span class="mord mathnormal">m</span><span class="mord mathnormal">es</span><span class="mord mathnormal">h</span><span class="mord cjk_fallback">将使用它来细化单元格。在本例中，</span><span class="mord mathnormal">t</span><span class="mord mathnormal">o</span><span class="mord mathnormal">p</span><span class="mord mathnormal">o</span><span class="mord mathnormal" style="margin-right:0.05764em;">S</span><span class="mord mathnormal">e</span><span class="mord mathnormal">t</span><span class="mord cjk_fallback">用于生成</span><span class="mord mathnormal">ce</span><span class="mord mathnormal" style="margin-right:0.05764em;">llS</span><span class="mord mathnormal">e</span><span class="mord mathnormal">t</span><span class="mord cjk_fallback">。这要求</span><span class="mord mathnormal">t</span><span class="mord mathnormal">o</span><span class="mord mathnormal">p</span><span class="mord mathnormal">o</span><span class="mord mathnormal" style="margin-right:0.05764em;">S</span><span class="mord mathnormal">e</span><span class="mord mathnormal">t</span><span class="mord mathnormal" style="margin-right:0.02778em;">D</span><span class="mord mathnormal">i</span><span class="mord mathnormal">c</span><span class="mord mathnormal">t</span><span class="mord cjk_fallback">存在于系统目录中并正确配置。因此，随后会复制并更改现有的文件。</span><span class="mord">‘‘‘</span><span class="mord mathnormal">ba</span><span class="mord mathnormal">s</span><span class="mord mathnormal">h</span><span class="mclose">?</span><span class="mspace" style="margin-right:0.2778em;"></span><span class="mrel">&gt;</span><span class="mspace" style="margin-right:0.2778em;"></span></span><span class="base"><span class="strut" style="height:0.625em;vertical-align:-0.1944em;"></span><span class="mord mathnormal">c</span><span class="mord mathnormal">p</span></span></span></span>FOAM_APP/utilities/mesh/manipulation/topoSet/topoSetDict system/

Toposetdict的示例操作子字典必须由以下内容替换：

actions
(
    {
        name c0;
        type cellSet;
        action new;
        source fieldToCell;
        sourceInfo
        {
            fieldName magGradalpha1;
            min 20;
            max 100;
        }
    }
);

现在可以根据System/Toposetdict中的定义生成cellset：

?>  topoSet
?>  refineHexMesh c0

当网格现在使用Paraview查看时，自由表面的区域现在应该有额外的分辨率。

2.4.2 变换点

在OpenFOAM网格格式中，与网格的比例和位置有关的唯一信息位于点位置向量中。如前所述，所有剩余的存储的网格信息是纯粹基于连通性的。也就是说，网格大小、位置和方向可以通过单独变换点位置来改变。为此，transformPoints网格实用程序随OpenFOAM提供。由于此实用程序相对简单，因此只显示了所需的语法。变换网格时最常用的选项是-rotate、-translate和-scale选项。

也可以在括号周围使用双引号而不是单引号。 执行任务的顺序是硬编码的，用户不能更改。 如果要确保在转换之前执行缩放，请运行TransformPoints两次。

• scale。按指定的标量量在任意或所有基本方向上缩放网格的点。-scale '（1.0 1.0 1.0）'不会变更点位置，而-scale '（2.0 2.0 2.0）'会在所有方向上将点位置均匀加倍。任何非均匀缩放都将沿给定方向拉伸或压缩网格。
• translate。按给定向量移动网格，有效地将此向量添加到网格中的每个点位置向量。
• rotate。旋转网格。旋转由输入向量定义。网格将经历使用第二个向量定向第一个向量所需的旋转。旋转网格时，也可以通过添加-rotateFields选项来旋转任何初始或边界向量与张量值。

这三个点转换的语法如下所示。

?>  transformPoints -scale '(x y z)'
?>  transformPoints -translate '(x y z)'
?>  transformPoints -rotateFields -rotate '( (x0 y0 z0) (x1 y1 z1) )'

2.4.3 镜像网格

有时，用对称平面生成网格，然后执行镜面反射比在一个步骤中网格化整个几何结构更容易。 mirrorMesh实用程序正是这样做的。 关于镜像进程本身的所有参数都是从位于system/mirrorMeshDict中的字典中读取的。

为了成功镜像网格，镜像平面必须是平面。在本例中，四分之一网格被镜像到全域中。首先，下面的固析例必须复制到选择的目录中并重新命名，以防以后混淆。mirrorMeshDict需要从现有案例复制到案例系统目录中。

?>  cp -r <span class="katex"><span class="katex-html" aria-hidden="true"><span class="base"><span class="strut" style="height:1em;vertical-align:-0.25em;"></span><span class="mord mathnormal" style="margin-right:0.02778em;">FO</span><span class="mord mathnormal">A</span><span class="mord"><span class="mord mathnormal" style="margin-right:0.10903em;">M</span><span class="msupsub"><span class="vlist-t vlist-t2"><span class="vlist-r"><span class="vlist" style="height:0.3283em;"><span style="top:-2.55em;margin-left:-0.109em;margin-right:0.05em;"><span class="pstrut" style="height:2.7em;"></span><span class="sizing reset-size6 size3 mtight"><span class="mord mathnormal mtight" style="margin-right:0.13889em;">T</span></span></span></span><span class="vlist-s">​</span></span><span class="vlist-r"><span class="vlist" style="height:0.15em;"><span></span></span></span></span></span></span><span class="mord mathnormal" style="margin-right:0.10903em;">U</span><span class="mord mathnormal" style="margin-right:0.00773em;">TOR</span><span class="mord mathnormal" style="margin-right:0.07847em;">I</span><span class="mord mathnormal">A</span><span class="mord mathnormal">L</span><span class="mord mathnormal" style="margin-right:0.05764em;">S</span><span class="mord">/</span><span class="mord mathnormal">s</span><span class="mord mathnormal">t</span><span class="mord mathnormal">ress</span><span class="mord mathnormal">A</span><span class="mord mathnormal">na</span><span class="mord mathnormal" style="margin-right:0.01968em;">l</span><span class="mord mathnormal">ys</span><span class="mord mathnormal">i</span><span class="mord mathnormal">s</span><span class="mord">/</span><span class="mord mathnormal">so</span><span class="mord mathnormal" style="margin-right:0.01968em;">l</span><span class="mord mathnormal">i</span><span class="mord mathnormal">d</span><span class="mord mathnormal" style="margin-right:0.02778em;">D</span><span class="mord mathnormal">i</span><span class="mord mathnormal">s</span><span class="mord mathnormal" style="margin-right:0.01968em;">pl</span><span class="mord mathnormal">a</span><span class="mord mathnormal">ce</span><span class="mord mathnormal">m</span><span class="mord mathnormal">e</span><span class="mord mathnormal">n</span><span class="mord mathnormal" style="margin-right:0.13889em;">tF</span><span class="mord mathnormal">o</span><span class="mord mathnormal">am</span><span class="mord">/</span><span class="mord mathnormal" style="margin-right:0.01968em;">pl</span><span class="mord mathnormal">a</span><span class="mord mathnormal">t</span><span class="mord mathnormal">eHo</span><span class="mord mathnormal" style="margin-right:0.01968em;">l</span><span class="mord mathnormal">e</span><span class="mord">.</span><span class="mclose">?</span><span class="mspace" style="margin-right:0.2778em;"></span><span class="mrel">&gt;</span><span class="mspace" style="margin-right:0.2778em;"></span></span><span class="base"><span class="strut" style="height:0.8889em;vertical-align:-0.1944em;"></span><span class="mord mathnormal">m</span><span class="mord mathnormal" style="margin-right:0.03588em;">v</span><span class="mord mathnormal" style="margin-right:0.01968em;">pl</span><span class="mord mathnormal">a</span><span class="mord mathnormal">t</span><span class="mord mathnormal">eHo</span><span class="mord mathnormal" style="margin-right:0.01968em;">l</span><span class="mord mathnormal">e</span><span class="mord mathnormal">mi</span><span class="mord mathnormal" style="margin-right:0.02778em;">rror</span><span class="mord mathnormal" style="margin-right:0.10903em;">M</span><span class="mord mathnormal">es</span><span class="mord mathnormal">h</span><span class="mord mathnormal" style="margin-right:0.05764em;">E</span><span class="mord mathnormal">x</span><span class="mord mathnormal">am</span><span class="mord mathnormal" style="margin-right:0.01968em;">pl</span><span class="mord mathnormal">e</span><span class="mclose">?</span><span class="mspace" style="margin-right:0.2778em;"></span><span class="mrel">&gt;</span><span class="mspace" style="margin-right:0.2778em;"></span></span><span class="base"><span class="strut" style="height:0.8889em;vertical-align:-0.1944em;"></span><span class="mord mathnormal">c</span><span class="mord mathnormal">d</span><span class="mord mathnormal">mi</span><span class="mord mathnormal" style="margin-right:0.02778em;">rror</span><span class="mord mathnormal" style="margin-right:0.10903em;">M</span><span class="mord mathnormal">es</span><span class="mord mathnormal">h</span><span class="mord mathnormal" style="margin-right:0.05764em;">E</span><span class="mord mathnormal">x</span><span class="mord mathnormal">am</span><span class="mord mathnormal" style="margin-right:0.01968em;">pl</span><span class="mord mathnormal">e</span><span class="mclose">?</span><span class="mspace" style="margin-right:0.2778em;"></span><span class="mrel">&gt;</span><span class="mspace" style="margin-right:0.2778em;"></span></span><span class="base"><span class="strut" style="height:0.7778em;vertical-align:-0.1944em;"></span><span class="mord mathnormal">c</span><span class="mord mathnormal">p</span><span class="mspace" style="margin-right:0.2222em;"></span><span class="mbin">−</span><span class="mspace" style="margin-right:0.2222em;"></span></span><span class="base"><span class="strut" style="height:0.4306em;"></span><span class="mord mathnormal" style="margin-right:0.02778em;">r</span></span></span></span>FOAM_APP/utilities/mesh/manipulation/mirrorMesh/mirrorMeshDict system/

下一步是定义将作为镜像平面的平面。这样的平面可以由原点和法向向量定义，在mirrorMeshDict中如下所示。发生镜像的patch会自动移除。

pointAndNormalDict
{
    basePoint    (0 0 0);
    normalVector (0 -1 0);
}

正确定义后，可以执行mirrorMesh，并检查网格是否存在错误：

?>  mirrorMesh
?>  checkMesh

这会产生半网格。对于第二次镜像，必须更改字典以说明不同的镜像平面：

pointAndNormalDict
{
    basePoint    (0 0 0);
    normalVector (-1 0 0);
}

再次运行mirrorMesh并使用checkMesh检查网格本身是否存在任何错误：

?>  mirrorMesh
?>  checkMesh

这将生成一个全域网格。

有很多开源应用程序可用于设计计算域几何形状以及生成网格。另一方面，也有用于相同目的的商业解决方案。几何图形定义和网格生成有时涉及冗长的工作流，直到几何图形减少到模拟可接受的级别，并且生成的网格具有足够的质量。对于简单的情况，例如验证，blockMesh仍然是许多用户首选的网格生成工具。自动网格生成解决方案（如snappyHexMesh、cfMesh、foamHexMesh和engrid）降低了网格生成的复杂性和持续时间，但仍涉及处理不同的网格生成参数。上述应用程序随其详细文档分发，因此本章旨在提供OpenFOAM中网格的组成方式信息，以及使用blockMesh和snappyHexMesh生成网格的附加信息，以扩展当前可用的文档。

1. Sept. 2014. Url: http://www.sourceflux.de/blog/adding-source-terms-equations-fvoptions/.
2. May 2015. Url: http://www.sourceflux.de/blog/a-changing-cell-set-in-openfoam/.
3. J. Höpken and T. Maric. Feature handling in snappyHexMesh. Url: www.sourceflux.de/blog/snappyhexmesh-snapping-edges.
4. OpenFOAM User Guide. OpenCFD limited. 2016.
5. E. de Villiers. 7 th OpenFOAM Workshop: snappyHexMesh Training. Url: http://www.openfoamworkshop.org/2012/downloads/Training/EugenedeVilliers/EugenedeVilliers-TrainingSlides.tgz (visited on 12/2013).

本章将介绍 OpenFOAM 模拟的结构以及边界条件和初始条件的定义。

注：本章的模拟测试用例可在第 3 章子目录下的示例案例库中找到。

OpenFOAM仿真文件是一个目录，其中包含存储在不同子目录中的一组文件。一些文件用于配置和控制仿真，其他文件用于存储生成的仿真数据。通常，OpenFOAM模拟案例的基于文件的组织相对简单易用：可以使用文本编辑器轻松编辑配置文件。使用基于文件的组织还有一个优点：可以轻松地参数化仿真。

使用icoFoam不可压缩Navier-Stokes求解器使用的方腔顶盖流教程案例解释了OpenFOAM案例(模拟)的标准组成。图3.1显示了型腔设置的示意图。

Case目录中存在icoFoam不需要的其他输入文件，它们由预处理实用程序使用。通常，所需输入字典的数量随着求解器复杂性的增加而增加。

下面的方腔案例的子目录列表显示了模拟文件的组织方式：

?>  cd $FOAM_TUTORIALS/incompressible/icoFoam/cavity
?>  ls *

0:
p U
constant:
polyMesh transportProperties
system:
controlDict fvSchemes fvSolution

0、constant和system目录是OpenFOAM模拟案例必须包含的标准目录。0目录保存应用于物理场的初始条件和边界条件。特定求解器使用的每个物理场都由一个以该物理场命名的文本文件表示。对方腔案例的问题，模拟中使用的场是压力场$p$和速度场$\mathbf{U}$，这些场的定义如第3.2节所示。

随着模拟的进行，求解器应用程序将结果模拟数据写入案例目录的新子目录。这些目录就是所谓的时间步目录，该目录中的文件是根据模拟时间值命名的。它们不仅包含已在0/目录中定义的那些物理场变量，还包含解算器的辅助变量，如体积通量phi。

尽管icoFoam求解器使用压力p场和速度U场的初始值开始模拟，但FVM在方程离散过程中会使用体积通量场(Phi)(有关详细信息，请参见第1章)。因此，时间步长目录将保存由求解器应用程序计算的体积通量场phi。

除了时间点目录之外，还可以写入其他模拟数据：

• 通过求解器应用程序（如体积通量phi）
• 通过通过与求解器一起运行的function object
• 在模拟完成之后，作为后处理的结果

Constant目录存储整个模拟过程中保持不变的模拟数据。这通常包含polyMesh子目录中的网格数据，以及各种配置文件：

• transportProperties。输运数据
• turbulenceProperties。湍流模型数据
• dynamicMeshDict。动网格控制数据

并非所有上述配置文件都存在于方腔实例的constant目录中。缺少的是icoFoam不需要的turbulenceProperties字典。需要提供哪些附加字典，这取决于选定的求解器应用程序。如果在缺少必要输入数据的情况下在案例目录中执行求解器，则会向用户提示错误消息，通知用户缺少或定义了错误的字典、字典参数或物理场变量。

如果求解器应用程序使用了动态网格功能，若网格节点位置或网格拓扑发生更改，则新的polyMesh文件夹也会写入每个时间点目录。第13章讨论了OpenFOAM的动态网格功能。

OpenFOAM中的配置文件通常称为字典文件，甚至简称为字典。这种命名是由于它们在源代码中的使用，源代码基于IODICTIONARY类。

这些方法见第1章。它还可能包含用于配置不同预处理和后处理应用程序（如设置字段）的字典。system目录中最重要的字典是controlDict，其控制与解算器运行时间相关的所有参数，以及将求解数据写入到case目录的频率。controlDict中定义的所有参数与用于模拟的求解器无关。有关使用controlDict控制模拟运行的更多信息，请参见第3.4节。[5]对controlDict中的参数进行了相当广泛的讨论，第3.4.1节对其中一些参数进行了解释。

0OpenFOAM simulation是一个包含用于配置模拟的不同子目录和文件的目录。这种OpenFOAM模拟的文件结构使得设置边界和初始条件非常简单。每个对模拟重要的物理量（压力、温度、速度场等）都有其文件，存储在0目录中：属于模拟第一个时间步的目录。本节介绍如何在这些文件中实际应用第1.3节中的边界条件和初始条件。根据场的张量秩（标量、向量、张量），使用稍微不同的语法设置它们各自的值。本次讨论仅限于案例配置；第1.3节和第10章分别讨论了边界条件的数值和设计细节。顾名思义，边界条件定义网格边界处的场值。初始条件指内部场变量的初始值。图1.12显示了内部值和边界值之间的差异示意图。

有关计算网格的更多信息，请参见第2.1章。第1章讨论了数值方法的基本原理。读者应该在一定程度上了解这两个主题。

在仔细查看如何定义边界条件文件之前，必须将用icoFoam模拟的方腔实例复制到选择的位置。要设置方腔仿真实例的基本边界条件，需要复制并重命名仿真实例目录：

?>  cp -r $FOAM_TUTORIALS/incompressible/icoFoam/cavity/cavity cavityOscillatingU
?>  cd cavityOscillatingU

列出0/目录表明定义了两个不同的物理场：p和U。这两个文件都可以直接从命令行中使用任何文本编辑器进行编辑。可以使用以下命令将压力边界条件文件的相关行打印到屏幕上：

?>  cat 0/p | tail -n 23 | head -21

dimensions [0 2 -2 0 0 0 0];
internalField uniform 0;
boundaryField
{
    movingWall
    {
        type zeroGradient;
    }
    fixedWalls
    {
        type zeroGradient;
    }
    frontAndBack
    {
        type empty;
    } 
}

内容显示边界条件文件的三个顶级对象：dimensions、internalField和boundaryField。第一个对象是一组标量维度(dimsionSet)，用于定义物理场的量纲。每个标量对应于特殊SI单位的幂，如维度集（dimensionSet.H）的声明源文件中所定义。

// 定义维度指数名称的枚举
enum dimensionType
{
    MASS, // 质量 kg
    LENGTH, // 长度 m
    TIME, // 时间 s
    TEMPERATURE, // 温度 K
    MOLES, // 物质的量 mol
    CURRENT, // 电流 A
    LUMINOUS_INTENSITY // 光强 Cd
};

另一个对象是internalField，其定义了物理场的初始条件。请注意，这不包括由最后一个子字典：boundaryField定义的边界。在本例中，所有网格值都设置为0。还可以在每个网格格的基础上定义初始值，这反过来要求用户为每个网格组成一个包含所需值的列表。该列表中的网格数量必须与网格中的单元数量相同，其组成如[5]所述。最后，在boundaryField子目录中定义边界条件。必须为每个patch和每个物理场指定边界条件。因此，每个patch的边界条件在boundaryField字典的子字典中定义。

3.2.1 设置边界条件

OpenFOAM的官方发布已经附带了许多边界条件。在案例目录中执行foamHelp命令时，会提供所有可用边界条件的列表。 它说明了场类型，只显示了与相应场兼容的边界条件。结果列表中显示的一些边界条件非常一般，而其他边界条件则非常特定于问题。对于速度场U，可以用以下方式调用：

?>  foamHelp boundary -field U

foamHelp要求在检查可用边界条件时生成网格。作为定义边界条件的一个小例子，使用了之前复制的cavityOscillatingU案例。在本例中，将速度场U的动壁边界条件更改为确保速度随时间振荡的边界条件（图3.1）。要设置振荡速度，我们将使用codedFixedValue边界条件。这种边界条件允许其用户编写C++代码，定义边界条件如何运行。这个C++代码然后在上面的编码边界条件的C++代码片段中，我们重载了OpenFOAM场赋值操作符，以将速度向量设置为cos(πt)(1，0，0)，使其每秒切换一次方向。除了上述对U场的改变外，还必须延长模拟时间，以便可以可视化几个速度振荡周期。为此，必须更改system/control Dict字典：将条目endTime从0.5修改为4.0。所有需要的调整都已完成，更改后的边界条件可以使用icoFoam解算器进行测试。在执行此操作之前，必须使用block Mesh工具生成网格：OpenFOAM编译成仿真实例中的库，存储在dynamicCode子文件夹中。该库动态链接到求解器应用程序，并在模拟中使用。由于它是实现边界条件的C++代码，边界条件的用户应该知道如何在OpenFOAM中用标量、向量或张量来执行它们各自物理场的基本操作。

本书第二部分介绍了C++编程知识。一般来说，关于边界条件、function object和其他OpenFOAM组件的使用信息可以在扩展代码指南中找到，也在第5章中介绍。或者在OpenFOAM源代码中读取边界条件的头文件通常更快。必须使用文本编辑器打开文件中codedFixedValue边界条件的声明。该文件位于:

$FOAM_SRC/finiteVolume/fields/fvPatchFields/derived/codedFixedValue/codedFixedValueFvPatchField.H

头文件顶部的注释部分的Description部分包含边界条件用法的描述，我们在这里针对movingWall边界进行修改：

movingWall
{
    type codedFixedValue;
    value uniform (0 0 0);
    name oscillatingFixedValue; // name of generated BC
    code
    #{
        operator==(Foam::cos(this-> db().time().value() * M_PI)*vector(1,0,0));
    #};
}

在上面的边界条件的C++代码片段中，我们重载了OpenFOAM场赋值操作符，以将速度向量设置为$\cos (\pi t)(1,0,0)$，使其每秒切换一次方向。除了上述对速度场的改变外，还必须延长模拟时间，以便可以观察几个速度振荡周期。为此，必须更改system/controlDict字典：将条目endTime从0.5修改为4.0。所有需要的调整都已完成，更改后的边界条件可以使用icoFoam解算器进行测试。在执行此操作之前，必须使用blockMesh工具生成网格：

?>  blockMesh
?>  icoFoam

可以使用paraView实现速度场的可视化。

注：要使用ParaView查看OpenFOAM的计算结果，可以创建一个扩展名为.foam的空文件，然后利用ParaView打开。

在振荡边界条件的驱动下，内部流动发生脉动。

在选择边界条件时总是需要仔细考虑，特别是在将结果与实验(验证)或精确解(验证)进行比较时。在 foam _ tutorial cases 教程中包含了使用不同边界条件的模拟示例。教程案例可以是一个使用不同的边界条件有用的信息来源。

3.2.2 设置初始条件

在这一节中，概述了如何设置模拟的初始条件。初始条件在模拟开始时定义内部物理场的值。有各种不同的工具用于用户指定初始化物理场。第八章提供有关如何开发新的前后处理应用的资料。

注：设定初始条件(IC)称为前处理。

初始条件可以是微不足道的，就像方腔案例一样：将初始内部速度和压力设置为统一的向量值0。在第一个时间步长之后，计算的流动解通常与初始设置的流动解有很大不同。对于定常的单相不可压缩流动，初始条件可以看作是加速收敛到定常的初始条件。如果这个猜测离最终结果，可能会出现计算发散。

在其他情况下，初始条件是至关重要的，因为它们决定了场如何从初始值在时间上演变。可压缩流动模拟在很大程度上依赖于初始压力、温度和/或密度来正确地计算状态方程。不可压缩的多相模拟需要非常精确的用于分离流体相的相指示器场的初始值。初始相位指示器中的较大误差会导致更大的曲率近似误差，从而导致很强的数值不稳定性和可能的灾难性失败。

如前所述，空腔示例情况对初始条件的定义相对宽容。对初始条件有特殊要求的情况是使用interFoam解算器模拟的damBreak情况。求解器interFoam是一个两相流动求解器，使用代数流体体积方法区分两个不可混溶和不可压缩的流体相。为此，引入了一个新的标量场：alpha.water。本例中的气体和液相的alpha.water值分别为0和1。对于第一个示例，将向damBreak案例添加一个水滴，如图3.2所示。

首先，创建damBreak教程案例的副本，并将案例重命名：

?>  run
?>  cp -r FOAM_TUTORIALS/multiphase/interFoam/laminar/damBreak/damBreak damBreakWithDrop

原始场值通常保存在子文件夹0中，或保存在0文件夹中的*.orig字段文件中。要重置alpha.water的初始值，请将0/alpha.water文件复制到0/alpha.water：

?>  cp 0/alpha.water 0/alpha.water