Computational Fluid Dynamics:Principles and Applications

第 1 章 绪论 计算流体力学（computational fluid dynamics，CFD）的历史始于 20 世纪 70 年代初。大约在那个时期，CFD 开始被用作一个缩写，指将物理学、数...

第 2 章 控制方程 2.1 流动及其数学描述 在开始推导描述流体行为的基本方程之前，有必要先澄清“流体动力学”这一术语的含义。它实际上研究的是大量个体粒子相互作用时的运动；在这里，这些粒子是分...

2.2 守恒定律 2.2.1 连续性方程 若只考虑单相流体，质量守恒定律表达的是：在这样的流体系统中，质量既不能产生，也不能消失。连续性方程中也不存在扩散通量的贡献，因为对于静止流体，任何质量变...

2.3 黏性应力 流体与微元表面之间的摩擦所产生的黏性应力，由应力张量 $ \overline{\overline{\tau}} $ 描述。在 笛卡尔（Cartesian）坐标系中，其一般形式为...

2.4 完整的 Navier–Stokes 方程组 前面几节分别推导了质量、动量和能量守恒定律。现在可将它们组合成一个方程组，以便更清楚地考察其中的各项。为此，回到式（2.2）所表达的矢量物理量...

参考文献 [1] Lax PD. 非线性双曲型方程的弱解及其数值计算（Weak solutions of nonlinear hyperbolic equations and their num...

第 3 章 控制方程的求解原理 上一章给出了完整的 Navier–Stokes／Euler 方程组，补充了完全气体的热力学关系，并为化学反应气体定义了附加输运方程。现在，我们可以求解整个控制方程...

3.2 时间离散 如本章开头所述，求解 Euler 方程和 Navier–Stokes 方程的绝大多数数值格式采用线方法，即分别进行空间与时间离散。这种方法具有最大的灵活性，因为可以根据待求问题...

3.3 湍流建模 对于无黏流或层流，求解控制方程（2.19）不存在根本性困难。然而，湍流模拟却是一个重大问题。尽管现代超级计算机性能很高，利用含时 Navier–Stokes 方程（2.19）直...

3.4 初始条件与边界条件 无论选择何种数值方法求解控制方程（2.19），都必须给出合适的初始条件与边界条件。初始条件确定 t = 0 时，或迭代格式第一步时的流体状态。显然，初始猜测越好（越接...

参考文献 [1] Frymier PD, Hassan HA, Salas MD. Navier-Stokes calculations using cartesian grids: I. La...

第 4 章 结构有限体积格式 如第 3 章引言所述，求解 Euler 和 Navier–Stokes 方程的绝大多数数值格式采用线方法，即分别进行空间和时间离散。因此，可以根据待求问题的需要，对...

4.2 一般离散方法 本章引言已提到控制体定义与流动变量布置的三种方法。本节将逐一详细介绍，并讨论其优缺点。 4.2.1 单元中心格式 当控制体与网格单元相同、流动变量定义在单元形心（如图 4....

4.3 对流通量的离散 前面各节讨论了空间离散的一般方法。本节将详细介绍如何近似对流通量。 如第 3.1.5 节所述，在有限体积方法中，基本可以选择： 中心格式； 通量矢量分裂格式； ...

4.4 黏性通量的离散 为获得相容的空间离散，黏性通量通常使用与对流通量相同的控制体。唯一例外是采用重叠控制体的单元顶点格式（第 4.2.2 节）；主要出于稳定性考虑 $ [97–99] $，黏...

参考文献 [1] Hirsch C. Numerical computation of internal and external flows, vols. 1 and 2. Chicheste...

第 5 章 非结构有限体积格式 如第 3 章引言所述，绝大多数求解 Euler 方程和 Navier–Stokes 方程的数值格式采用线方法，即分别进行空间离散和时间离散。其主要优点是可以为空间...

5.2 一般离散方法 本章开头已提到，定义控制体及确定流动变量位置的两种最常用方法，是单元中心格式和中位对偶格式。本节将更详细地介绍它们。 不过，首先简要讨论非结构网格流动求解器所需的基本数据结...

5.3 对流通量的离散 前几节讨论了空间离散方法的一般问题及所需数据结构。下文进一步介绍对流通量计算的具体实现。 如 3.1.5 节所述，在有限体积方法框架内，主要可选择以下格式： 中心格式...

