Ansys Fluent提供了组分PDF输运模型，用于模拟湍流火焰中的有限速率化学效应。有关使用组分PDF输运模型的信息，请参阅用户指南中的“组分PDF输运问题建模”。以下各节介绍了此模型背后的理论：

• 7.2.1. 概述
• 7.2.2. 组分PDF输运理论
• 7.2.3. 拉格朗日求解方法
• 7.2.4. 欧拉求解方法

7.2.1 概述

组分PDF输运模型，类似于层流有限速率模型（参见直接使用有限速率动力学（无TCI）（第243页））和EDC模型（参见涡耗散概念（EDC）模型（第249页）），应在您对模拟湍流反应流中的有限速率化学动力学效应感兴趣时使用。通过适当的化学机理，可以预测如CO和NOx等动力学控制组分，以及火焰熄灭和点火。PDF输运模拟计算成本较高，建议您从小网格开始建模，最好是在二维中进行。

Ansys Fluent对组分PDF输运方程有两种不同的离散化方法，即拉格朗日和欧拉方法。拉格朗日方法在严格意义上比欧拉方法更准确，但需要更长的运行时间才能收敛。

有关Ansys Fluent中组分PDF输运模式的限制，请参阅Fluent用户指南中的限制。

7.2.2 组分PDF输运理论

湍流燃烧由反应的纳维-斯托克斯方程控制。尽管这一方程组是准确的，但其直接求解（即解析所有湍流尺度）对于实际的湍流流动而言，计算成本过高。在“Species Transport and Finite-Rate Chemistry”（第237页）中，组分方程经过雷诺平均处理，导致湍流标量通量和平均反应率的未知项出现。在Ansys Fluent中，湍流标量通量通过梯度扩散模型来模拟，将湍流对流视为增强的扩散过程。平均反应率可以通过Laminar、Eddy-Dissipation或EDC Finite-Rate化学模型来模拟。由于反应率总是高度非线性的，因此在湍流流动中模拟平均反应率既困难又容易出错。

另一种替代方法是对组分和能量方程进行雷诺平均，转而推导它们单点联合概率密度函数（PDF）的输运方程。这个PDF，记作$P$，可以被视为代表流体在每个组分和温度状态下所花费时间的比例。$P$在$N$个组分和温度空间中具有$N+1$维。从PDF中，可以计算任何单点热化学矩（例如，平均或均方根温度、平均反应率）。成分PDF输运方程是从纳维-斯托克斯方程推导出来的，如[528]（第1087页）所示：

$$\begin{array}{rl}\frac{\partial }{\partial t}\left(\rho P\right)+\frac{\partial }{\partial x_i}\left(\rho u_{i}P\right)+\frac{\partial }{\partial {\psi }_k}\left(\rho S_{k}P\right)& =-\frac{\partial }{\partial x_i}\left[\rho ⟨u_i^{\prime \prime }|\psi ⟩P\right]\\ & +\frac{\partial }{\partial {\psi }_k}[\rho ⟨\frac{1}{\rho }\frac{\partial J_{ik}}{\partial x_i}|\psi ⟩P]\end{array}\text{\hspace{1em}(7.148)}$$

式中

$$P=\text{Favre joint PDF of composition}$$

$\rho =$ 流体平均密度

$u_i=$ Favre 平均流体速度矢量

$S_k=$ 组分 $k$ 的反应速率

$\psi =$ 组成空间矢量

$u_i^{\prime \prime }=$ 流体速度波动矢量

$$J_{i,k}=\text{molecular diffusion flux vector}$$

符号 $⟨\dots ⟩$ 表示期望值，而 $⟨A\mid B⟩$ 是事件 $A$ 在事件 $B$ 发生的条件下的条件概率。

在方程 7.148（第277页）中，左侧的项是封闭的，而右侧的项则不是，需要建模。左侧的第一个项是不稳定状态下PDF的变化率，第二个项是由于平均速度场对流引起的PDF变化，第三个项是由于化学反应引起的PDF变化。PDF输运方法的主要优势在于，高度非线性的反应项是完全封闭的，无需建模。右侧的两个项分别表示由于湍流引起的标量对流（湍流标量通量）和分子混合/扩散引起的PDF变化。

