大多数把 PyBaMM 当工具用的教程，都是从”跑一条放电曲线”开始的。这篇反过来：先看清 PyBaMM 内部是怎么把电化学 PDE 编译成可求解的方程组的。理由很实际——当求解器在半夜第 300 圈崩掉时，你需要知道崩在哪一层，而这只能靠理解它的分层架构。
顺便澄清一个搜索时的常见拼写问题：正式名称是 PyBaMM（Python Battery Mathematical Modelling），不是 PyBAMM。EIS 部分曾经在独立包 pybammeis 里，现已并入核心，直接用 pybamm.EISSimulation 即可，它在频域内对底层 DAE 系统做线性化。
PyBaMM 电池建模系列（共 4 篇）：① 架构与求解器 → ② EIS 标签 → ③ 数据集管线 → ④ 训练 SOH/RUL。本文是第 ① 篇。
一、把 PyBaMM 看成一条 DAE 编译管线
PyBaMM 架构上最关键的设计，是在任何数值离散化之前先把电化学抽象成符号表达式树。Simulation 类扮演的是编译器前端：把连续介质方程翻译成离散状态空间。
实例化顺序不是随意的，它就是编译的各个阶段：
- 模型：SPM、SPMe，或基于 P2D 多孔电极理论的 DFN。
- 子模型选项：SEI 生长、析锂过电势、颗粒破裂力学、LAM、表面积公式。这一步必须在构造模型时完成——事后往参数里塞值不会启用任何机制。
- 参数化：注入 OCP 曲线等非线性函数和标量属性，把符号映射为物理量。
- 实验协议：电流密度输入、截止电压、CCCV 循环、边界条件。
- 离散化与求解：有限体积法做空间离散，PDE 转 DAE，交给能处理极端刚性的 BDF 求解器。
理解这个顺序有直接的调试价值：报错信息出现在哪一阶段，决定了你该去查什么。参数缺失在第 3 阶段报，网格问题在第 5 阶段报，而”模型跑通了但没有退化”根本不报错——那是第 2 阶段的沉默失败。
二、DFN 的形式化与刚性从哪来
DFN 求解固相与电解液相耦合的守恒定律。颗粒内固相锂浓度 $c_s(r,x,t)$ 在球坐标下受菲克第二定律支配：
$$ \frac{\partial c_s}{\partial t} = \frac{1}{r^2} \frac{\partial}{\partial r} \left( r^2 D_s(c_s) \frac{\partial c_s}{\partial r} \right) $$
通过 Butler-Volmer 动力学与颗粒-电解液界面耦合，决定体积转移电流密度 $j(x,t)$：
$$ j = a_s i_0 \left[ \exp\left(\frac{\alpha_a F \eta}{R T}\right) – \exp\left(-\frac{\alpha_c F \eta}{R T}\right) \right] $$
交换电流密度 $i_0$ 依赖 $c_s$、电解液浓度 $c_e$ 与局部过电势 $\eta = \phi_s – \phi_e – U_{OCP}(c_s)$。
刚性的来源可以说得更具体：这个系统的特征时间常数横跨约十个数量级——双电层电容响应在微秒级，固相扩散在小时级。刚性比就是这两端的比值。显式方法（比如 RK4）受 CFL 条件约束，步长必须小到能分辨最快的那个过程，于是要模拟一小时就得走上亿步。这不是”慢”，是根本不可行。
三、符号树与自动微分
PyBaMM 不硬编码稀疏矩阵，而是用符号节点（pybamm.Variable、pybamm.grad、pybamm.div）构造方程。下面这段是可以直接运行的，用来看清表达式树和符号求导：
import pybamm
x = pybamm.Variable("x")
expr = 3 * x**2 + 2 * x
print(expr) # 打印表达式树本身
print(expr.diff(x)) # 符号求导，得到新的表达式树而非数值
真实模型里的扩散算子形式如下。注意 grad/div 必须作用在带 domain 的变量上——空间算子需要知道自己在哪个几何域上离散，缺了 domain 会直接报错：
c_s = pybamm.Variable("Solid concentration", domain="negative particle")
D_s = pybamm.Parameter("Solid diffusion coefficient")
N_s = -D_s * pybamm.grad(c_s) # 扩散通量
dcdt = -pybamm.div(N_s) # 浓度变化率
符号图的价值在于：它让 PyBaMM 可以无缝替换参数集、在任意网格上离散，并通过 CasADi 做自动微分，解析地构造精确雅可比矩阵。精确雅可比直接减少隐式求解器每步的牛顿迭代次数——这不是锦上添花，对刚性系统而言它决定了能不能收敛。
四、求解器怎么选
空间离散后 DFN 产生一个庞大的刚性 ODE+代数方程系统。PyBaMM 把符号树编译成 CasADi 的 SX/MX 计算图，生成优化过的 C 代码来求值右端项和雅可比，再交给 SUNDIALS（IDA/IDAS 或 CVODES）的变阶变步长 BDF 方法。求解器自适应步长：电流阶跃瞬间用极小步长，静置阶段可以拉到数分钟。
| 求解器 | 适用 | 注意 |
|---|---|---|
CasadiSolver |
默认，通用性最好 | 长循环实验上明显慢于 IDAKLU |
IDAKLUSolver |
长时间退化模拟的首选 | 稀疏线性代数，规模越大优势越明显 |
ScipySolver |
纯 ODE 模型的快速试验 | 处理不了带代数约束的 DAE |
容差是另一个必须显式管理的旋钮。默认容差在健康电池上没问题，但在深度退化状态下会引入与电化学现象无法区分的数值噪声——你会以为看到了某种退化特征，实际是求解器在抖。做退化研究时，把容差和求解器版本一起写进样本元数据。
五、三个模型不是精度递进，是不同的状态空间
对做机器学习的人来说，这一节最关键：SPM、SPMe、DFN 的区别不只是”准不准”，它们表达的物理量根本不同。用错模型生成的数据，标签本身就没有对应的物理含义。
- SPM：假设电解液电导率无穷大、浓度均匀。适合宏观 SOH 或简单 RUL 预测，前提是电解液动力学不构成限制。高倍率下它系统性偏乐观，不是随机误差。
- SPMe：重新引入电解液浓度梯度的解析近似。覆盖了最主要的一阶效应，成本远低于 DFN，是大范围参数扫描的实用选择。
- DFN：解析负极、隔膜、正极的空间分布。如果目标是高倍率极化、局部析锂，或频域中固液耦合主导的 EIS 弧段，DFN 不可替代。
实用组合：SPMe 扫参数空间，DFN 在筛出的子集上复算。前提是把用了哪个模型记进元数据，否则两批数据混在一起就再也分不开了。
六、AI 样本的 Schema 设计
数据集要为可证伪性而设计。一条有审计价值的样本应包含：
- 边界与初始条件：初始 $c_s$ 分布（SOC）、温度场、倍率限制、频率激励网格。
- 观测特征：终端电压统计量、$Z_{re}$ 与 $Z_{im}$、高频截距、Warburg 尾系数。
- 内嵌物理标签：SEI 消耗的锂摩尔数（LLI）、颗粒破裂导致的活性物质体积分数损失（LAM，正负极分开）、局部 ECM 参数。
- 求解器指纹：求解器类型、容差、版本、失败标志、步长统计。
七、样本审计表
用这张表判断一条样本是否具备研究价值，而不只是”能跑出一条曲线”。
| 审计维度 | 必须记录 | 研究价值 | 不可接受的信号 |
|---|---|---|---|
| 模型结构 | SPM/SPMe/DFN、启用的退化子模型、参数集版本 | 决定样本能表达哪些物理机制 | 只存电压曲线，不知来自哪个模型或参数集 |
| 实验协议 | C-rate、温度、截止电压、静置时间、EIS 频率网格 | 决定观测量与内部状态的可辨识性 | 训练与测试共享同一条轨迹的相邻 cycle |
| 求解器状态 | 求解器类型与版本、容差、失败标志、步长统计 | 区分真实物理变化与数值误差 | 不收敛的高度退化样本仍被写入训练集 |
| 轨迹完整性 | 实际完成圈数 vs 请求圈数 | 识别被提前终止的截断轨迹 | 不同长度的轨迹混在一起当等价样本 |
| 标签提取 | SOH、LLI、LAM、ECM 参数的提取公式与单位 | 保证学到的是可解释的物理量 | 标签只是 health_score，无单位、无公式、无来源 |
八、三个方法论陷阱
- 把仿真数据当成 ground truth。DFN 生成的数据只是特定数学理论的投影。必须用非线性最小二乘与真实充放电设备数据校准。
- 追求样本数量而忽视物理多样性。朴素蒙特卡洛采样出 10,000 条几乎相同的曲线，会造成数据流形坍缩。应在正交的物理维度上做拉丁超立方采样。
- 忽视容差设置。从深度退化状态提取阻抗特征时，容差不当会引入无法与电化学现象区分的数值噪声。
还有一个更根本的：合成数据最大的风险不是不准，而是太自洽。真实电池有制造离散性、温度不均、传感器噪声，PDE 解一样都没有。在仿真上漂亮的指标，到真实电池上可能直接崩掉。
References
- Python Battery Mathematical Modelling (PyBaMM), Journal of Open Research Software
- PyBaMM EIS Simulation documentation
- PyBaMM coupled degradation notebook
- pybamm-eis（已归档，功能并入核心）
本文开头那个”求解器在半夜第 300 圈崩掉”的场景，单独写了一篇排查指南：PyBaMM 求解器收敛失败：先分清数值问题还是物理问题——包含 IDAKLU 全部选项的默认值，以及一个设了却静默不生效的陷阱。