実在気体効果
理论与物理
概述
老师，“真实气体效应”就是指理想气体假设不成立的情况对吧？在什么情况下会成为问题呢？
大致分为两种情况。一种是高温下气体分子发生振动激发、解离、电离的情况（如极超音速再入、等离子体等）。另一种是高压·低温下分子间力和分子体积不可忽略的情况（如超临界CO₂循环、LNG工艺等）。两者本质上都是对 $pv = RT$ 理想气体状态方程的偏离。
是用压缩因子 $Z$ 来测量这个偏离对吧？
没错。压缩因子定义为
，理想气体时 $Z = 1$。高压天然气中 $Z \approx 0.8$，超临界CO₂中 $Z$ 有时会低至0.2-0.5。
状态方程（EOS）
那么用什么状态方程来代替理想气体呢？
我们来看看几个代表性的EOS（状态方程）。
范德华方程:
这里 $a$ 表示分子间引力，$b$ 表示排除体积。这在历史上很重要，但精度有限。实际应用中广泛使用的是Peng-Robinson（PR）或Soave-Redlich-Kwong（SRK）方程。
Peng-Robinson EOS:
这里 $\kappa = 0.37464 + 1.54226\omega - 0.26992\omega^2$，$\omega$ 是偏心因子。
只要有临界温度 $T_c$、临界压力 $p_c$ 和偏心因子 $\omega$ 就可以用于任何物质了是吧。
是的。对于CO₂，$T_c = 304.1$ K, $p_c = 7.38$ MPa, $\omega = 0.225$。不过PR-EOS倾向于低估液相密度5-15%。有时会加入Peneloux体积修正项来修正这一点。
高温空气的热化学模型
高温侧的真实气体效应是如何建模的呢？
对于空气，随着温度升高，依次发生以下现象。
| 温度范围 | 现象 | 影响 |
|---|---|---|
| < 800 K | 理想气体行为 | $\gamma \approx 1.4$ |
| 800-2500 K | O₂振动激发 | $\gamma$ 下降 |
| 2500-4000 K | O₂解离 | O原子生成 |
| 4000-9000 K | N₂解离 | N原子生成 |
| > 9000 K | 电离 | e⁻, N⁺, O⁺ 生成 |
假设化学平衡时，通过吉布斯自由能最小化来求各化学物种的组成。有限速率反应模型则需要指定各反应的阿伦尼乌斯参数。
比热比 $\gamma$ 随温度变化，这影响似乎很大啊。
正是如此。激波后温度超过3000 K时，$\gamma$ 会从1.4下降到1.1-1.2。这会导致激波角和密度比与理想气体的预测值有很大差异。
范德华在1873年意识到的事情——分子有“大小”
为真实气体理论奠定基础的是荷兰物理学家Johannes van der Waals。在1873年的博士论文中，他提出了“气体分子具有有限大小，且分子间存在引力”这一在当时具有革命性的想法。由此推导出的范德华状态方程，是第一个能够描述高压·低温区域气体行为的实用模型。现代工程中使用的Peng-Robinson或Redlich-Kwong方程，也都是对范德华思想的精炼。150年前的博士论文，直接影响了今天可压缩CFD中状态方程的选择。
各项的物理含义
- 时间项 $\partial(\rho\phi)/\partial t$：想象一下打开水龙头的瞬间。一开始水流会不稳定地喷溅，过一会儿才变成稳定的水流，对吧？描述这个“变化过程中”的就是时间项。心脏搏动导致血流脉动，发动机阀门每次开闭引起流动变化，这些都是非定常现象。那么定常分析是什么？就是只看“经过足够时间流动稳定之后”——也就是令此项为零。计算成本会大幅下降，因此先用定常分析求解是CFD的基本策略。
- 对流项 $\nabla \cdot (\rho \mathbf{u} \phi)$：把落叶丢进河里会怎样？会被水流带着往下游漂，对吧？这就是“对流”——流体的运动携带物质的效果。暖风的暖气能到达房间另一端，也是因为空气这个“搬运工”通过对流输送热量。这里有趣的是——这项包含“速度×速度”，因此是非线性的。也就是说，流速变快这项会急剧增强，变得难以控制。这就是湍流的根本原因。常见的误解：“对流和传导差不多”→ 完全不一样！对流是流动携带，传导是分子传递。效率有天壤之别。
- 扩散项 $\nabla \cdot (\Gamma \nabla \phi)$：有过在咖啡里倒入牛奶后放置不管的经历吗？即使不搅拌，过一会儿也会自然混合。那就是分子扩散。那么下一个问题——蜂蜜和水，哪个更容易流动？当然是水，对吧？因为蜂蜜的粘度（$\mu$）高，所以难流动。粘度越大扩散项越强，流体的运动就变得“粘稠”。雷诺数小的流动（缓慢、粘稠）中扩散占主导。相反，Re数大的流动中对流占压倒性优势，扩散则成为配角。
- 压力项 $-\nabla p$：注射器的活塞一推，液体就从针头有力地射出，对吧？为什么？因为活塞侧压力高，针头侧压力低——这个压力差就是推动流体的力。水坝放水也是同样原理。天气图上等压线密集的地方呢？没错，会刮强风。“有压力差的地方就会产生流动”——这就是纳维-斯托克斯方程压力项的物理含义。这里的误解点：CFD中的“压力”多指表压而非绝对压力。切换到可压缩分析时结果突然变得奇怪，原因可能就是混淆了绝对压力/表压。
- 源项 $S_\phi$：受热的空气会上升——为什么？因为比周围空气轻（密度低），被浮力推上去了。这个浮力作为源项添加到方程中。其他例子还有，燃气灶的火焰产生化学反应热，工厂的电磁泵对金属熔液施加洛伦兹力…这些都是“从外部向流体注入能量或力”的作用，都用源项表示。忘记源项会怎样？自然对流分析中忘记加入浮力，流体就完全不动——就像冬天房间里开了暖气，暖空气却不上浮，这种物理上不可能的结果。
假设条件与适用范围
- 连续介质假设：克努森数 Kn < 0.01（分子平均自由程 ≪ 特征长度）时成立
- 牛顿流体假设：剪切应力与应变速率呈线性关系（非牛顿流体需要粘度模型）
- 不可压缩假设（Ma < 0.3时）：将密度视为常数。马赫数0.3以上需考虑可压缩性效应
- Boussinesq近似（自然对流）：仅在浮力项中考虑密度变化，其他项使用恒定密度
- 不适用的情形：稀薄气体（Kn > 0.1）、超音速·极超音速流动（需要激波捕捉）、自由表面流动（需要VOF/Level Set等方法）
量纲分析与单位制
| 变量 | SI单位 | 注意点·换算备忘 |
|---|---|---|
| 速度 $u$ | m/s | 入口条件从体积流量换算时，注意截面积单位 |
| 压力 $p$ | Pa | 区分表压与绝对压力。可压缩分析使用绝对压力 |
| 密度 $\rho$ | kg/m³ | 空气: 约1.225 kg/m³@20°C、水: 约998 kg/m³@20°C |
| 粘性系数 $\mu$ | Pa·s | 注意与运动粘度系数 $\nu = \mu/\rho$ [m²/s] 混淆 |
| 雷诺数 $Re$ | 无量纲 | $Re = \rho u L / \mu$。层流/湍流转换的判定指标 |
| CFL数 | 无量纲 | $CFL = u \Delta t / \Delta x$。直接关系到时间步长的稳定性 |
数值解法与实现
EOS的数值实现
将真实气体EOS集成到CFD中时，计算上有什么问题吗？
有几个重要的点。理想气体可以从守恒变量（$\rho, \rho\mathbf{u}, \rho E$）解析地转换到温度或压力，但PR-EOS等则需要迭代计算。
具体来说，当给定内能 $e$ 和密度 $\rho$ 时，需要通过牛顿法从下式求解温度 $T$：
这里 $e_{departure}$ 是从EOS导出的偏离函数（departure function），对于PR-EOS为：
这个迭代计算需要针对单元数×时间步长进行，因此计算成本会增加。
相比理想气体，计算会慢多少呢？
一般是2-5倍左右。由于EOS的求值次数占主导，实际工作中通常用查表法（look-up table）来加速。预先计算好温度和压力的2D表格，运行时只需插值即可。
超临界流体的数值困难
超临界状态在数值上有什么难点呢？
临界点附近热力学性质会发生急剧变化。定压比热 $c_p$ 在拟临界温度（pseudo-critical temperature）附近出现峰值，密度也急剧变化。这种陡峭的变化是导致数值振荡或发散的原因。
例如CO₂的超临界条件（p = 8 MPa, T = 305 K附近），$c_p$ 会跃升到通常的10倍以上。密度仅变化几K就会从700 kg/m³急剧变化到200 kg/m³。
变化这么剧烈啊。数值格式需要特别的处理呢。
黎曼求解器的扩展
真实气体也能用Roe格式或HLLC吗？
可以用，但需要修改。计算Roe平均状态时必须去掉理想气体假设。真实气体的Roe格式中，平均音速使用广义的 $\Gamma = v/(c_p)(\partial p / \partial T)_v$ 按如下方式计算（如Vinokur的公式化）：
HLLC对EOS的依赖较少，实现更容易，因此在真实气体计算中更常被选用。
考虑到实现的麻烦程度，HLLC更实用呢。
是的。OpenFOAM的rhoCentralFoam基于KNP格式，不依赖于EOS，因此只需替换EOS就能进行真实气体计算，这是它的优点。
真实气体查找表的“网格点不足”问题
实时计算真实气体的状态方程（范德华法、Peng-Robinson法等）会使CFD速度慢几十倍。因此在实际应用中，使用在温度·压力网格点上预先计算热力学量的查找表（LUT）。问题在于“表格网格点太少则插值误差大，太多则内存和读取时间增加”的权衡。将密集网格集中在冲击波通过的高温高压区域，在平缓区域使用粗网格的“自适应表格”是高效的，但其设计本身已成为一项技术诀窍。
迎风格式（Upwind）
一阶迎风: 数值扩散大但稳定。二阶迎风: 精度提高但有振荡风险。高雷诺数流动中必不可少。
中心差分（Central Differencing）
二阶精度，但Pe数 > 2时会产生数值振荡。适用于低雷诺数的扩散主导流动。
TVD格式（MUSCL、QUICK等）
通过限制器函数抑制数值振荡同时保持高精度。对捕捉激波或陡峭梯度有效。
有限体积法 vs 有限元法
FVM: 自然满足守恒律。CFD的主流。FEM: 对复杂形状·多物理场有利。SPH等无网格法也在发展中。
CFL条件（库朗数）
显式方法: CFL ≤ 1为稳定条件。隐式方法: CFL > 1也稳定，但影响精度和迭代次数。LES: 推荐 CFL ≈ 1。物理意义: 一个时间步内信息传播不超过一个网格。
残差监控
连续性方程·动量·能量的各项残差下降3~4个数量级可判断为收敛。质量守恒的残差尤其重要。
松弛因子
压力: 0.2~0.3、速度: 0.5~0.7是常见的初始值。发散时降低松弛因子。收敛后可提高以加速。
非定常计算的内部迭代
在每个时间步内迭代直到收敛到定常解。内部迭代次数: 5~20次为参考值。残差在时间步间波动时需重新审视时间步长。
SIMPLE法的比喻
SIMPLE法是“交替调整”的方法。先假设求解速度（预测步），然后根据该速度修正压力以满足质量守恒（修正步），再用修正后的压力修正速度——重复这个抛接球过程以逼近正确答案。类似于两人调整架子使其水平的作业。
なった
詳しく
報告