湍流标量通量项是未封闭的，在Ansys Fluent中通过梯度扩散假设进行建模。

$$-\frac{\partial }{\partial x_i}\left[\rho ⟨u_i^{\prime \prime }\mid \psi ⟩P\right]=\frac{\partial }{\partial x_i}\left(\frac{\rho {\mu }_t}{S{c}_t}\frac{\partial P}{\partial x_i}\right)\text{\hspace{1em}(7.149)}$$

其中，${\mu }_t$ 表示湍流粘度，而 $S{c}_t$ 是湍流Schmidt数。正如《湍流》(第41页)中所述，进行成分PDF输运模拟需要一个湍流模型，该模型决定了 ${\mu }_t$。

由于描述的是单点PDF，因此缺乏相邻点的信息，所有梯度项（如分子混合）都是未封闭的，必须进行建模。混合模型至关重要，因为在反应物和热量共同扩散的最小分子尺度上会发生燃烧。在PDF方法中对混合进行建模并不简单，这也是PDF输运方法中最薄弱的环节。关于混合模型的描述，请参见《粒子混合》(第279页)。

7.2.3 拉格朗日求解方法

采用拉格朗日蒙特卡罗方法来求解 $N$ +1 维的PDF输运方程。蒙特卡罗方法适用于高维方程，因为计算成本随维数线性增加。其缺点是引入了统计误差，这些误差必须得到仔细控制。

为了求解建模后的PDF输运方程，我们将其比拟为一个具有相同解的随机微分方程(SDE)。蒙特卡罗算法涉及的虚拟粒子因粒子对流而在物理空间中随机移动，同时也因分子混合和反应而在成分空间中移动。这些粒子具有质量，并且在平均意义上，一个单元内粒子质量的总和等于该单元的质量（单元密度乘以单元体积）。由于实际网格单元体积变化很大，因此需要调整粒子质量，以控制单元内的粒子数量大致恒定且均匀。

对流、扩散和反应过程被分解为若干步骤进行处理，具体细节将在后续章节中阐述。关于分数步长法的详细信息，请参阅[92]（第1062页）。

该方法的介绍将在以下章节中展开：

• 7.2.3.1. 粒子对流
• 7.2.3.2. 粒子混合
• 7.2.3.3. 粒子反应

7.2.3.1 粒子对流

在Ansys Fluent中，采用了一种空间上二阶精度的拉格朗日方法，该方法包含两个步骤。在第一步对流阶段，粒子被推进到新的位置。

$$x_i^{1/2}=x_i^0+\frac{1}{2}{u}_i^0\Delta t\text{\hspace{1em}(7.150)}$$

其中

$x_i=$ 粒子位置矢量

$u_i=$ 粒子位置处的流体速度矢量的Favre平均值

$\Delta t=$ 粒子时间步长

对于非定常流动，粒子时间步长即为物理时间步长。对于稳态流动，则为每个单元计算局部时间步长。

$$\Delta t=\min \left(\Delta t_{\text{conv }},\Delta t_{\text{diff }},\Delta t_{\text{mix }}\right)\text{\hspace{1em}(7.151)}$$

其中

$\Delta t_{\text{conv }}=$ 对流数 $\times \Delta x/\text{(单元流体速度)}$

$\Delta t_{\text{diff }}=$ 扩散数 $\times {\left(\Delta x\right)}^2/\left(\text{ 单元湍流扩散率 }\right)$

$\Delta t_{\text{mix }}=$ 混合数 $\times \text{湍流时间尺度}$

$\Delta x=$ 特征单元长度 =${\mathrm{体积}}^{1/D}$ 其中 $D$ 是问题维度

在第一次对流步骤之后，包括扩散和反应在内的所有其他子过程都会被处理。

