12.9.1 简介

Ansys Fluent 针对模拟液体液滴撞击壁面并形成薄膜的过程，提供了一种特定的边界条件。该模型特别适用于内燃机的模拟。Lagrangian 壁面薄膜模型采用瞬态 Lagrangian 粒子来离散物理效应，因此非常适合快速瞬态效应。此外，还存在一种能够处理稳定壁面薄膜的 Eulerian 壁面薄膜模型（Eulerian 壁面薄膜（第837页））。

喷雾与壁面相互作用是进气道燃油喷射（PFI）发动机混合物形成过程中的重要环节。燃油喷雾撞击表面，通常在进气门附近的进气道以及进气门本身，液滴在此处溅射并随后蒸发。蒸发的混合物被卷入发动机气缸，与新鲜充气和气缸内的任何残留气体混合。被压缩和燃烧的混合物最终通过排气口排出。这一过程每分钟重复200至8000次，具体取决于发动机。

由于壁面上形成了油膜，进气道中仍残留着几个周期的燃油。这使得油膜对于进气道燃油喷射（PFI）发动机的碳氢化合物排放变得重要起来。此外，直喷式（DI）发动机的燃烧室内也可能形成油膜。在直喷发动机中，燃油直接喷入燃烧室，如果喷射事件发生在循环的早期或晚期，喷雾可能会撞击活塞。直喷发动机内壁面油膜的建模，尤其是在柴油发动机中，由于发动机表面存在积碳而变得复杂。这些积碳在撞击时会吸收液态油膜。人们认为积碳会在循环后期吸附燃油，然而这一现象非常复杂，尚未被充分理解。

此外，壁面油膜的形成是选择性催化还原（SCR）后处理过程的一个重要因素。利用氨对氮氧化物进行选择性催化还原是一项成熟的技术，用于满足全球严格的氮氧化物（NOx）排放标准。SCR系统设计者面临多项技术挑战，如氨滑脱、尿素结晶、低温沉积以及其他潜在问题。如果不将壁面油膜形成的物理过程纳入SCR预测模型，可能会在NOx转化预测中引入较大误差。

最后，壁面油膜建模在燃气轮机运行和电子设备冷却中也扮演着重要角色。

DPM粒子用于模拟壁面薄膜。Ansys Fluent中的壁面薄膜模型允许液滴撞击边界表面并形成薄膜。该模型可以细分为四个主要子主题：与壁面边界初始撞击时的相互作用、随后在表面上的追踪、薄膜变量的计算以及与气相的耦合。图12.7：拉格朗日壁面薄膜的溅射、动量、热和质量传递机制（第498页）示意性地展示了拉格朗日壁面薄膜模型考虑的基本机制。

图12.7：拉格朗日壁面薄膜的溅射、动量、热和质量传递机制

拉格朗日壁面薄膜模型的主要假设和限制如下：

• 薄膜层很薄，厚度小于500微米。这一限制是由于薄膜中线性速度分布的假设。

• 薄膜粒子的温度变化相对缓慢，因此可以利用解析积分方案。

• 假设薄膜粒子直接接触壁面表面，壁面到薄膜的热传递通过传导进行。

• 薄膜温度永远不会超过液体的沸腾温度。

如果你想模拟喷雾撞击非常热的壁面，壁面射流模型可能更合适，因为在壁面射流撞击模型中，假设液滴下方存在一层蒸汽层，使它们无法直接接触边界表面。对于典型操作条件下的柴油发动机计算，这可能是一个更准确的假设。

有关在Ansys Fluent中使用拉格朗日壁面薄膜模型的更多信息，请参阅《Fluent用户指南》中的“拉格朗日壁面薄膜建模”。

12.9.2 莱顿弗罗斯特温度考虑因素

在拉格朗日壁膜（LWF）模型中，壁膜蒸发速率方程（扩散控制模型的方程12.275（第519页），对流-扩散控制模型的方程12.282（第520页），以及壁膜边界层模型的方程12.291（第523页）和方程12.292（第523页））假设液体膜完全接触壁面。当壁膜温度达到沸点时，计算沸点速率；然而，假设仍然是液体膜完全覆盖壁面。

LWF沸点速率方程（恒温壁和热流壁分别为方程12.283（第520页）和方程12.284（第521页））给出了沸点速率对壁温的线性依赖关系。然而，液体膜沸点现象要复杂得多，随着壁面与液体膜之间的温度差增加，从热壁到膜的热传递速率偏离线性依赖关系。这种偏离是由在热壁表面和沸点膜之间形成的蒸汽口袋引起的，这显著影响了膜沸点过程的热传递和流体动力学特性。绝热蒸汽层完全覆盖热壁的极限条件被称为莱顿弗罗斯特点。

Fluent在多相模型背景下处理壁面沸点现象，特别是在壁面沸点模型（第733页）中。然而，在LWF模型的背景下，预测莱顿弗罗斯特点对于评估LWF沸点速率的充分性并避免非物理结果至关重要。

12.9.2.1 默认壁温限制器

在Ansys Fluent中，默认的壁面温度限制$T_{w\_lim}$是液膜沸点$T_{\text{boil }}$与壁面和液膜之间温度差$\Delta T$之和：

$$T_{w\_\text{lim }}=\left(T_{\text{boil }}+\Delta T\right)\text{\hspace{1em}(12.205)}$$

默认的 $\Delta T$ 为 $50\mathrm{K}$。

12.9.2.2 莱顿弗罗斯特温度报告

在进行LWF计算过程中，Ansys Fluent会计算莱顿弗罗斯特温度并将其与壁面温度进行比较。如果壁面温度高于莱顿弗罗斯特温度，Fluent会发出警告。需要注意的是，莱顿弗罗斯特点的计算仅用于信息提示，并不用于限制计算中的壁面温度。

对于等温表面，莱顿弗罗斯特温度 $T_{\text{Leid, iso }}$ 可以近似为液体临界温度 $T_{crit}$ 和液体表面张力 $\sigma$ 的函数 [51]（第1059页）：

$$T_{Leid,iso}=\frac{27}{32}{T}_{crit}\left(1-\exp \left(-0.52{\left[\frac{{10}^4{\left(\rho /A\right)}^{4/3}}{\sigma }\right]}^{1/3}\right)\right)\text{\hspace{1em}(12.206)}$$

其中，$\rho$ 表示固体的密度，$A$ 是固体的原子序数。方程 12.206（第 500 页）中的最后一项引入了对固-液表面张力效应的修正。

了解 $T_{\text{Leid, iso }}$ 以及固体表面的传热特性后，莱顿弗罗斯特温度 $T_{\text{Leid }}$ 可按以下公式计算：

$$\frac{T_{Leid,iso}-T_f}{T_{Leid}-T_f}=\exp \left(0.00175\beta \right)\mathrm{erfc}\left(0.042\sqrt{\beta }\right)\text{\hspace{1em}(12.207)}$$

其中，$T_f$ 表示液膜温度（K），erfc 是互补误差函数，

$$\beta =\frac{1}{k\rho c_p}$$

其中：

$k=$ 固体导热系数 $\left(\mathrm{W}/\mathrm{m}-\mathrm{K}\right)$

$\rho =$ 固体密度 $\left(\mathrm{kg}/{\mathrm{m}}^3\right)$

$c_p=$ 固体比热容 $\left(\mathrm{J}/\mathrm{kg}-\mathrm{K}\right)$

12.9.3 边界撞击期间的相互作用

预测液滴与墙壁的撞击并不直接，因为涉及大量参数。任何撞击的结果都由液滴液体（粘度、表面张力、密度）和撞击条件（液滴直径和液滴速度）的性质决定。另外两个对撞击结果有重要贡献的因素是表面粗糙度和表面液膜的厚度。此外，墙壁温度也可能影响撞击动力学。

有几种模型可以描述撞击物理过程并预测撞击结果。Ansys Fluent 提供了以下模型：

• Stanton 和 Rutland 模型 [629] (p. 1093), [631] (p. 1093) 使用墙壁温度和基于韦伯数的准则来定义撞击状态（粘附、反弹、扩散和飞溅）之间的固定值阈值。在 Fluent 实现中，原始的 Stanton-Rutland 模型已通过 O'Rourke 和 Amsden [495] (p. 1086) 的概念进行了扩展。更多详情请参见 Stanton-Rutland 模型 (p. 501)。

• Kuhnke 模型 [324] (p. 1076) 考虑了干燥和湿润的墙壁，并额外考虑了液滴属性和撞击条件、表面粗糙度和膜厚度的影响。Kuhnke 模型区分以下状态：反弹、扩散、飞溅和干燥飞溅（热破碎）。更多信息请参见 Kuhnke 模型 (p. 508)。

• Stochastic Kuhnke 模型 [606] (p. 1092), [37] (p. 1059) 源自 Kuhnke 模型，具有以下特点：

• 使用为选择性催化还原（SCR）建模应用开发和调整的不同状态转换准则。

• 在关键温度转换过程中引入随机效应。

• 引入了蒸发溅射机制，采用“部分蒸发”概念，即在特定条件下假设部分撞击液滴完全蒸发。

详情参见随机Kuhnke模型（第514页）。