5.4 黏性通量的离散 计算式（5.2）的扩散通量 $ \bar{F}_{ v} $，需要知道控制体面上的流动量及其一阶导数。为方便，黏性通量采用与对流通量相同的控制体，以获得相容的空间离散并简...

参考文献 [1] Hirsch C. Numerical computation of internal and external flows, vols. 1 and 2. Chicheste...

第 6 章 时间离散 采用线方法，即对控制方程（2.19）分别进行空间离散和时间离散，并为每个控制体写出方程，得到以下关于时间的耦合常微分方程组： $$ \frac{\mathrm{d}\lef...

6.2 隐式时间推进格式 在式（6.2）中取 $ \beta\neq0 $，可得到各种隐式时间积分格式。对定常流动问题，取 $ \omega=0 $ 的隐式格式最适合；非定常流动见 6.3 节。...

6.3 非定常流动的求解方法 非定常流动模拟如今在许多工程领域十分重要，例如叶轮机械静止与旋转部件的相互作用、活塞发动机、流固耦合、直升机空气动力学、气动声学、湍流的 DNS 或 LES、爆轰等...

参考文献 [1] Jameson A, Schmidt W, Turkel E. Numerical solutions of the Euler equations by finite vol...

第 7 章 湍流建模 湍流与层流相比，最突出的特征是流体分子沿复杂而不规则的路径作混沌运动。强烈的混沌运动使流体各层之间发生剧烈混合。由于分子之间以及分子与固体壁面之间的动量和能量交换增强，在相...

7.2 一阶封闭模型 一阶封闭模型是近似计算 Reynolds／Favre 平均 Navier–Stokes 方程中 Reynolds 应力的最简便方法。它们基于 7.1.5、7.1.6 节分别...

7.3 大涡模拟 Smagorinsky 早在 1963 年就在气象学（大气环流）中采用了所谓的 LES 方法 $ [98] $。Deardorff 于 1970 年给出了 LES 的首个工程应...

参考文献 [1] Liu C, Liu Z. Multigrid methods and high order finite difference for flow in transition....

第 8 章 边界条件 任何数值模拟都只能考虑真实物理区域或系统的一部分。区域截断会引入人工边界，必须在这些边界上指定某些物理量的值。此外，与流体接触的壁面构成物理区域的自然边界。边界条件的数值处...

8.2 固体壁面 8.2.1 无黏流动 在无黏流动中，流体沿表面滑动。由于没有摩擦力，速度矢量必须与表面相切，等价于表面法向没有流动，即 $$ \vec{ v}\cdot\vec{n}=0\qq...

8.3 远场 对翼型、机翼、汽车及其他构型周围的外部流动进行数值模拟，必须采用有限区域，因此需要人工远场边界条件。远场边界条件的数值实现必须满足两个基本要求：首先，与无限区域相比，截断区域不应对...

8.4 入口／出口边界 人们为 Navier–Stokes 方程的入口，尤其是出口（又称开放）边界条件，提出了多种数值实现方法 $ [22–26] $。这里重点介绍为叶轮机械应用发展的方法。适宜...

8.5 喷注边界 喷注速度和其他流动变量由给定的质量流量 $ \dot{m} $ 及喷注温度 $ T_{inj} $ 计算。这里假设质量沿边界法向喷入。对流通量采用式（2.21）计算，逆变速度 ...

8.6 对称面 若流动关于某条直线或某个平面对称，首先必须满足没有穿过边界的通量，即对称边界上的法向速度为零。此外，以下梯度必须为零： 标量沿边界法向的梯度； 切向速度沿边界法向的梯度；...

8.7 坐标切口 这种边界条件仅出现在结构网格中。坐标切口是人工边界，并非物理边界；它由计算坐标不同、物理位置相同的网格点构成一条线，三维时为一个面。这意味着网格经过折叠后与自身接触。11.1....

8.8 周期边界 某些实际应用中的流场沿一个或多个坐标方向具有周期性，此时只需在其中一个重复区域内模拟流动，并通过周期边界条件保证其与其余物理区域之间的正确相互作用。 周期边界有两种基本类型。第...

8.9 网格块之间的接口 在 3.1 节讨论结构网格空间离散时，已经指出通常无法在几何形状复杂的区域中生成单个网格（见图 3.4）。当时介绍了两种解决方法：多块方法和 Chimera 方法。下面...

8.10 非结构网格边界处的流动梯度 5.3.4 节已指出，中位对偶格式计算流动梯度时需要谨慎处理。如果在三角形或四面体网格上采用式（5.50）的 Green–Gauss 方法计算梯度，区域边界...

参考文献 [1] Kroll N, Jain RK. Solution of two-dimensional Euler equations—experience with a finite v...

第 9 章 加速技术 为加速定常问题中控制方程（2.19）的求解，人们发展了多种方法。这些加速技术也适用于双时间推进格式（6.3 节）的内迭代。本章讨论以下数值方法： 局部时间推进； 焓...

9.2 焓阻尼 某些情况下，总焓 H（式（2.12））在整个流场中保持常数，例如由 Euler 方程描述、没有热源及外力的外部流动。可以利用这一特点降低计算成本。 图 9.1 非结构网格无黏...

9.3 残差光顺 通过优化级系数 $ [5-7] $，可在一定程度上改变显式多级时间推进格式（6.1.1、6.1.2 节）的最大 CFL 数和收敛性质。Jameson 和 Baker $ [8]...

9.4 多重网格 多重网格是一种强有力的加速技术（参见图 3.7）。它在一系列逐渐变粗的网格上求解控制方程，然后将粗网格的解更新量合并并加到最细网格的解上。Brandt $ [19] $ 最初为...

9.5 低马赫数预处理 在低亚声速马赫数范围，当流速远小于声速时，控制方程（2.19）的对流项出现刚性。以下示例可以说明这一点。三维情况下，五个特征值为 $$ \begin{aligned}(\...

9.6 并行化 如今，几乎所有计算机——台式机、笔记本电脑和平板电脑——的中央处理器（CPU）都由多个核心组成。此外，现代图形处理器（GPU）还提供数百乃至数千个可进行通用编程的处理单元。然而，...

参考文献 [1] Kroll N, Jain RK. Solution of two-dimensional Euler equations—experience with a finite v...

第 10 章 相容性、精度与稳定性 在有限体积方法中，需要采用适当方法，近似计算每个控制体上 Navier–Stokes 方程（2.19）的面积分。第4章和第5章介绍了多种空间离散方案。由于对穿...

10.2 离散格式的精度 离散格式的精度与其截断误差有关。例如，若截断误差的主项与 $ \Delta x $ 成正比，则称为空间一阶精度格式；若主项随 $ (\Delta x)^{2} $ 变化...

10.3 von Neumann 稳定性分析 在实现新的离散方法之前，至少应大致了解它对数值格式稳定性与收敛行为的影响。已有大量实践证实，von Neumann 稳定性分析能够可靠地评估求解格式...

参考文献 [1] Roache PJ. Quantification of uncertainty in computational fluid dynamics. Ann Rev Fluid ...

第 11 章 网格生成原理 在数值求解控制方程之前，必须离散所有边界表面，并在流动区域内部生成体网格。如第3章开头的3.1节所述，基本上可选择： 结构网格； 非结构网格。 图11.1...

11.2 非结构网格 非结构网格通常由二维三角形或三维四面体组成。但如今，用多种单元类型构造非结构网格越来越常见。例如，以六面体或棱柱体离散边界层，流动区域其余部分用四面体填充，并在六面体或棱柱...

参考文献 [1] Brodersen O, Hepperle M, Ronzheimer A, Rossow CC, Schöning B. The parametric grid genera...

第 12 章 软件应用 软件包包含若干一维、二维流动求解器及网格生成器的源代码和可执行程序，并提供多种一维、二维流动算例的输入数据与网格。此外，还包含两个用于显式、隐式时间推进格式 von Ne...

12.2 一维结构网格生成器 grid_1d 目录包含生成一维结构网格的程序。每个网格节点对应一个截面积，沿 x 轴的面积分布对应 Laval 喷管，按以下关系计算： $$ \begin{gat...

12.3 二维结构网格生成器 grid_struct_2d 目录包含三个生成二维外流与内流结构网格的程序。第一个的源代码位于 cgrid 子目录，用于生成翼型周围 C 型网格，见11.1.1节及...

12.4 结构网格转非结构网格程序 grid_unstr_2d 目录中的程序，将二维结构网格转换为非结构三角形网格。它连接距离较短的一对对角节点，将结构网格的每个四边形单元划分为两个三角形。所得...

12.5 准一维 Euler 求解器 structured_1d 目录包含求解准一维 Euler 方程的程序，描述喷管或管道中的无黏流动。方程的守恒微分形式为 $ [3] $ $$ \frac{...

12.6 二维结构 Euler／Navier–Stokes 求解器 structured_2d 目录包含在结构贴体网格上求解二维 Euler 方程和 Navier–Stokes 方程的程序。空间...

12.7 二维非结构 Euler／Navier–Stokes 求解器 unstructured_2d 和 C++/unstructured_2d 目录包含在非结构三角形网格上求解二维 Euler...

12.8 并行化 Parallelization 目录包含采用异步代理库、CUDA、MPI 和 OpenMP 实现的并行程序示例，见9.6节。待并行化的程序求解二维 Laplace 方程， $$...

参考文献 [1] Farin GE. Curves and surfaces for computer aided geometric design—a practical guide. 3rd...

附录 A.1 控制方程的微分形式 假设对流通量与黏性通量连续，则先应用 Gauss 定理，可将控制方程（2.19）由积分形式转为微分形式，得到 $$ \frac{\partial}{\parti...

A.2 Euler 方程的拟线性形式 这种形式有助于理解控制方程的数学性质，也便于发展和分析数值格式。 拟线性形式由控制方程的微分形式（A.4）推导得到。这里仅考虑 Euler 方程，即 $$ ...

A.3 控制方程的数学类型 偏微分方程（PDE）的数学类型，可用经典的二阶拟线性方程说明： $$ a\frac{\partial^{2}U}{\partial x^{2}}+b\frac{\pa...

A.4 旋转参考系中的 Navier–Stokes 方程 在某些情况下，例如叶轮机械、螺旋桨或地球物理应用中，计算域绕某条轴匀速旋转。此时，将 Navier–Stokes 方程变换到旋转参考系更...

A.5 动网格上的 Navier–Stokes 方程 在研究流固耦合或模拟外挂物分离等问题时，必须在移动、甚至变形的网格上求解控制方程。最常用的两种方法，是任意 Lagrange–Euler（A...

A.6 薄剪切层近似 高 Reynolds 数流动中，黏性应力仅在物体周围的狭窄区域影响流动。因此，可合理假设壁面法向梯度占主导，忽略其他方向的梯度 $ [21, 22] $。这就是 Navie...

A.7 PNS 方程 当满足以下三个条件时： 流动为定常，即 $\partial\vec W/\partial t=0$； 流体主要沿一个方向运动，且不得发生边界层分离； 横流分量可...

A.8 Navier–Stokes 方程的轴对称形式 模拟圆柱形物体绕流或圆柱形腔内流动时，通常可假设流动关于纵轴对称。这意味着没有旋流，且流动量沿周向的导数为零。在这些条件下，可从 Navie...

A.9 对流通量 Jacobian 矩阵 对流通量 Jacobian 表示对流通量对守恒变量的梯度，即 $$ \bar{A}_{c}=\frac{\partial\vec{F}_{c}}{\pa...

A.10 黏性通量 Jacobian 矩阵 黏性通量 Jacobian 表示黏性通量对守恒变量的梯度，即 $$ \bar{A}_{ v}=\frac{\partial\vec{F}_{\math...

A.11 从守恒变量到特征变量的变换 对流通量 Jacobian（A.9节）具有实特征值和完整的特征向量组，因此可对角化为 $ [25] $ $$ \bar{A}_{c}=\bar{T}\bar...

A.12 GMRES 算法 考虑线性方程组 $$ \bar{A}\vec{x}=\vec{b}.\tag{A.90} $$ 寻求以下形式的近似解： $$ \vec{x}=\vec{x}_{0}...

A.13 张量记号 坐标 $ (x_{i}) $ 或速度分量 $ ( v_{i}) $ 这样的表达式表示一阶张量，具有三个分量，对应矢量。因此， $$ \begin{aligned}x_{i}&...

参考文献 [1] Sommerfeld AJW. Partial differential equations in physics. New York: Academic Press; 194...