最后，第二次对流步骤计算为 $$x_i^1=x_i^{1/2}+\Delta t\left(u_i^{1/2}-\frac{1}{2}{u}_i^0+\frac{1}{\rho S{c}_t}\frac{\partial {\mu }_{eff}}{\partial x_i}+{\xi }_i\sqrt{\frac{2{\mu }_{eff}}{\rho \Delta tS{c}_t}}\right)\text{\hspace{1em}(7.152)}$$

其中

$\rho =$ 平均单元流体密度

$u_i=$ 颗粒位置处的平均流体速度矢量

${\mu }_{eff}=$ 有效粘度

$S{c}_t=$ 湍流施密特数

${\xi }_i=$ 标准化正态随机向量

7.2.3.2 粒子混合

组分和热的分子混合必须被建模，并且通常是 PDF 输运方法中最大的建模误差来源。Ansys Fluent 提供了三种分子扩散模型：修正的 Curl 模型 [278] (p. 1073), [480] (p. 1085)，IEM 模型（有时称为 LSME 模型）[146] (p. 1065) 和 EMST 模型 [635] (p. 1093)。

7.2.3.2.1 修正的 Curl 模型

对于修正的Curl模型，从单元格中的所有粒子中随机选择几对粒子，并将它们的个体成分向其平均成分移动。在粒子质量相等的特殊情况下，所选粒子对的数量计算为： $$N_{\text{pair }}=\frac{1.5C_{\phi }N\Delta t}{{\tau }_t}\text{\hspace{1em}(7.153)}$$

其中 $N=$ 单元格中的粒子总数

$C_{\phi }=$ 混合常数（默认值 = 2）

${\tau }_t=$ 湍流时间尺度（对于 $k-\epsilon$ 模型，这是 $k/\epsilon$）

对于粒子质量可变的一般情况，使用 [480]（第 1085 页）中的算法。

对于每一对粒子，选择一个均匀随机数 $\xi$，每个粒子的成分 $\phi$ 向这对粒子的平均成分移动，移动因子与 $\xi$ 成正比： $${\phi }_i^1=\left(1-\xi \right){\phi }_i^0+\xi \frac{\left({\phi }_i^0{m}_i+{\phi }_j^0{m}_j\right)}{\left(m_i+m_j\right)}\text{\hspace{1em}(7.154)}$$

$${\phi }_j^1=\left(1-\xi \right){\phi }_j^0+\xi \frac{\left({\phi }_i^0{m}_i+{\phi }_j^0{m}_j\right)}{\left(m_i+m_j\right)}$$

其中 ${\phi }_i$ 和 ${\phi }_j$ 分别是粒子 $i$ 和 $j$ 的组成向量，$m_i$ 和 $m_j$ 分别是粒子 $i$ 和 $j$ 的质量。

7.2.3.2.2 IEM 模型

对于交互通过均值交换（IEM）模型，单元内所有颗粒的组成向均值方向移动一小段距离： $${\phi }^1={\phi }^0-\left(1-e^{-0.5C_{\phi }\Delta t/{\tau }_t}\right)\left({\phi }^0-\tilde{\phi }\right)\text{\hspace{1em}(7.155)}$$

其中 ${\phi }^0$ 是混合前的组成，${\phi }^1$ 是混合后的组成，$\bar{\phi }$ 是颗粒位置处的 Favre 平均组成向量。

7.2.3.2.3 EMST 模型

从物理角度来看，混合发生在彼此相邻的流体颗粒之间。改进的 Curl 和 IEM 混合模型没有考虑这种局部性，这可能是误差的来源。欧几里得最小生成树（EMST）模型混合在组成空间中彼此接近的颗粒对。由于标量场在局部是平滑的，因此在组成空间中接近的颗粒在物理空间中也很可能接近。颗粒配对由欧几里得最小生成树决定，这是将一个颗粒连接到至少另一个颗粒的边的集合的最小长度。EMST 混合模型比改进的 Curl 和 IEM 混合模型更准确，但会稍微增加计算成本。有关 EMST 模型的详细信息，请参阅参考文献 [635]（第 1093 页）。

7.2.3.2.4 液体反应