• 用户自定义的撞击模型可通过使用DEFINE_IMPINGEMENT、DEFINE_FILM_REGIME和DEFINE_SPLASHING_DISTRIBUTION用户自定义函数来实现，具体描述见Fluent用户指南和Fluent定制手册中的用户自定义函数部分。

注意：

• 对于能量方程未求解的情况，拉格朗日壁面薄膜模型适用，其撞击/溅射模型基于以下假设：

• 壁面温度为$300\mathrm{K}$

• 液膜的饱和温度等于$5000\mathrm{K}$

• 后续章节中提及的颗粒速度相对于壁面坐标系。除非另有说明，无量纲参数$Re$和$We$基于垂直于壁面的颗粒速度。

12.9.3.1 Stanton-Rutland模型

12.9.3.1.1 机制定义

壁面相互作用基于Stanton [629]（第1093页）和O’Rourke [496]（第1086页）的工作，根据局部信息计算液滴与壁面相互作用的机制。四种机制（粘附、反弹、扩散和溅射）基于撞击能量和壁面温度。以下图表有助于展示界限。

图12.8：壁面相互作用准则的简化决策图

低于液体的临界转变温度时，撞击液滴可能粘附、扩散或溅射；高于临界转变温度时，颗粒可能反弹或溅射。

机制划分的标准基于撞击能量和液体的临界转变温度。撞击能量定义为...

$$E^2=\frac{\rho V_{pn}^2{d}_p}{\sigma }\left(\frac{1}{\min \left(h_0/d_p,1\right)+{\delta }_{bl}/d_p}\right)\text{\hspace{1em}(12.208)}$$

其中，$V_{pn}$ 表示粒子垂直于壁面的速度，$\rho$ 是液体密度，$d_p$ 是液滴直径，$\sigma$ 是液体表面张力，$h_0$ 是液膜高度。这里，${\delta }_{bl}$ 表示边界层厚度，其定义为

$${\delta }_{bl}=\frac{d_p}{\sqrt{Re}}\text{\hspace{1em}(12.209)}$$

雷诺数定义为 $Re=\rho V_{\mathrm{pn}}{d}_p/\mu$。通过如第12.208式（第502页）那样定义能量，壁面上薄膜的存在抑制了飞溅现象，但当薄膜高度趋近于零时，并不会产生非物理结果。

临界过渡温度定义为：

$${\mathrm{T}}_{\text{crit }}={\mathrm{T}}_{\text{crit }}^{\ast }{\mathrm{T}}_{\text{sat }}\text{\hspace{1em}(12.210)}$$

其中，${\mathrm{T}}_{\text{sat }}$ 表示饱和温度，而 ${\mathrm{T}}_{\text{crit }}^{\ast }$ 是临界温度因子，其值根据应用场景在1至1.5之间变化。${\mathrm{T}}_{\text{crit }}^{\ast }$ 的默认值为1；因此，默认的临界转变温度 ${\mathrm{T}}_{\text{crit }}$ 将等于液滴的饱和温度 ${\mathrm{T}}_{\text{sat }}$。

作为拉格朗日壁面薄膜模型的替代方案，临界转变温度 $T_{\text{crit }}$ 可以通过在液体饱和温度 $T_{sat}$ 上增加一个恒定偏移量 $T_{cal}$ 来计算：

$$T_{\text{crit }}=T_{\text{sat }}+T_{\text{cal }}\text{\hspace{1em}(12.211)}$$

当无量纲能量 $E$ 小于 16 时，应用附着状态，此时颗粒速度设定为壁面速度。在扩散状态中，颗粒的初始方向和速度根据壁面射流模型设定，其中液滴沿表面特定方向的概率通过与无粘性液体射流的类比给出，该射流的动量通量具有经验定义的径向依赖性。

如果壁面温度高于液体的临界转变温度，低于临界冲击能量 $\left(E_{cr}\right)$ 的冲击事件会导致颗粒从壁面反弹。

当冲击能量超过一个临界能量阈值，定义为 $E_{cr}$ =57.7 时，会发生溅射。溅射液滴包的数量在 Wall 边界条件对话框中设定，默认数量为 4，但您可以选择 0 到 10 之间的数字。溅射算法遵循 Stanton [629]（第 1093 页）的描述，并在溅射（第 503 页）中详细说明。

12.9.3.1.2 反弹

在反弹状态中，反弹颗粒的速度由以下公式给出：

$$V_{rt}=\frac{5}{7}{V}_{pt}\text{\hspace{1em}(12.212)}$$

$$V_{rn}=e_n{V}_{pn}\text{\hspace{1em}(12.213)}$$

其中，$V_{\mathrm{rt}}$ 和 $V_{\mathrm{rn}}$ 分别是反弹速度的切向和法向分量，$V_{\mathrm{pt}}$ 是冲击速度的切向分量，而 $e_n$ 则是由以下公式给出的法向恢复系数：

$$e_n=0.993-1.76{\Theta }_I+1.56{\Theta }_I^2-0.49{\Theta }_I^3\text{\hspace{1em}(12.214)}$$

${\Theta }_I$ 是从壁面测量的冲击角度，以弧度为单位。

反弹粒子的偏转角 $\psi$ 在入射方向的 -90° 到 90° 之间随机取值。

12.9.3.1.3 飞溅

图 12.9：Stanton-Rutland 模型：冲击与飞溅（第 504 页）图形化地描述了 Stanton-Rutland 飞溅模型中的以下飞溅液滴变量。

• 冲击粒子的直径，$d_p$

• 冲击粒子的速度，$V_p$ 及其变化性

• 飞溅的二次粒子的直径，$d_i$

• 飞溅的二次粒子的速度，$V_i$

• 从壁面法线测量的冲击角度 ${\Theta }_I$

• 从壁面法线测量的飞溅粒子的反射角 ${\theta }_s$

• 从入射方向的偏转角 $\psi$

图 12.9：Stanton-Rutland 模型：冲击与飞溅

如果冲击表面的粒子具有足够高的能量，粒子会发生飞溅并产生多个新粒子。您可以在 DPM 标签或 Wall Film 标签的 Wall 边界条件对话框中明确设置每次冲击产生的粒子数量。飞溅包的最小数量为三个。飞溅包的属性（直径、大小和方向）从以下章节描述的实验获得的分布函数中随机抽样。在拉格朗日壁膜模型中将飞溅包的数量设置为零会关闭拉格朗日模型的飞溅计算。请记住，每个飞溅包可以被视为分布曲线的离散样本，并且在 Wall 边界条件对话框中选择飞溅液滴的数量并不限制飞溅液滴的数量，只是代表这些液滴的包的数量。

因此，对于每个溅射的包裹，通过从累积概率分布函数（CPDF）中采样，也称为F，该函数来自威布尔分布函数，并根据Mundo等人[465]（第1084页）的数据进行拟合，从而获得不同的直径。方程为

$$f\left(d\right)=2\frac{d}{D^2}\exp \left[-{\left(\frac{d}{D}\right)}^2\right]\text{\hspace{1em}(12.215)}$$

该分布代表了在飞溅液滴样本中找到直径为 $d_i$ 的液滴的概率。此分布与 O’Rourke [496]（第1086页）所使用的 Nakamura-Tanasawa 分布函数相似，其分布参数 $D=\sqrt{2}{d}_{\max }$。为确保分布函数在韦伯数增加时产生符合物理规律的结果，采用了 O’Rourke [496]（第1086页）提供的以下 $d_{\max }$ 表达式。飞溅尺寸分布峰值处的粒径与撞击粒径的比值由下式给出：

$$d_{\max }/d_p=\max \left(\frac{E_{\text{crit }}^2}{E^2},\frac{6.4}{We},0.06\right)\text{\hspace{1em}(12.216)}$$

能量表达式由方程12.208（第502页）给出。低韦伯数冲击由方程12.216（第504页）中的第二项描述，并且在O’Rourke [496]（第1086页）分析的任何实验中，对于极高能量冲击，飞溅峰值直径比绝不会小于0.06。方程12.216（第504页）中的韦伯数是使用粒子法向速度和直径定义的：

$$We=\frac{\rho V_{pn}^2{d}_p}{\sigma }\text{\hspace{1em}(12.217)}$$

累积概率分布函数（CPDF）是必需的，以便能够从实验数据中采样直径。通过积分公式12.215（第504页）可以得到CPDF。

$$F\left(d\right)=1-\exp \left[-{\left(\frac{d}{D}\right)}^2\right]\text{\hspace{1em}(12.218)}$$

其范围在零和一之间。通过反转公式12.218（第505页）并在零和一之间采样CPDF，得到了直径的表达式（该直径是撞击韦伯数$We$和撞击能量的函数$D$）。因此，第$i$个溅射液滴包的直径表达式为$d_i=D\sqrt{-\ln \left(1-c_i\right)}$，其中$c_i$是第$i$个随机样本。一旦确定了溅射液滴的直径，通过在给定直径处评估公式12.215（第504页），就可以确定在该样本中找到该液滴的概率。

为了防止溅射颗粒直径变得过大，Ansys Fluent采用了分布限制方法：

• 尾部限制方法（拉格朗日壁膜模型的默认设置）限制了

$$d_{99.9}/d_p<1.0$$

