浸入边界法（IBM）
理论与物理
IBM的基本概念
老师，浸没边界法（Immersed Boundary Method）和ALE法有什么区别？
ALE法中网格与结构界面是贴合（body-conforming）的，而IBM中是在固定的正交网格（Cartesian grid）上“嵌入”结构。结构的边界与网格无关地定义，通过体积力（forcing term）或插值来施加边界条件。
这是Peskin（1972）为心脏瓣膜血流模拟而构思的方法，原始公式如下。
其中 $\mathbf{F}$ 是拉格朗日界面上的力密度，$\delta$ 是狄拉克δ函数，$\mathbf{X}(t)$ 是界面位置。
狄拉克δ函数是如何离散化的？
使用Peskin的正则化δ函数。四点离散δ函数是标准方法，
其中 $h$ 是网格间距，$\phi$ 是支撑宽度为4的光滑核函数。
IBM的分类
IBM也有好几种类型吗？
大致可分为两大体系。
| 分类 | 方法 | 原理 | 精度 |
|---|---|---|---|
| 连续力法 | Peskin型 | 在欧拉方程中添加源项 | 一阶（δ函数的扩散） |
| 离散力法 | Fadlun型、幽灵单元法 | 直接修正离散方程 | 可达二阶精度 |
| 切割单元法 | 笛卡尔切割单元 | 用界面切割单元 | 二阶精度 |
连续力法（Peskin型）实现简单，但界面会因δ函数的支撑宽度而模糊。壁面 $y^+$ 控制困难，因此不适用于高雷诺数的湍流壁面流动。
离散力法（幽灵单元法或直接力法）可以锐利地表现界面，壁面边界层的分辨率也能达到接近ALE法的精度。但实现复杂，尤其是移动界面时的切割单元稳定处理困难。
切割单元法和IBM有什么区别？
严格来说，切割单元法是IBM的一种，在界面穿过单元的位置对单元进行几何切割，并在切割单元上满足守恒律。守恒性高、精度好，但需要处理切割单元体积变得非常小的问题（小单元问题）。通过单元合并或通量再分配等方法应对。
催生IBM的是“无法分析心脏瓣膜”的烦恼
浸没边界法（IBM）的原型诞生于20世纪70年代，源于Charles Peskin在纽约大学挑战心脏二尖瓣模拟。瓣膜形状复杂运动，传统的边界贴合网格需要每步重新划分网格，不切实际。于是产生了“在固定的正交网格上分散瓣膜影响力”的想法。为医疗应用而发明的IBM，如今也用于水下机器人和风力涡轮机设计，这向我们展示了基础研究影响的深远。
各项的物理意义
- 时间项 $\partial(\rho\phi)/\partial t$：想象一下拧开水龙头的瞬间。最初水流不稳定地喷溅，过一会儿才变成稳定的水流，对吧？描述这个“正在变化的过程”的就是时间项。心脏搏动导致血流脉动，发动机阀门每次开闭引起流动变化，这些都是非定常现象。那么定常分析是什么？就是“经过足够长时间流动稳定后”的状态——也就是令此项为零。计算成本大幅降低，因此先用定常求解是CFD的基本策略。
- 对流项 $\nabla \cdot (\rho \mathbf{u} \phi)$：把落叶扔进河里会怎样？会被水流带着向下游移动，对吧？这就是“对流”——流体运动携带物质的效果。暖气的热风能到达房间另一端，也是因为空气这个“搬运工”通过对流输送热量。这里有趣的是——此项包含“速度×速度”，因此是非线性的。也就是说，流速变快时此项会急剧增强，变得难以控制。这就是湍流的根本原因。常见误解：“对流和传导差不多”→ 完全不一样！对流是流动携带，传导是分子传递。效率有天壤之别。
- 扩散项 $\nabla \cdot (\Gamma \nabla \phi)$：有过在咖啡里倒入牛奶后放置的经历吗？即使不搅拌，过一会儿也会自然混合。那就是分子扩散。那么下一个问题——蜂蜜和水，哪个更容易流动？当然是水。因为蜂蜜的粘度（$\mu$）高，所以不易流动。粘度越大，扩散项越强，流体运动就越“粘稠”。雷诺数小的流动（缓慢、粘稠）中扩散占主导。相反，雷诺数大的流动中对流占压倒性优势，扩散则成为配角。
- 压力项 $-\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$。直接关系到时间步长的稳定性 |
数值解法与实现
直接力法
请告诉我实用的IBM实现方法。
直接力法（Mohd-Yusof, 1997; Fadlun et al., 2000）应用最广泛。在界面附近的欧拉单元中，将速度强制为所需的边界条件值。
其中 $\mathbf{u}^*$ 是无强制力的中间速度，$\mathbf{u}_{BC}$ 是界面上应满足的速度（无滑移条件则为结构的速度）。
相当简洁呢。精度如何？
界面位于单元中间时为一阶精度。改进版使用基于界面距离的插值达到二阶精度。幽灵单元法在界面内侧设置虚拟单元（ghost cell），通过反射插值施加边界条件。
IBM-FSI耦合
IBM如何与结构进行耦合？
IBM-FSI耦合流程如下。
1. 在欧拉网格上求解流体（包含IBM强制力）
2. 将界面上的流体力插值并传递到拉格朗日结构
3. 推进结构时间步（更新位移·速度）
4. 根据新的界面位置更新IBM掩码
5. 返回步骤1
与ALE法相比，IBM-FSI的优缺点如下。
| 比较项目 | ALE-FSI | IBM-FSI |
|---|---|---|
| 网格变形 | 需要 | 不需要 |
| 大变形 | 困难（需要重新划分网格） | 容易 |
| 接触·碰撞 | 非常困难 | 可以应对 |
| 壁面边界层精度 | 高 | 稍低 |
| 实现的容易度 | 中等 | 中等～高 |
| 守恒性 | 高 | 取决于方法 |
IBM对大变形适应性强是最大的优点呢。
没错。心脏瓣膜开闭、降落伞展开、旗帜飘扬等结构大幅变形甚至拓扑结构改变的问题，正是IBM的用武之地。
商用软件中的IBM
商用CFD软件中实现了IBM吗？
有限但存在。
| 软件 | IBM功能 | 备注 |
|---|---|---|
| STAR-CCM+ | 重叠网格（与IBM类似） | 重叠网格实质上起到IBM的作用 |
| Ansys Fluent | 无（用重叠网格替代） | 可通过UDF实现直接力法 |
| OpenFOAM | immersedBoundary（ESI版） | 幽灵单元型IBM |
| Palabos | 标准搭载 | 基于格子玻尔兹曼法 |
严格需要IBM时，使用研究代码（Nek5000、CaNS、PeleLM、AFiD等）或在OpenFOAM中自定义实现是比较现实的选择。
用IBM分析鱼群游动——由此发现的“节能秘诀”
IBM擅长处理多个物体复杂运动的问题，被广泛用于鱼群游动模拟。斯坦福大学等研究通过数值模拟证实，后方的鱼巧妙利用前方鱼产生的涡流，可以削减高达约50%的推进能耗。若非IBM，很难计算多条鱼交替游动，这个发现或许也不会诞生。自然界的节能策略反哺工程学——生物流体力学是个低调而有趣的领域。
迎风格式（Upwind）
一阶迎风：数值扩散大但稳定。二阶迎风：精度提高但有振荡风险。高雷诺数流动中必不可少。
中心差分（Central Differencing）
二阶精度，但佩克莱特数 Pe > 2 时会发生数值振荡。适用于低雷诺数的扩散主导流动。
TVD格式（MUSCL、QUICK等）
通过限制器函数抑制数值振荡同时保持高精度。对捕捉激波或陡峭梯度有效。
有限体积法 vs 有限元法
FVM：自然满足守恒律。CFD的主流。FEM：对复杂形状·多物理场有利。SPH等无网格法也在发展中。
CFL条件（库朗数）
显式法：CFL ≤ 1 为稳定条件。隐式法：CFL > 1 也稳定，但影响精度和迭代次数。LES：推荐 CFL ≈ 1。物理意义：一个时间步内信息传播不超过一个单元。
残差监控
连续性方程·动量·能量的各项残差下降3~4个数量级可判断为收敛。质量守恒残差尤其重要。
松弛因子
なった
詳しく
報告