液体中的反应通常在低湍流水平（小 Re）下发生，反应物的扩散性较低（大 Sc）。对于此类流动，默认的混合常数 $C_{\phi }=2$ 会高估混合速率。液体微混合选项从模型湍流 [530]（第 1088 页）和标量 [187]（第 1067 页）谱中插值 $C_{\phi }$。

7.2.3.3 粒子反应

粒子组成向量表示为 $$\phi =\left(Y_1,Y_2,\dots ,Y_N,T\right)\text{\hspace{1em}(7.156)}$$

其中 $Y_k$ 是第 $k$ 种物质的质分数，$T$ 是温度。

对于反应分数步长，反应源项被积分如下： $${\phi }^1={\phi }^0+{\int }_0^{\Delta t}Sdt\text{\hspace{1em}(7.157)}$$

其中 $S$ 是化学源项。大多数真实的化学机理包含数十种物质和数百种反应。通常，反应在达到点火温度之前不会发生，但一旦发生，就会非常迅速地进行，直到反应物被消耗完。因此，一些反应的时间尺度非常快，大约为 ${10}^{-10}\mathrm{s}$，而其他反应的时间尺度则慢得多，大约为 $1\mathrm{s}$。这种时间尺度的差异导致数值刚性，这意味着需要大量的计算工作来积分公式 7.157（第 281 页）中的化学源项。在 Ansys Fluent 中，反应步骤（即计算 ${\phi }^1$）可以通过直接积分或现场自适应制表（ISAT）来执行，如下文所述。

在 Ansys Fluent 中，典型的稳态 PDF 传输模拟可能包含 50000 个单元，每个单元有 20 个粒子，并且需要 1000 次迭代才能收敛。因此，至少需要 ${10}^{\circ }$ 次刚性 ODE 积分。由于每次积分通常需要几十或几百毫秒，平均而言，化学的直接积分对 CPU 的要求极高。

对于给定的反应机理，方程 7.157（第 281 页）可以视为一种映射。在初始组成向量 ${\phi }^0$ 的情况下，最终状态 ${\phi }^1$ 仅取决于 ${\phi }^0$ 和映射时间 $\Delta t$。理论上，如果在模拟前能够构建一张涵盖所有可实现 ${\phi }^0$ 状态和时间步长的表格，那么通过查表就可以避免积分运算。但实际上，这种先验的表格构建是不可行的，因为需要一个完整的 $N+3$ 维表格（$N$ 种物质、温度、压力和时间步长）。为了说明这一点，考虑一个每个维度有 $M$ 个点的结构化表格。所需的表格大小为 $M^{N+3}$，而对于 $M=10$ 个离散点和 $N=7$ 种物质的保守估计，表格将包含 ${10}^{10}$ 条目。

进一步审视，存储整个可实现空间的完整数据是非常浪费的，因为大多数区域从未被访问。例如，在一个实际的燃烧器中找到 $Y_{OH}=1$ 和 $T=300K$ 的组成是不现实的。事实上，对于稳态的三维层流模拟，化学反应可以通过空间位置向量参数化。因此，映射必须位于 $N+3$ 维组成空间内的一个三维流形上。因此，只需对组成空间的这一访问区域进行表格化就足够了。

然而，访问的区域取决于特定的化学机制、分子传输特性、流动几何形状以及边界条件。因此，在模拟之前无法得知访问的区域，表格也无法预处理。相反，表格必须在模拟过程中建立，这被称为原位表格化。Ansys Fluent采用ISAT [526]（第1087页）动态地表格化化学映射关系，并加速求解时间。ISAT是一种在飞行中（原位）带误差控制的（自适应表格化）表格化访问组分空间区域的方法。当正确使用ISAT时，通常可以实现两到三个数量级的加速。然而，为了最优地使用ISAT，理解其工作原理至关重要。

7.2.4 欧拉求解方法