其中，$d_{99.9}$ 表示粒径，对应于累积概率密度函数（CPDF）的 $99.9\%$ 处。

• 峰值限制法（欧拉壁膜模型的默认设置）限制了粒径比。

$$d_{\max }/d_p<P_{\max }$$

其中，$P_{\max }$ 是一个限制值，根据 Mundo [465]（第1084页）从实验数据中确定的范围为0.2到0.3（默认值为0.22）。

可以通过文本用户界面（TUI）命令选择限制方法和参数 $P_{\max }$。

每个包裹中的液滴数量与该液滴大小的PDF值成正比。然后，总粒子数由溅射的总质量决定。

溅射粒子的质量分数 $y_s$ 是从Yarin和Weiss [724]（第1099页）的实验数据中获得的。O’Rourke等人 [496]（第1086页）给出的溅射质量分数表达式为：

$$y_s=\{\begin{array}{l}1.8\times {10}^{-4}\left(E^2-E_{crit}^2\right),E_{crit}^2<E^2<7500\\ 0.75,E^2>7500\end{array}\text{\hspace{1em}(12.219)}$$

Stanton等人[629]（第1093页）提出的飞溅质量分数表达式如下：

$$y_s=\{\begin{array}{ll}0.0& K_u\le 18.0\\ -27.2+3.15K_u-0.116K_u^2+1.40\times {10}^{-3}{K}_u^3& 18.0<K_u<33.0\\ 0.75& K_u\ge 33.0\end{array}\text{\hspace{1em}(12.220)}$$

其中，$K_u$ 是无量纲冲击速度，定义为：

$$K_u=G\sqrt{K}\text{\hspace{1em}(12.221)}$$

其中，$K$ 由公式 12.236（第 508 页）定义，常数 $G$ 取值为 2.37。

飞溅质量分数方法可通过 TUI 命令选择，对于拉格朗日壁面薄膜（LWF），默认采用 O’Rourke 表达式；对于欧拉壁面薄膜（EWF），默认采用 Stanton 表达式。

作者（O’Rourke 等人[496]（第 1086 页））指出，对于典型的柴油喷雾，几乎所有撞击都远高于上限，因此飞溅事件几乎总是会弹出撞击液滴质量的 75%。为了得到液滴总数量的表达式，我们注意到质量守恒要求飞溅包裹的总质量之和必须等于飞溅质量分数，即

$$\frac{\rho \pi }{6}{N}_{\text{tot }}\sum _{i=1}^{N_{\text{parcels }}}\left(f_i{d}_i^3\right)=y_s{m}_p\text{\hspace{1em}(12.222)}$$

其中，$m_p$表示撞击气块的总质量。总溅射液滴数量的表达式为

$$N_{\text{tot }}=y_s{m}_p/\left(\frac{\rho \pi }{6}\sum _{i=1}^{N_{\text{parcels }}}\left(f_i{d}_i^3\right)\right)\text{\hspace{1em}(12.223)}$$

每个包裹溅射出的液滴数量由 $N_{tot}$ 乘以该液滴尺寸对应的概率密度函数（PDF）值来确定。

为了计算溅射液滴离开表面的速度，还需对速度的法向分量进行额外的相关性采样。采用拟合Mundo [465]（第1084页）数据的高斯威布尔函数作为法向分量的概率密度函数（PDF）。概率密度表示为

$$f\left(\frac{V_{in}}{V_{pn}}\right)=\left[\frac{b_v}{{\Theta }_v}{\left(\frac{V_{in}/V_{pn}}{{\Theta }_v}\right)}^{b_v-1}\right]\exp \left[-{\left(\frac{V_{in}/V_{pn}}{{\Theta }_v}\right)}^{b_v}\right]\text{\hspace{1em}(12.224)}$$

其中，索引 $i$ 指的是第 $i^{\text{th}}$ 个溅射的液滴，而下标 $p$ 指的是撞击的液滴。$b_v$ 和 ${\theta }_v$ 是 ${\theta }_l$ 的函数，${\theta }_l$ 是从表面法线测量的撞击角度：

$$b_v=\{\begin{array}{cc}2.1,& {\theta }_I\le {50}^{\circ }\\ 1.10+0.02{\theta }_I,& {\theta }_I>{50}^{\circ }\end{array}$$

$${\theta }_v=0.158e^{0.017{\theta }_I}$$

速度的切向分量是通过反射角表达式 ${\Theta }_s$ 获得的：

$${\theta }_s=65.4+0.226{\theta }_l\text{\hspace{1em}(12.225)}$$

结合

$$V_{it}=V_{in}\tan \left({\Theta }_s\right)\text{\hspace{1em}(12.226)}$$

请注意，反射角 ${\theta }_s$ 是从表面法线测量的，并且 ${\theta }_s$ 和 ${\theta }_l$ 都以度为单位。

根据 [630]（第 1093 页）的计算，飞溅液滴的方位角 $\psi$ 假设遵循分布 $H\left(\psi \right)$（公式 12.227，第 507 页）。该处理基于 Naber 和 Reitz [468]（第 1084 页）开发的方法。

$$H\left(\psi \right)=\{\begin{array}{ll}1/2\pi & \beta =0\\ \frac{\beta /\pi }{2\left(e^{\beta }-1\right)}\exp \left[\beta \left(1-\left(\left|\psi \right|/\pi \right)\right)\right]& \beta >0\end{array}\text{\hspace{1em}(12.227)}$$

其中，$\beta$ 是一个与冲击角度 ${\theta }_I$ 相关的参数，其关系为：

$$\sin {\theta }_I=\frac{e^{\beta }+1}{e^{\beta }-1}\frac{{\beta }^2}{{\beta }^2+{\pi }^2}\text{\hspace{1em}(12.228)}$$

$\psi$ 是通过对公式12.227（第507页）和公式12.228（第507页）进行积分来确定的，具体如下：

• 对于 $\beta =0$ 的情况：

$$\psi =2\pi p\text{\hspace{1em}(12.229)}$$

对于 $\beta >0$ 的情况：

$$\psi =\{\begin{array}{ll}\left(\pi /\beta \right)\log \left[e^{-\beta }+2p\left(1-e^{-\beta }\right)\right]& \text{ for }p<0.5\\ \left(\pi /\beta \right)\log \left[2-e^{-\beta }-2p\left(1-e^{-\beta }\right)\right]& \text{ for }p>0.5\end{array}\text{\hspace{1em}(12.230)}$$

在上式中，$p$ 是一个在0到1之间均匀分布的随机数。

$$p=\mathrm{rnd}\left(0,1\right)$$

最后，对新形成的微粒进行能量平衡计算，以确保新水滴的动能和表面能之和不超过旧水滴的能量总和。能量平衡的计算公式如下：

$$\frac{1}{2}\sum _{i=1}^{N_{\text{parcel }}}\left(m_i{V}_i^2\right)+\pi \sigma \sum _{i=1}^{N_{\text{parcel }}}\left(N_i{d}_i^2\right)=\frac{1}{2}{m}_p{V}_{pn}^2+\pi \sigma \left(N_p{d}_p^2\right)-E_s\text{\hspace{1em}(12.231)}$$

其中，$E_s$ 是发生溅射的阈值能量。为了确保能量守恒，计算以下修正因子：

$$E_k=\frac{\frac{1}{2}{m}_p{V}_{\mathrm{pn}}^2+\pi \sigma \left(N_p{d}_p^2\right)-E_s-\pi \sigma \sum _{i=1}^{N_{\text{parcel }}}\left(N_i{d}_i^2\right)}{\frac{1}{2}\sum _{i=1}^{N_{\text{parcel }}}\left(m_i{V}_i^2\right)}\text{\hspace{1em}(12.232)}$$

阈值能量 $E_s$ 计算如下：

$$E_s=\frac{1}{2}{m}_p{V}_s^2\text{\hspace{1em}(12.233)}$$

其中，$V_s$ 表示从临界韦伯数计算得出的溅射阈值速度。

$$V_s=\sqrt{W{e}_{cr}\frac{\sigma }{\rho d_p}}\text{\hspace{1em}(12.234)}$$

临界韦伯数则由以下公式给出：

$$W{e}_{cr}={18}^2{d}_p{\left(\frac{\rho }{\sigma }\right)}^{1/2}{\left(\frac{\mu }{\rho }\right)}^{1/4}{\left(\frac{V_{pn}}{d_p}\right)}^{3/4}\text{\hspace{1em}(12.235)}$$

由于飞溅液滴速度的采样点数量相对较少，因此需要这一修正因子（详见Stanton [631]（第1093页））。在方程12.232（第507页）中，飞溅包的各组成部分乘以$E_k$的平方根，以确保能量守恒。因此，飞溅包的法向和切向速度分量可表示为：

$$V_{in}^{\prime }=\sqrt{E_K}{V}_{in}\text{ and }{V}_{it}^{\prime }=\sqrt{E_K}{V}_{it}$$

需要注意的是，飞溅事件本质上是瞬态的，因此飞溅子模型仅在 Ansys Fluent 中进行非稳态追踪时可用。飞溅还可能导致发生飞溅的单元中的源项大幅增加，这可能会导致时间步长之间收敛困难。因此，在启用飞溅功能进行模拟时，可能需要使用更小的时间步长。

12.9.3.2 Kuhnke 模型

12.9.3.2.1 状态定义

图 12.10：Kuhnke 冲击模型状态（第 508 页）展示了 Kuhnke [324]（第 1076 页）冲击模型的冲击图。该模型通过基于无量纲变量 $K$ 和 $T^{\ast }$ 将相关冲击现象分类为四种状态，考虑了所有相关的冲击现象：扩散（或沉积）、反弹、飞溅和干燥飞溅（热破碎）。沉积和飞溅状态会导致壁膜的形成。

图 12.10：Kuhnke 冲击模型状态

无量纲变量 $K$ 定义为：

$$K=W{e}^{5/8}L{a}^{1/8}\text{\hspace{1em}(12.236)}$$

其中：

$$We=\frac{\rho V_{pn}^2{d}_p}{\sigma }\text{\hspace{1em}(12.237)}$$

$$La=\frac{\sigma \rho d_p}{{\mu }^2}\text{\hspace{1em}(12.239)}$$

无量纲变量 $T^{\ast }$ 定义如下：

$$T^{\ast }=\frac{T_w}{T_{\text{sat }}}\text{\hspace{1em}(12.239)}$$

在上述方程中，使用了以下符号表示：

$We=$ 基于法向冲击速度的韦伯数

$La=$ 拉普拉斯数

$\rho =$ 液滴密度 $\left(\mathrm{kg}/{\mathrm{m}}^3\right)$

$d_p=$ 液滴直径

$V_{pn}=$ 法向冲击速度 $\left(\mathrm{m}/\mathrm{s}\right)$

$\sigma =$ 液滴表面张力 $\left(\mathrm{N}/\mathrm{m}\right)$

$\mu =$ 液滴粘度 $\left(\mathrm{kg}/\mathrm{m}-\mathrm{s}\right)$

$T_w=$ 壁面温度 (K)

$T_{\text{sat }}=$ 液滴饱和温度 (K)

当壁面温度 $T_w$ 达到临界过渡温度 $T_{\text{crit }}$，标志着膜沉积的开始，使用方程 12.239 (第 509 页) 来定义临界温度因子 $T_{\text{crit }}^{\ast }$：

$$T_{\text{crit }}^{\ast }=T_{\text{crit }}/T_{\text{sat }}\text{\hspace{1em}(12.240)}$$

更具体地说，Kuhnke模型考虑了以下几种情况：

• 扩散：在以下任一条件下，将撞击的液滴加入壁面液膜中：

$$\text{-}K<K_{\text{crit }}\text{and}{T}^{\ast }<T_{\text{crit }}^{\ast }$$

• 当$K<K_{\text{crit }}$且墙面湿润时

• 回弹：在以下任一条件下，撞击的液滴会原封不动地反弹，尺寸不变：

$$-K<K_{\text{crit }}\text{and}{T}^{\ast }>T_{\text{crit }}^{\ast }$$

• $我们<5$ 且墙面潮湿

• 飞溅效果：若满足以下条件，则将水滴雾化为更小的水滴：

$$-K>K_{\text{crit }}$$

在飞溅状态下，根据壁面温度和壁面表面条件，部分撞击的液滴质量可能会沉积形成液膜（溅射），或者液滴可能被完全雾化（干溅射或热破碎）。

临界转变参数 $K_{\text{crit }}$ 和 $T_{\text{crit }}^{\ast }$ 计算如下。

与壁面温度相关的转变由临界温度因子 $T_{\text{crit }}^{\ast }$ 控制，其值根据应用在 1.1 到 1.5 之间变化。$T_{\text{crit }}^{\ast }$ 的默认值为 1；因此，默认的临界转变温度 ${\mathrm{T}}_{\text{crit }}$ 将等于液滴液体的饱和温度 ${\mathrm{T}}_{\text{sat }}$。Kuhnke [324]（第 1076 页）建议单个液滴撞击时取 1.1，液滴链撞击时取 1.16。对于大多数尿素-水悬浮液，$T_{\text{crit }}^{\ast }$ 设定为 1.4 [66]（第 1060 页）。

作为拉格朗日壁面液膜模型的替代方案，临界转变温度 $T_{\text{crit }}$ 可以通过液体的饱和温度 $T_{\text{sat }}$ 加上一个恒定偏移量 $T_{\text{cal }}$ 来计算：

$$T_{\text{crit }}=T_{\text{sat }}+T_{\text{cal }}\text{\hspace{1em}(12.241)}$$

对于SCR悬浮液，推荐的$T_{cal}$值范围在167至${207}^{\circ }\mathrm{C}$之间。

参数$K$的两个关键值（一个对应于高温下反弹与热破碎之间的转变，即$K_{\text{hot }}$，另一个对应于低温下沉积与溅射之间的转变，即$K_{\text{cold }}$）取决于多个参数和壁面条件（湿或干），具体描述如下。

• 湿壁

$$\begin{array}{rl}{K}_{\text{crit }}\left(\delta ,La\right)=& \phi \left(\delta ,{\delta }_0,125,58.7,-100\right)+\\ & \left(408.4{\delta }_0-283.6\right)\cdot \mathrm{Weibull}\left(\delta ,{\delta }_0,3\right)\end{array}\text{\hspace{1em}(12.242)}$$

无量纲薄膜高度 $\delta$ 的计算公式如下：

$$\delta =h_f/d_p$$

其中，$h_f$ 表示液膜高度 $\left[\mathrm{m}\right]$，$d_p$ 表示液滴直径 $\left[\mathrm{m}\right]$。

公式 12.242（第 510 页）中的其他参数表述如下：

$$\phi \left(x,x_o,y_l,y_r,s\right)=\frac{y_r-y_l}{\pi }\arctan \left(\pi s\frac{x-x_o}{y_r-y_l}\right)+\frac{y_r+y_l}{2}\text{\hspace{1em}(12.243)}$$

$${\delta }_0=\phi \left(La,600,0.85,1.08,0.003\right)\text{\hspace{1em}(12.244)}$$

$$\begin{array}{c}\text{Weibull}\left(x,x_0,q\right)=\frac{q}{x_0}{\left(\frac{x}{x_0}\right)}^{q-1}e–{\left(x/x_0\right)}^q\end{array}\text{\hspace{1em}(12.245)}$$

这里，$x$、$x_o$、$y_{l^{\prime }}$、$y_{r^{\prime }}$、$s$ 和 $q$ 是常数，其值如公式 12.242（第510页）所示。

• 干墙

$$K_{\text{crit }}\left(T^{\ast },{\epsilon }_a\right)=\left(1-\lambda \right)K_{\text{cold }}+\lambda K_{\text{hot }}\text{\hspace{1em}(12.246)}$$

其中

${\epsilon }_a=$ 无量纲表面粗糙度 $R_a/d_p$

$R_a=$ 表面粗糙度 (m) （平均绝对距离与均值的比值）

以及

$$\lambda \left(T^{\ast }\right)=\arctan \left(p\left(T^{\ast }-T_{\text{crit }}^{\ast }\right)\right)/\pi +1/2\text{\hspace{1em}(12.247)}$$

将下面的文本翻译成中文：与

$$p=\tan \left(\pi \left(\omega -0.5\right)/2\right)/d{T}^{\ast }$$

其中，$\omega$ 是混合函数的阈值转换值，设定为 0.99，而 $d{T}^{\ast }$ 则代表混合区域的宽度。

参数 $K_{\text{hot }}$ 在 20 至 40 之间随机取样。

$$K_{\text{hot }}=20+20p\text{\hspace{1em}(12.248)}$$

其中，$p$ 是一个在 0 到 1 之间均匀分布的随机数：

$$p=\mathrm{rnd}\left(0,1\right)$$

参数 $K_{\text{cold}}$ 是通过以下方程计算的：

$$K_{\text{cold }}\left(T^{\ast },{\epsilon }_a\right)=K_{\min }\left(T^{\ast }\right)-12.75\left(\ln \left({\epsilon }_a\right)-\chi \right)+12.75\sqrt{{\left(\ln \left({\epsilon }_a\right)-\chi \right)}^2+5}\text{\hspace{1em}(12.249)}$$

将下面的文本翻译成中文：随着

$$\chi =-\frac{K_{\min }\left(T^{\ast }\right)+43.6}{25.5}$$

$$K_{\min }\left(T^{\ast }\right)=\{\begin{array}{ll}{K}_{T_{\min }}+\left(K_{T_{\max }}-K_{T_{\min }}\right)e^{c_T\left(T^{\ast }-1\right)}& T^{\ast }\le 1\\ K_{T_{\max }}& 1<T^{\ast }<T_{\text{trans }}^{\ast }\end{array}\text{\hspace{1em}(12.250)}$$

式中

$$K_{T_{\min }}=54$$

$$K_{T_{\max }}=130$$

$$c_T=13$$

飞溅质量分数确定如下：

$$y_s=\min \left[1,\left(T^{\ast }-0.8\right)\left(1-B\right)/\left(T_{\text{crit }}^{\ast }-0.8\right)+B\right]\text{\hspace{1em}(12.251)}$$

式中

$$B=\{\begin{array}{ll}0.2+0.6p& \text{ for dry wall }\\ 0.2+0.9p& \text{ for wet wall }\end{array}\text{\hspace{1em}(12.252)}$$

12.9.3.2.2 反弹

在反弹状态下，反弹粒子的速度由以下公式给出：

其中，$p=\mathrm{rnd}\left(0,1\right)$

$$V_{rt}=V_{pt}\text{\hspace{1em}(12.253)}$$

$$V_{rn}=e_n{V}_{pn}\text{\hspace{1em}(12.254)}$$

其中，$V_{rt}$ 和 $V_{rn}$ 分别表示反弹速度的切向和法向分量，$V_{pt}$ 是入射速度的切向分量，而 $e_n$ 则是通过以下公式计算得到的法向恢复系数：

$$e_n=\sqrt{0.678e^{-{0.04415We}_0}}\text{\hspace{1em}(12.255)}$$

其中，$W{e}_0$ 是冲击前的正常冲击韦伯数。

反弹颗粒的偏差角 $\psi$，从入射方向随机取值在 -90° 到 90° 之间。

12.9.3.2.3 飞溅

图 12.11：Kuhnke 模型：冲击液滴（第 513 页）图形化地描述了 Kuhnke 飞溅模型中的以下飞溅液滴变量：

• 冲击颗粒的直径，$d_p$

• 冲击颗粒的速度，$V_p$ 加上变异性

• 飞溅的次级颗粒的直径，$d_i$

• 飞溅颗粒的平均速度，$V_i$ 加上变异性

• 冲击角度 ${\Theta }_I$，从壁面测量

• 飞溅颗粒的反射角度 ${\theta }_s$，从壁面测量

• 从入射方向的偏差角 $\psi$ 注意：

注意，通常情况下，冲击和反射角度在 Stanton/Rutland 和 Kuhnke 模型框架中测量方式不同：

• Stanton 和 Rutland 模型：从壁面的法线测量

• Kuhnke 模型：从壁面测量

图 12.11：Kuhnke 模型：冲击液滴

飞溅颗粒的属性计算如下。

首先，飞溅液滴的平均直径计算如下：

$$d_{10}={\gamma }_{10}{d}_{\mathrm{p}}\text{\hspace{1em}(12.256)}$$

将下面的文本翻译成中文：随着

$${\gamma }_{10}=\{\begin{array}{ll}3.3e3.6{\left({\theta }_l/\pi \right)}^{2}W{e}^{-0.65}& \text{ dry wall }\\ 2.2e3.6{\left({\theta }_l/\pi \right)}^{2}W{e}^{-0.36}& \text{ wet wall }\end{array}\text{\hspace{1em}(12.257)}$$

其中

$d_{\mathrm{p}}=$主液滴直径 $\left[\mathrm{m}\right]$

$We=$冲击颗粒的正常冲击韦伯数

${\theta }_I=$冲击角度，单位为弧度

接下来，根据Stanton和Rutland [629]（第1093页），[631]（第1093页）的原始提议，飞溅液滴直径 $d_{si}$ 从连续Weibull概率密度函数${weibull}_{cpdfi}$ 的逆函数中随机采样，该函数描述如下：

$${\text{ weibull }}_{\text{cpdfi }}\left(d_i,d_{10},q\right)=d_{10}{\left(-\log \left(1-d_i\right)\right)}^{1/q}\text{\hspace{1em}(12.258)}$$

当 $q=2$ 时，也采用了如《飞溅》（第503页）中所述的极限方法。

每个飞溅包裹中的液滴数量根据总体质量平衡，通过方程12.222（第506页）确定。

平均飞溅液滴速度 $V_1$ 是根据次级韦伯数 $W{e}_{a1}$ 和平均液滴尺寸 $d_{10}$ 计算得出的：

$$V_1=\sqrt{\frac{\sigma W{e}_{a1}}{\rho d_{10}}}\text{\hspace{1em}(12.259)}$$

式中

$$W{e}_{a1}=\{\begin{array}{ll}{\gamma }_{10}\left(W{e}_{a0}\left(1-\eta {\sin }^2{\theta }_I\right)+12\right)-12/v_{32}& \text{ dry wall }\\ \max \left(51-7.1e^{3.4\left({\theta }_I/\pi \right)},c_{a}W{e}_{a0}\right)& \text{ wet wall }\end{array}\text{\hspace{1em}(12.260)}$$

其中，$W{e}_{a0}$ 和 $W{e}_{a1}$ 是基于绝对颗粒速度的，并且

$${\eta }^{\prime }=0.85$$

$$c_a=-0.378{\left({\Theta }_I/\pi \right)}^2-0.123\left({\Theta }_I/\pi \right)+0.156$$

$$v_{32}=2$$

与溅射液滴直径类似，溅射液滴速度 $V_i$ 从连续Weibull概率密度函数的倒数中随机抽样，其分布平均值为 $V_1$。

喷射角度 ${\theta }_S$ 则通过假设一个逻辑分布，从平均值 ${\bar{\theta }}_S$ 计算得出：

$${\Theta }_S=\{\begin{array}{ll}9.3+\ln \left({\epsilon }_z\right)\left(2.7-0.03{\Theta }_I\right)+0.22{\Theta }_I& \text{ cold dry wal }\\ 0.225{\Theta }_{I}e{\left(0.017{\Theta }_I-0.937\right)}^2& \text{ wet wall }\\ 0.96{\Theta }_I{e}^{-0.0045We}& \text{ hot wall }\end{array}\text{\hspace{1em}(12.261)}$$

角度（以度为单位）按照逻辑PDF logistic $\left({\theta }_S,{\bar{\theta }}_S,{\sigma }_l\right)$ 分布，其中 ${\sigma }_l=4$。

此处使用的符号如下：

${\epsilon }_z=$ 无量纲表面粗糙度 $R_z/d_p$

$R_z=$ 表面粗糙度 $\left(\mathrm{m}\right)$（最大高度差的平均值）

连续逻辑PDF的逆函数由以下公式给出：

$${\mathrm{logistic}}_{\text{cpdfi }}\left({\theta }_S,{\bar{\theta }}_S,{\sigma }_l\right)={\sigma }_l\log \left[{\theta }_S/\left(1-{\theta }_S\right)\right]+{\bar{\theta }}_S\text{\hspace{1em}(12.262)}$$

偏差角 $\psi$ 的确定方法如下：

$$\psi =-\frac{\pi }{\beta }\ln \left(1-p\left(1-e^{-\beta }\right)\right),p=\mathrm{rnd}\left(0,1\right)\text{\hspace{1em}(12.263)}$$

式中

$$\beta =\{\begin{array}{ll}\sqrt{\left(1+8.872\cos \left(1.152{\theta }_l\right)\right)/\left(1-\cos {\theta }_l\right)}& {\theta }_l\le {80}^{\circ }\\ \left({\pi }^2/2\right)\cos {\theta }_l& {\theta }_l\ge {80}^{\circ }\end{array}$$

12.9.3.3 随机Kuhnke模型

随机Kuhnke模型的状态图如图12.12所示：随机Kuhnke模型的状态图（第515页）。

图12.12：随机Kuhnke模型的状态图

该模型通过根据临界温度和韦伯数转换标准将所有相关的冲击现象分类为以下四个核心状态，考虑了所有相关的冲击现象：

• 反弹

• 扩散

• 飞溅

• 蒸发飞溅

扩散和飞溅状态导致壁膜的形成。此外，引入了过渡状态以避免在膜沉积过程中发生突变。

状态转换标准定义如下。定义了上下沉积极限$T_{DL,u}$和$T_{DL,l}$之间的临界温度范围：

$$T_{DL,u}=T_{sat}+\Delta T_{DL}\text{\hspace{1em}(12.264)}$$

$$T_{DL,l}=T_{DL,u}-\Delta T_{lim}\text{\hspace{1em}(12.265)}$$

其中，

$T_{\text{sat }}=$ 液膜饱和温度 [K]

$\Delta T_{DL}=$ 上沉积极限偏移 $\left[\mathrm{K}\right]$

$\Delta T_{\text{lim }}=$ 沉积范围温度差 $\left[\mathrm{K}\right]$

韦伯数准则在湿壁和干壁上的定义有所不同：

• 湿壁上溅射的临界韦伯数：

$$W{e}_{\text{crit },w}=A_{w}L{a}^{-0.18}\text{\hspace{1em}(12.266)}$$

• 干壁溅射临界韦伯数：

$$W{e}_{crit,D}=A\cdot L{a}^{-0.18}\text{\hspace{1em}(12.267)}$$

在公式12.266（第515页）和公式12.267（第516页）中：

$La=$ 拉普拉斯数

$A_w=$ 拉普拉斯数常数

$A=$ 根据表12.2插值得到的壁面粗糙度$r_z$的函数：参数$A$作为壁面粗糙度的函数（第516页）

表12.2：参数$A$作为壁面粗糙度的函数

拉普拉斯数表示为：

$$La=\rho \sigma d/{\mu }^2\text{\hspace{1em}(12.268)}$$

其中，

$\rho =$ 液滴密度

$\sigma =$ 液滴表面张力

$d=$ 液滴直径

$\mu =$ 液滴粘度

最后，Stochastic Kuhnke 模型根据液滴的韦伯数（参见 Ansys Fluent 理论指南中“Kuhnke 模型”部分的韦伯数定义）和壁面条件，为每次撞击事件分配适当的模式。这些模式包括：

• 反弹

$$\text{Regime boundaries:}We<5\text{and}\{\begin{array}{l}{T}_{wall}<T_{DL,l}\text{ and Wet wall }\\ \text{ or }\\ T_{wall}>T_{DL,l}\end{array}$$

冲击后的液滴在大小不变的情况下被喷射出去。根据Ansys Fluent理论指南中关于Kuhnke模型的“反弹”部分，计算了反弹颗粒的速度分量。

• 扩展：

在流体动力学模拟中，液滴的冲击与反弹行为是研究的重要环节。本段描述了液滴在受到冲击后，其尺寸保持不变，并被重新发射出去的现象。这一过程的速度分量计算，是依据Ansys Fluent软件中的Kuhnke模型进行的，该模型详细规定了液滴反弹时的物理行为，确保了模拟结果的准确性与可靠性。通过这样的模拟，科研人员能够更好地理解液滴在复杂环境中的动态响应，进而为相关工程应用提供理论支持。

$$\text{Regime boundaries:}{T}_{\text{wall }}<T_{DL,l}\text{and}\{\begin{array}{l}5<\mathrm{We}<{\mathrm{We}}_{\mathrm{SW}}\text{ and Wet wall }\\ \text{ or }\\ \mathrm{We}<{\mathrm{We}}_{\mathrm{SD}}\text{ and Dry wall }\end{array}$$

撞击的液滴被加入到壁膜中。

• 飞溅

$$\text{Regime boundaries:}{T}_{\text{wall }}<T_{DL,l}\text{and}\{\begin{array}{l}\mathrm{We}>{\mathrm{We}}_{\mathrm{SW}}\text{ and Wet wall }\\ \text{ or }\\ \mathrm{We}>{\mathrm{We}}_{\mathrm{SD}}\text{ and Dry wall }\end{array}$$

冲击液滴的一部分质量沉积形成薄膜，其余部分则被雾化。飞溅质量分数和飞溅液滴变量根据Ansys Fluent理论指南中Kuhnke模型的“飞溅”部分进行计算。

• 蒸发飞溅

区域边界：$T_{\text{wall }}>T_{DL,u}$ 且 $We>5$

假设一小部分冲击的液体喷雾在撞击时立即蒸发，而其余部分则完全雾化，释放出较小尺寸的液滴。此区域内发生的蒸发事件由随机数采样控制。与飞溅区域类似，飞溅液滴变量如Ansys Fluent理论指南中Kuhnke模型的“飞溅”部分所述进行计算。

• 过渡区域

区域边界：$T_{DL,l}<T_{\text{wall }}<T_{DL,u}$ 且 $We>5$

在这些条件下，冲击液滴可能被分配到任意一个区域（扩散、飞溅或蒸发飞溅）。这取决于随机数与从公式12.269（第517页）计算出的沉积比$dr$​的比较，结合液滴韦伯数和湿壁或干壁标准。

$$dr=\frac{T_{\text{wall }}-T_{DL,l}}{T_{DL,u}-T_{DL,l}}\text{\hspace{1em}(12.269)}$$

12.9.4 分离与剥离子模型

关于薄膜分离模型的详细描述，请参见薄膜分离（第840页）。有关薄膜剥离模型的详细信息，请参见薄膜剥离（第842页）。

12.9.5 壁面薄膜颗粒的守恒方程

以下描述了壁面薄膜中单个颗粒的动量、质量和能量守恒方程。基于颗粒的方法最初由O'Rourke [495]（第1086页）提出，后续推导大部分基于该工作。

12.9.5.1 动量

薄膜中颗粒的动量方程为： $$\rho h\frac{d\overrightarrow{u_p}}{dt}={\tau }_g\overrightarrow{t}+\overrightarrow{{\tau }_w}+\overrightarrow{F}+\rho h\left(\overrightarrow{g}-\overrightarrow{{\alpha }_w}\right)\text{\hspace{1em}(12.270)}$$ $\overrightarrow{u_p}=$ 表示薄膜颗粒速度

${\tau }_g=$ 表示气体流对薄膜表面剪应力的大小

$\overrightarrow{t}=$ 表示薄膜表面速度方向的单位向量

$\overrightarrow{{\tau }_w}=$ 表示壁面对薄膜施加的应力

$\overrightarrow{F}=$ 表示保持薄膜在表面所需每单位面积的力

$h=$ 表示颗粒位置处的当前薄膜高度

$\rho h\left(\overrightarrow{g}-\overrightarrow{{\alpha }_w}\right)=$ 表示体积力项

尽管薄膜厚度值很小，但由于在移动边界模拟中观察到极高的加速度率，体积力项可能变得非常重要。

壁面应力 $\overrightarrow{{\tau }_w}$ 的表达式为： $$\overrightarrow{{\tau }_w}=-2\frac{{\mu }_l}{h}\left(\overrightarrow{u_p}-\overrightarrow{u_w}\right)\text{\hspace{1em}(12.271)}$$ 其中，${\mu }_l$ 代表液膜粘度，而 $\overrightarrow{u_w}$ 则是壁面速度。

力 $\overrightarrow{F}$ 虽然为了完整性被包含在第12.270节（第517页）的方程中，但在Fluent求解的方程组中并未对其进行建模。取而代之的是，通过明确施加以下条件来考虑其影响： $${\overrightarrow{u}}_p\cdot {\hat{n}}_{\alpha }=0\text{\hspace{1em}(12.272)}$$ 在每个时间步长，对于所有代表壁膜的粒子，Ansys Fluent 求解一个形式如下的粒子位置方程： $$\frac{d{\overrightarrow{u}}_p}{dt}=\alpha -\beta {\overrightarrow{u}}_p\text{\hspace{1em}(12.273)}$$ 在重新排列方程12.270（第517页）并代入方程12.271（第518页）后，薄膜粒子的加速度可表示为： $$\frac{d\overrightarrow{u_p}}{dt}=\frac{{\tau }_g\overrightarrow{t}}{\rho h}+\frac{2{\mu }_l\overrightarrow{u_w}}{\rho h^2}+\left(\overrightarrow{g}-\overrightarrow{{\alpha }_w}\right)-\frac{2{\mu }_l}{\rho h^2}\overrightarrow{u_p}\text{\hspace{1em}(12.274)}$$

12.9.5.2 薄膜传质

12.9.5.2.1 薄膜蒸发与沸腾

当薄膜颗粒温度高于蒸发温度 $T_{vap}$ 时，应用薄膜蒸发定律。在 Ansys Fluent 中，蒸发速率采用以下模型进行建模：

• 扩散控制模型

• 对流/扩散控制模型

• 气体侧边界层模型

• 薄膜沸腾模型

12.9.5.2.1.1 扩散控制模型

薄膜蒸发速率受表面暴露于气相的梯度扩散控制。薄膜表面与气相之间的蒸汽浓度梯度决定了蒸发速率。 $${\dot{N}}_i=k_{film}\left(C_{i,s}-C_{i,\infty }\right)\text{\hspace{1em}(12.275)}$$ 其中，$N_i$ 表示蒸气的摩尔通量 $\left(\mathrm{kmol}/{\mathrm{m}}^2-\mathrm{s}\right)$，$k_{\text{film }}$ 是传质系数 $\left(\mathrm{m}/\mathrm{s}\right)$，而 $C_{i,s}$ 和 $C_{i,\infty }$ 分别代表薄膜表面和主体气体中的蒸气浓度 $\left(\mathrm{kmol}/{\mathrm{m}}^3\right)$。

薄膜表面的蒸气浓度通过薄膜温度下的饱和蒸气压来评估，而主体气体浓度则从流动场解中获得。蒸发的速率对饱和蒸气压非常敏感，这一点与液滴蒸发类似。

薄膜传质系数通过以下方式获得： $$k_{\text{film }}=Nu\cdot D_i/\text{ len }\text{\hspace{1em}(12.276)}$$ 其中，$D_i$ 表示物质 $i$ 的二元扩散系数，$Nu$ 是努塞尔特相关系数，而 len 则是特征长度（以米为单位），计算方法为薄膜包覆区域 $A_p$ 的平方根。 $$\text{ len }=\sqrt{A_p}\text{\hspace{1em}(12.277)}$$ $A_p$ 被视为墙体表面积 $A_w$ 的质量加权百分比： $$A_p=A_w\cdot m_p/m_{\text{film }}\text{\hspace{1em}(12.278)}$$ 在计算质量传递系数的公式12.276（第519页）中，采用努塞尔特相关性，并将普朗特数替换为施密特数，具体如下：

• 层流流动： $$Nu=0.332{\mathrm{Re}}^{1/2}S{c}^{1/3}\text{\hspace{1em}(12.279)}$$
• 湍流： $$Nu=0.0296{\mathrm{Re}}^{4/5}S{c}^{1/3}\text{\hspace{1em}(12.280)}$$ 雷诺数 $Re$ 是根据颗粒膜的特征长度 $len$ 以及与壁面平行的颗粒相对速度分量计算得出的。

对于多组分蒸发，每个组分的扩散性用于方程 12.276（第 519 页）和方程 12.279（第 519 页）中，以获得液膜混合物中每个组分的蒸发速率。

颗粒包的质量 $m_p$ 按照以下方式减少： $$m_p\left(t+\Delta t\right)=m_p\left(t\right)-N_i{A}_p{M}_{w,i}\Delta t\text{\hspace{1em}(12.281)}$$ 其中，$M_{w,i}$ 表示气相物种 $i$ 的分子量，该物种是从液体中蒸发的蒸汽所加入的。为了考虑单个包裹中的质量损失，薄膜颗粒的直径会减小。这保持了包裹中液滴数量的恒定。

12.9.5.2.1.2 对流/扩散控制模型

薄膜蒸发速率由以下因素决定： $${\dot{N}}_i=k_{\text{film }}{\rho }_{\infty }\ln \left(1+B_m\right)\text{\hspace{1em}(12.282)}$$

$$B_m=\frac{Y_{i,\text{ surf }}-Y_{i,\infty }}{1-Y_{i,\text{ surf }}}$$ 其中

${\dot{N}}_i=$ 物质 $i$ 的蒸气质量通量 $\left(\mathrm{kg}/{\mathrm{m}}^2/\mathrm{s}\right)$

$Y_{i,\text{ surf }}$ 和 $Y_{i,\infty }=$ 物质 $i$ 在液膜表面和主体中的质量分数

$k_{\text{film }}=$ 液膜传质系数 $\left(\mathrm{m}/\mathrm{s}\right)$，根据公式 12.276计算

${\rho }_{\infty }=$ 气体主体的密度 $\left(\mathrm{kg}/{\mathrm{m}}^3\right)$

12.9.5.2.1.3 热解模型

单速率热解模型使用公式 12.87 (第 468 页)，恒速模型使用公式 12.88 (第 468 页)。如果考虑拉格朗日壁面液膜，可以应用允许区分自由流粒子和液膜粒子热解速率的次级速率模型。该模型使用公式 12.87 (第 468 页) 分别计算从自由粒子到气体主体相的质量传递速率（使用粒子前置指数因子 $A_{\text{particle }}$ 和粒子活化能 $E_{\text{particle }}$）以及从液膜到气体主体相的质量传递速率（使用液膜前置指数因子 $A_{\text{film }}$ 和液膜活化能 $E_{\text{film }}$）。

12.9.5.2.1.4 气体侧边界层模型

当主体中物质的分压低于液膜表面的蒸气压时，Fluent 使用公式 12.291 (第 523 页) 和公式 12.292 (第 523 页) 分别计算层流和湍流的蒸发速率。这些公式使用壁函数来确定传质系数。详情请参见“液膜冷凝”(第 522 页)。

12.9.5.2.1.5 带有导热传热模型的液膜沸腾模型

当液膜温度达到沸点时，液膜沸腾速率通过将能量平衡方程 (公式 12.301 (第 526 页)) 的右侧设为 0 并求解液膜沸腾速率 ${\dot{m}}_p$ 来计算。

得到的沸腾速率表达式如下：

恒温墙体 $${\dot{m}}_p=\frac{A_p}{h_{fg}}\left[-k_p\left(1+\frac{h_{\text{film }}}{k_p+h_{\text{film }}}\right)T_p+k_p{T}_w+h_{\text{film }}\left(1-\frac{h_{\text{film }}}{k_p+h_{\text{film }}}\right)T_{\infty }\right]\text{\hspace{1em}(12.283)}$$ 热流墙 $${\dot{m}}_p=\frac{A_p}{h_{fg}}\left[-\left(\frac{k_p{h}_{film}}{k_p+h_{film}}\right)T_p+{\dot{q}}_w+h_{film}\left(1-\frac{h_{film}}{k_p+h_{film}}\right)T_{\infty }\right]\text{\hspace{1em}(12.284)}$$ 包含： $$k_p=k_l/\left(h/2\right)\text{\hspace{1em}(12.285)}$$

$$h_{\text{film }}=\frac{k_{\infty }}{\text{ len }}Nu$$ 这里 $${\dot{m}}_p=\text{film parcel boiling rate}\left(\mathrm{kg}/\mathrm{s}\right)$$ $A_p=$ 薄膜包裹面积 $\left({\mathrm{m}}^2\right)$

$h_{fg}=$ 潜热 $\left(\mathrm{J}/\mathrm{kg}\right)$

$T_p=$ 薄膜包裹温度 (K)

$T_w=$ 壁面温度 (K)

$T_{\infty }=$ 主流相温度 (K)

${\dot{q}}_w=$ 壁面热流 $\left(\mathrm{W}/{\mathrm{m}}^2\right)$

$h=$ 薄膜高度 $\left(\mathrm{m}\right)$

$k_l=$ 液态薄膜的热导率 $\left(\mathrm{W}/\mathrm{m}/\mathrm{K}\right)$

如果启用了气体侧边界层模型，则传热系数 $h_{\text{film }}$ 的计算方法如下。

• 层流流动：$h_{film}=k_{\infty }/l$，其中长度 $l$ 是从薄膜表面到主流的垂直距离。

• 湍流流动：$h_{\text{film }}$ 根据壁面函数确定，如《薄膜能量传递》(p. 524)中所述。

在其他情况下，Fluent 使用努塞尔数表达式： $$h_{\text{film }}=k_{\infty }Nu/\text{len}\text{\hspace{1em}(12.286)}$$ 其中

$k_{\infty }$ = 连续相的热导率 $\left(\mathrm{W}/\mathrm{m}/\mathrm{K}\right)$

$\text{len}$ = 特征长度 $\left(\mathrm{m}\right)$

努塞尔数 $Nu$ 根据公式 12.299 (第 526 页) 计算。如果薄膜通过与冷壁接触而被冷却至沸点以下，蒸发速率计算将恢复到蒸发方程 (公式 12.275 (第 519 页) 或公式 12.282 (第 520 页))。

12.9.5.2.1.6 带有对流换热模型的薄膜沸腾模型

对流换热模型的沸腾速率表达式如下：

对于恒温壁： $$\begin{array}{c}\frac{d{m}_p}{dt}=[-\frac{1}{1-f\left(\mathrm{Pr}\right)}[h_{wall}+h_{\infty }]T_p+h_{wall}{T}_w+f[\mathrm{Pr}\left]\right[\frac{h_{wall}+h_{\infty }}{1-f\left(\mathrm{Pr}\right)}]T_w+h_{\infty }{T}_{\infty }]\frac{A_p}{h_{lat}}]\\ =\left[-\frac{1}{1-f\left(\mathrm{Pr}\right)}\left(h_{wall}+h_{\infty }\right)T_p+\left(\frac{h_{wall}+f\left(\mathrm{Pr}\right)h_{\infty }}{1-f\left(\mathrm{Pr}\right)}\right)T_w+h_{\infty }{T}_{\infty }\right]\frac{A_p}{h_{lat}}\end{array}\text{\hspace{1em}(12.287)}$$ 对于热流墙： $$\frac{d{m}_p}{dt}=\left[-h_{\infty }{T}_p+\frac{Q_{wall}}{A_p}\left(1+f\left(\mathrm{Pr}\right)\frac{h_{\infty }}{h_{wall}}\right)+h_{\infty }{T}_{\infty }\right]\frac{A_p}{h_{lat}}\text{\hspace{1em}(12.288)}$$ 随着 $$f\left(\mathrm{Pr}\right)=\frac{3}{8}\stackrel{-1/3}{\mathrm{Pr}}\text{\hspace{1em}(12.289)}$$ 在上述方程中，

$m_p=$ 薄膜颗粒的质量（千克）

$T_p=$ 薄膜颗粒的温度（开尔文）

$T_{\infty }=$ 整体温度（开尔文）

$T_w=$ 壁面温度（开尔文）

$h_{\text{lat }}=$ 蒸发潜热 $\left(J/kg\right)$

$h_{\text{wall }}=$ 从壁面到液膜的热传递系数 $\left(W/K\cdot {\mathrm{m}}^2\right)$

$h_{\infty }=$ 从液膜到整体相的热传递系数 $\left(W/K\cdot {\mathrm{m}}^2\right)$

$A_p=$ 液膜包裹面积，取为壁面面积的质量加权百分比 $\left({\mathrm{m}}^2\right)$

$Q_{\text{wall }}=$ 施加在壁面上的热通量（瓦特）

$\mathrm{Pr}=$ 液体的普朗特数 $\mathrm{Pr}=\frac{c_{p,l}{\mu }_l}{k_l}$

从液膜到整体相的热传递系数 $\left(h_{\infty }\right)$ 的获取方法如传导热传递模型中的薄膜沸腾模型所述（第520页）。

从壁面到液膜的热传递系数 $h_{\text{wall }}$ 的获取方法如方程12.312所述（第529页）。

如果液膜通过接触冷壁而冷却到沸点以下，蒸发速率的计算将返回到蒸发方程（方程12.275（第519页）或方程12.282（第520页））。

12.9.5.2.2 薄膜凝结

当整体中物种 $i$ 的分压超过其在薄膜表面的蒸气压时，凝结模型被触发。FLUENT拉格朗日壁膜模型在层流和湍流条件下使用两种不同的表达式来描述壁膜凝结。推导遵循[739]（第1100页），假设气相中物种 $i$ 的凝结发生在非冷凝物的存在下。

在Ansys Fluent拉格朗日壁膜冷凝模型中，组分$i$在存在非冷凝组分的情况下，会在壁膜上发生冷凝。对于组分$i$通过距离膜表面$l$处的平面的冷凝速率，菲克定律可以表述为： $$\dot{m}=-D_i{\rho }_{\infty }/\left(1-y_i\right)\frac{\partial y_i}{\partial l}\text{\hspace{1em}(12.290)}$$ $\dot{m}=$ 物种蒸气的质量通量 $\left(\mathrm{kg}/{\mathrm{m}}^2-\mathrm{s}\right)$

${\rho }_{\infty }=$ 主体相的密度 $\left(\mathrm{kg}/{\mathrm{m}}^3\right)$

$D_i=$ 物种 i 的扩散系数 $\left({\mathrm{m}}^2/\mathrm{s}\right)$

$y_i=$ 物种 i 的质量分数

对于层流，物种 $i$ 的质量分数从薄膜表面到主体气体条件进行积分，得到对数项，而在湍流情况下，直接采用有限差分。最终表达式为：

• 层流边界层： $$\dot{m}=k_{\text{film }}{\rho }_{\infty }\ln \left[\frac{1-Y_{i\infty }}{1-Y_{is}}\right]\text{\hspace{1em}(12.291)}$$
• 湍流边界层： $$\dot{m}=-k_{\text{film }}{\rho }_{\infty }\frac{Y_{i\infty }-Y_{is}}{1-Y_{is}}\text{\hspace{1em}(12.292)}$$ 其中，

$k_{\text{film}}$ = 质量传递系数 $\left(\mathrm{m}/\mathrm{s}\right)$

$Y_{is}$ 和 $Y_{i\infty }$ = 物种 $i$ 在薄膜表面和主体中的质量分数，分别

对于层流流动，方程 12.291（第 523 页）中的质量传递系数 $k_{film}$ 计算如下： $$k_{film}=\frac{D_i}{l}\text{\hspace{1em}(12.293)}$$ 其中，$l$ 是从薄膜表面到本体区域的法向距离。$l$ 计算为从相邻壁面单元中心到壁面的法向距离。

对于湍流流动，方程 12.292（第 523 页）中的湍流质量传递系数 $k_{film}$ 由壁面函数确定。当前实现使用方程 4.348（第 141 页）所述的壁面函数。当在湍流模型中启用增强壁面处理时，实现遵循增强壁面处理 $\epsilon$ 方程（EWT-ε）（第 145 页）中的计算。

对于多组分颗粒，方程 12.291（第 523 页）或方程 12.292（第 523 页）分别应用于每个液体组分。

有关如何使用薄膜冷凝模型的更多信息，请参阅《Fluent 用户指南》中的薄膜冷凝模型部分。

12.9.5.2.3 多组分薄膜模型

对于多组分拉格朗日壁面薄膜，蒸发速率计算为各组分蒸发速率的总和。液体薄膜中每个组分的蒸发速率从以下公式计算：

• 扩散控制模型：方程 12.275（第 519 页）和方程 12.276（第 519 页）

• 对流/扩散控制模型：方程 12.282（第 520 页）

如果在薄膜中为任何组分启用了热解模型，Ansys Fluent 应用方程 12.164（第 486 页）以获得单速率下的各组分热解速率，或方程 12.165（第 486 页）以获得恒定速率。

如果考虑拉格朗日壁面薄膜，可以使用次级速率模型。该模型使用方程 12.164（第 486 页）分别计算自由流颗粒的蒸发速率（使用颗粒预指数因子 $A_{\text{particle}}$ 和颗粒活化能 $E_{\text{particle}}$）和薄膜的蒸发速率（使用薄膜预指数因子 $A_{film}$ 和薄膜活化能 $E_{\text{film}}$）。

当多组分膜表面的总蒸气压超过细胞压时，Ansys Fluent采用多组分膜沸腾方程，根据壁面施加的热边界条件，使用方程12.283（第520页）或方程12.284（第521页）的平均形式。平均潜热$h_{fg}$根据方程12.170（第487页）计算，基于各液体组分的潜热和膜的总体沸腾速率，通过对方程12.283（第520页）或方程12.284（第521页）计算的传导传热模型速率进行平均，以及对方程12.287（第522页）或方程12.288（第522页）计算的对流传热模型速率进行平均，并结合方程12.171（第487页）计算的分馏蒸发量，来计算每个液体组分的总体沸腾速率。

12.9.5.3 膜中的能量传递

12.9.5.3.1 传导传热模型

为了获得膜中温度的方程，必须考虑来自气体侧和壁侧的能量通量。在没有膜蒸发/沸腾的情况下，从/到膜粒子的能量传递由以下公式给出： $$\frac{d}{dt}\left(m_p{c}_p{T}_p\right)=Q_{p,\text{ cond }}+Q_{p,\text{ conv }}\text{\hspace{1em}(12.294)}$$ 其中，

$m_p,c_p$, 和 $T_p$ 分别表示薄膜包裹的质量、比热容和温度。

$Q_{p,\text{ cond }}$ 表示从壁面传导的热量。

$Q_{p,\text{ conv }}$ 表示薄膜表面处的对流热通量。

图 12.13：薄膜中双线性温度分布的假设

当壁面上的热边界条件设定为恒定温度时，$Q_{p,\text{ cond }}$ 计算如下： $$Q_{p,\text{ cond }}=\frac{k_l{A}_p}{\left(h/2\right)}\left(T_w-T_p\right)\text{\hspace{1em}(12.295)}$$ 其中

$k_l=$ 液膜的热导率

$A_p=$ 根据方程12.278（第519页）计算的液膜包面积

$h=$ 液膜高度

对于热通量由边界条件给定的壁面，${\dot{q}}_w$，方程12.295（第525页）被替换为： $$Q_{p,\text{ cond }}={\dot{q}}_w{A}_p\text{\hspace{1em}(12.296)}$$ 对流热通量，即 $Q_{p,\text{ conv }}$，其计算公式如下： $$Q_{p,\text{ conv }}=A_p{h}_{\text{film }}\left(T_{\infty }-T_s\right)\text{\hspace{1em}(12.297)}$$ 其中，

$h_{\text{film }}=$ 传热系数 $\left(\mathrm{W}/\mathrm{K}-{\mathrm{m}}^2\right)$

$T_s=$ 薄膜表面温度

如果启用了薄膜冷凝或气体侧边界层模型，传热系数 $h_{\text{film }}$ 的计算方法如下：

• 层流： $$h_{\text{film }}=\frac{k_{\infty }}{l}\text{\hspace{1em}(12.298)}$$ 其中，$l$ 表示从薄膜表面到体相内部的法向距离。$l$ 的计算方法是从相邻壁面单元的中心到壁面的法向距离。

• 湍流流动：

$h_{\text{film}}$ 的确定依赖于壁面函数。当前的实现采用了标准壁面函数（第138页）中描述的壁面函数，具体见方程4.343（第140页）。当湍流模型中启用了增强壁面处理时，实现遵循增强壁面处理 $\epsilon$ 方程（EWT- $\epsilon$）（第145页）中的计算方法。

在其他情况下，Fluent 使用努塞尔数表达式来计算 $h_{\text{film}}$： $$h_{\text{film }}=k_{\infty }Nu/\text{ len }$$ 式中 $$k_{\infty }=\text{thermal conductivity of the continuous phase}\left(\mathrm{W}/\mathrm{m}/\mathrm{K}\right)$$

$$\text{len}=\text{the characteristic length}\left(\mathrm{m}\right)\text{computed from Equation 12.277 (p. 519)}$$ 努塞尔数 $Nu$ 的计算公式如下： $$\text{Laminar Flow:}Nu=0.332{\mathrm{Re}}^{1/2}\mathrm{Pr}1/3\text{\hspace{1em}(12.299)}$$

$$\text{Turbulent Flow:}Nu=0.0296{\mathrm{Re}}^{4/5}\mathrm{Pr}1/3$$ 在公式12.297（第525页）中，薄膜表面温度 ( T_s ) 是根据界面守恒关系确定的： $$h_{\text{film }}\left(T_{\infty }-T_s\right)=\frac{k_l}{\left(h/2\right)}\left(T_s-T_p\right)\text{\hspace{1em}(12.300)}$$ 当粒子在蒸发过程中质量发生变化时，需在方程12.294（第524页）中增加一项，以考虑蒸发热的影响： $$\frac{d}{dt}\left(m_p{c}_p{T}_p\right)=Q_{p,\text{ cond }}+Q_{p,\text{ conv }}-{\dot{m}}_p{h}_{fg}\text{\hspace{1em}(12.301)}$$ 其中，${\dot{m}}_p$ 是液膜包裹体的蒸发速率 $\left(\mathrm{kg}/\mathrm{s}\right)$，而 $h_{fg}$ 是蒸发潜热 (J/kg)。