磁気流体力学（MHD）
理论与物理
磁流体力学基础
老师，磁流体力学（MHD）是流体力学和电磁学结合起来的学科吗？
正是如此。MHD是研究导电性流体（液态金属、等离子体、电解质溶液等）中电磁场与流体运动相互作用的领域。当磁场作用于有电流流过的流体时，会产生洛伦兹力，从而改变流动。反之，导电性流体的运动也会感应出磁场。
常见的应用实例：
- 连续铸造: 钢水的电磁制动、电磁搅拌
- 铝电解精炼: Hall-Herault电解槽内熔融铝的流动控制
- 核聚变反应堆: 等离子体的磁场约束
- MHD泵: 无运动部件输送液态金属
- 空间推进: MPD推进器
MHD控制方程
MHD是纳维-斯托克斯方程与麦克斯韦方程组耦合而成的体系。
修正的纳维-斯托克斯方程（添加洛伦兹力）：
其中 $\mathbf{J} \times \mathbf{B}$ 是洛伦兹力（体积力）。
欧姆定律（在运动的导电性流体中）：
$\sigma$ 是电导率，$\mathbf{E}$ 是电场，$\mathbf{u} \times \mathbf{B}$ 是流体运动感应的电动势。
磁场感应方程：
右边第一项是磁场的对流（冻结），第二项是磁场的扩散。
流体的运动方程加入了电磁力，同时磁场的方程也包含了流速。这是完全的双向耦合呢。
MHD无量纲数
我们来整理一下MHD特有的无量纲数。
| 无量纲数 | 定义 | 物理意义 |
|---|---|---|
| 哈特曼数 Ha | $BL\sqrt{\sigma/\mu}$ | 电磁力/粘性力 |
| 磁雷诺数 Rm | $\mu_0 \sigma U L$ | 磁场对流/磁场扩散 |
| 斯图尔特数（相互作用参数）N | $\sigma B^2 L / (\rho U) = \text{Ha}^2/\text{Re}$ | 电磁力/惯性力 |
| 磁普朗特数 Pm | $\mu_0 \sigma \nu$ | 动量扩散/磁场扩散 |
工业用液态金属（钢水、熔融铝等）的 $\text{Pm} \sim 10^{-6}$ 极小。这意味着磁场的扩散远快于流体的动量扩散，因此 $\text{Rm} \ll 1$。在这种情况下，磁场几乎由外部施加的场决定，只需考虑其对流体的影响即可（低Rm近似）。
哈特曼数越大，磁场的影响就越强，对吧。
是的。Ha = 0 时为普通流体力学，Ha → ∞ 时流动完全受磁场约束。连续铸造中 Ha ∼ 100-1000 是典型范围。
核聚变反应堆的等离子体约束——MHD守护着人类的梦想
在ITER（国际热核聚变实验反应堆）中，使用托卡马克型磁场约束1亿摄氏度的等离子体。等离子体稳定性分析的核心正是MHD理论。等离子体是否会从磁场中"逃逸"？是否会发生"扭曲不稳定性"或"气球不稳定性"？通过将这些作为MHD方程的特征值问题来求解，可以优化线圈设计和电流分布。核聚变之所以一直被称为"30年后实现"，正说明了MHD稳定性分析的难度。
各项的物理意义
- 时间项 $\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以上需考虑可压缩性效应
- 布西内斯克近似（自然对流）：仅在浮力项中考虑密度变化，其他项使用恒定密度
- 不适用的情形：稀薄气体（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$。直接关系到时间步长的稳定性 |
数值解法与实现
MHD数值解法
MHD的耦合方程在CFD中是怎么求解的呢？
低Rm近似（工业液态金属的大多数情况）和高Rm问题（等离子体、天体物理）所用的方法不同。
低Rm近似解法
当 $\text{Rm} \ll 1$ 时，感应磁场与施加磁场相比可以忽略，磁场可以作为已知的外部场 $\mathbf{B}_0$ 来处理。电场的控制方程为：
$\phi$ 是电势。求解这个泊松方程得到 $\mathbf{E} = -\nabla\phi$，然后计算 $\mathbf{J} = \sigma(-\nabla\phi + \mathbf{u} \times \mathbf{B}_0)$ 和 $\mathbf{F}_L = \mathbf{J} \times \mathbf{B}_0$，将其作为源项添加到N-S方程中。
计算步骤：
1. 根据流场计算 $\mathbf{u} \times \mathbf{B}_0$
2. 求解电势的泊松方程
3. 计算电流密度 $\mathbf{J}$
4. 将洛伦兹力 $\mathbf{F}_L$ 加入N-S方程更新流场
5. 重复步骤1-4直至收敛
只是增加了一个电势方程，基本上还是在N-S的延长线上呢。
高Rm问题解法
对于 $\text{Rm} \gg 1$（等离子体物理、天体物理）的情况，需要完整求解磁场感应方程。
主要的数值方法：
- 约束输运 (CT) 法：通过单元面上的磁通守恒严格保持 $\nabla \cdot \mathbf{B} = 0$
- 散度清除法：添加修正方程来阻尼 $\nabla \cdot \mathbf{B}$ 的误差
- 矢量势法：令 $\mathbf{B} = \nabla \times \mathbf{A}$，自动保证 $\nabla \cdot \mathbf{B} = 0$
保持 $\nabla \cdot \mathbf{B} = 0$ 是数值MHD的根本性课题。如果这一点被破坏，就会产生非物理的磁单极力，导致计算崩溃。
哈特曼流动 — MHD基本验证问题
MHD最基本的解析解是哈特曼流动。在平行平板间的流动上施加垂直于壁面的磁场 $B_0$。
Ha = 0 时为抛物线分布（泊肃叶流动），Ha → ∞ 时，除了壁面附近的哈特曼层（厚度 $\delta_H \sim H/\text{Ha}$）外，分布变得平坦。
| Ha | 速度分布 | 哈特曼层厚度 |
|---|---|---|
| 0 | 抛物线 | 无 |
| 10 | 略微平坦化 | $H/10$ |
| 100 | 中心区域几乎平坦 | $H/100$ |
| 1000 | 仅在哈特曼层有梯度 | $H/1000$ |
Ha = 1000 时，哈特曼层只有通道宽度的1/1000啊。对网格的要求似乎很严格。
正是如此。要解析哈特曼层，需要在壁面附近布置与Ha数成比例的网格层数。这是MHD计算的主要成本之一。
MHD数值解法的难关——保持"磁场散度为零"的艰辛
MHD分析的数值实现中最令人头疼的就是维持 ∇·B=0（磁通连续性）。用有限体积法离散化磁场时，数值误差会累积导致 ∇·B≠0，从而产生非物理的磁力，使解发生畸变。为了防止这种情况，开发了"散度清除（divergence cleaning）"和"约束输运法（constrained transport）"等方法。特别是在太阳风的大规模MHD模拟中，如果散度误差超过1%，就会产生数值上的"磁单极子"，导致等离子体向错误的方向飞散。磁场散度为零是物理的要求，也是数值解法的战场。