拉格朗日求解方法通过随机跟踪拉格朗日粒子穿过域来求解组分PDF输运方程。由于需要大量的粒子来表示PDF，并且需要大量的迭代来减少统计误差并显式地输运粒子通过域，因此计算成本较高。欧拉PDF输运模型通过假设PDF的形状来克服这些限制，从而可以推导出欧拉输运方程。消除了随机误差，并且输运方程隐式求解，计算上更为经济。多维PDF形状假设为δ函数的乘积。与拉格朗日PDF模型一样，高度非线性的化学源项被封闭。然而，湍流标量通量和分子混合项必须建模，并分别用梯度扩散和IEM模型封闭。

$N_s+1$ 维度的组成 PDF（$N_s$ 种组分和焓）被表示为 $N_e$ 个 delta 函数（或模式）的集合。这个假定的 PDF 具有以下形式： $$P\left(\psi ;\overrightarrow{x},t\right)=\sum _{n=1}^{N_e}{p}_n\left(\overrightarrow{x},t\right)\prod _{k=1}^{N_s}\delta \left[{\psi }_k-<{\phi }_k{>}_n\left(\overrightarrow{x},t\right)\right]\text{\hspace{1em}(7.158)}$$ 其中，$p_n$ 是每个模式中的概率，$<{\phi }_k{>}_n$ 是第 $n^{\text{th}}$ 模式中组分 $k$ 的条件平均组成，${\psi }_k$ 是组分 $k$ 的组成空间变量，而 $\delta \left(\dots \right)$ 是 delta 函数。

欧拉 PDF 输运方程是通过将方程 7.158（第 282 页）代入封闭的组成 PDF 输运方程（方程 7.148（第 277 页）与方程 7.149（第 278 页）和方程 7.155（第 280 页））推导出来的。未知项 $p_n$ 和 $<{\phi }_k{>}_n$ 是通过强制这个输运 PDF 的较低阶矩与 RANS 较低阶矩输运方程匹配来确定的，使用直接矩量法（DQMOM）方法 [187]（第 1067 页），[412]（第 1081 页）。由此得到的输运方程为：

• 概率（第 $n^{\text{th}}$ delta 函数的幅度）： $$\frac{\partial \rho p_n}{\partial t}+\frac{\partial }{\partial x_i}\left(\rho u_i{p}_n\right)=\nabla \left(\rho \Gamma \nabla p_n\right)\text{\hspace{1em}(7.159)}$$

• 组分 $k$ 的概率加权条件平均值： $$\frac{\partial \rho s_{k,n}}{\partial t}+\frac{\partial }{\partial x_i}\left(\rho u_i{s}_{k,n}\right)=\nabla \left(\rho \Gamma \nabla s_{k,n}\right)+\rho \left(M_{k,n}+S_{k,n}+C_{k,n}\right)\text{\hspace{1em}(7.160)}$$

其中，$p_n$ 是第 $n^{\text{th }}$ 模式的概率，$s_{k,n}=p_n<{\phi }_k{>}_n$ 是第 $n^{\text{th }}$ 模式的第 $k^{\text{th }}$ 组分概率加权条件平均组成。$\Gamma ={\mu }_l+{\mu }_t/S{c}_t$ 是有效湍流扩散系数。项 $M_{k,n}$、$S_{k,n}$ 和 $C_{k,n}$ 分别表示混合、反应和修正项。注意，只求解 $N_e-1$ 个概率方程，第 $N^{\text{th }}$ 概率计算为 1 减去 $N_e-1$ 个已求解概率的总和。

更多信息，请参见以下章节：

• 7.2.4.1. 反应
• 7.2.4.2. 混合
• 7.2.4.3. 修正
• 7.2.4.4. 组成平均值和方差的计算

7.2.4.1 反应

方程 7.160（第 282 页）中第 $k^{\text{ th }}$ 组成和第 $n^{\text{ th }}$ 模式的反应源项 $S_{k,n}$ 计算如下： $$S_{k,n}=p_{n}S{\left(<{\phi }_k{>}_n\right)}_k\text{\hspace{1em}(7.161)}$$

其中，$S{\left(\right)}_k$ 是第 $k^{\text{th }}$ 组分的净反应速率。

7.2.4.2 混合

微混合项 $M_{k,n}$ 采用 IEM 混合模型进行建模： $$M_{k,n}=\frac{C_{\phi }}{\tau }\left(<{\phi }_k{>}_n-{\psi }_k\right)\text{\hspace{1em}(7.162)}$$ 其中，$\tau$ 是湍流时间尺度，$C_{\phi }$ 是混合常数。

因此，对于两模式 DQMOM-IEM 模型，组分 $k$ 的混合项为： $$M_{k,1}=\frac{C_{\phi }}{\tau }\left(p_1{s}_{k,2}-p_2{s}_{k,1}\right)\text{\hspace{1em}(7.163)}$$

$$\begin{array}{ccc}{M}_{k,2}& =& \frac{C_{\phi }}{\tau }\left(p_{{}_2}{s}_{k,1}-p_{{}_1}{s}_{k,2}\right)\end{array}$$ The default value of $C_{\phi }$ is 2, which is appropriate for gas-phase combustion. For reactions in liquids, where the diffusivities are much smaller than gases, the Liquid Micro-Mixing option interpolates $C_{\phi }$ from model turbulence [530] (p. 1088) and scalar [187] (p. 1067) spectra.

7.2.4.3 Correction

Using assumptions to ensure realizability and boundedness, the correction terms $C_{k,n}$ in Equa-

tion 7.160 (p. 282) for the $k^{\text{ th }}$ composition are determined from the linear system, $$\sum _{n=1}^{N_e}<{\phi }_k{>}_n^{m_k-1}{C}_{k,n}=\sum _{n=1}^{N_e}\left(m_k-1\right)<{\phi }_k{>}_n^{m_k-2}{p}_n{c}_{k,n}\text{\hspace{1em}(7.164)}$$

7.2.4.3 修正

利用假设确保可实现性和有界性，方程 7.160（第 282 页）中第 $k$ 个组分的修正项 $C_{k,n}$ 由线性系统确定， $$\sum _{n=1}^{N_e}<{\phi }_k{>}_n^{m_k-1}{C}_{k,n}=\sum _{n=1}^{N_e}\left(m_k-1\right)<{\phi }_k{>}_n^{m_k-2}{p}_n{c}_{k,n}\text{\hspace{1em}(7.164)}$$

其中 $m_k$ 是每个组分 $k$ 的非负整数低阶矩 $\left(1-N_e\right)$。注意，随着 $m_k$ 的增加，矩阵的条件数降低，这会减少高阶模态模拟的稳定性。

方程 7.164（第 283 页）中的耗散项 $c_{k,n}$ 计算如下， $$c_{k,n}=\Gamma \nabla \left(<{\phi }_k{>}_n\right)\cdot \nabla \left(<{\phi }_k{>}_n\right)\text{\hspace{1em}(7.165)}$$

对于双模态 DQMOM-IEM 模型，第 $k$ 个组分的修正项为， $$C_{k,1}=\frac{\Gamma }{<{\phi }_k{>}_1-<{\phi }_k{>}_2}\left(p_1{\left(\nabla <{\phi }_k{>}_1\right)}^2+p_2{\left(\nabla <{\phi }_k{>}_2\right)}^2\right)\text{\hspace{1em}(7.166)}$$

$$C_{k,2}=-C_{k,1}$$

7.2.4.4 成分均值和方差的计算

平均成分（组分 $k$ 或能量）的计算方法如下： $$\overline{{\phi }_k}=\sum _{n=1}^{N_e}{s}_{k,n}\text{\hspace{1em}(7.167)}$$

平均成分（组分 $k$ 或能量）的计算公式为： $$\overline{{\phi }_k}=\sum _{n=1}^{N_e}{s}_{k,n}\text{\hspace{1em}(7.167)}$$

其方差的计算公式为： $$\overline{{\phi }_k^{\prime 2}}=\left(\sum _{n=1}^{N_e}{s}_{k,n}<{\phi }_k{>}_n\right)-{\bar{\phi }}_k^2\text{\hspace{1em}(7.168)}$$