阻抗谱不是电池的一个附加特征。EIS 把从超快电子传输到宏观固相扩散的整条时间尺度谱系压缩进频域：高频截距界定欧姆电阻，中频圆弧参数化电荷转移动力学与 SEI 双电层电容，低频 Warburg 尾支配有限边界扩散极限。一次测量同时看到了三个物理层次。
但把它变成可用的机器学习标签，难点不在于画出一张形状合理的 Nyquist 图。难点在于阻抗谱天然是欠定的——截然不同的物理参数组合可以产生几乎相同的曲线。这篇讲怎么打破这种简并，以及怎么记录才能让半年后的你还知道每条曲线是什么意思。
PyBaMM 电池建模系列（共 4 篇）：① 架构与求解器 → ② EIS 标签 → ③ 数据集管线 → ④ 训练 SOH/RUL。本文是第 ② 篇。
一、EIS 的复变函数基础
小信号 EIS 响应是在非线性稳态工作点附近的线性扰动。基础对象是复阻抗传递函数：
$$ Z(\omega) = \frac{\tilde{V}(\omega)}{\tilde{I}(\omega)} = Z_{re}(\omega) + j Z_{im}(\omega) $$
等效电路模型中，$n$ 阶 RC 网络耦合 Warburg 元件 $Z_W$：
$$ Z_{ECM}(s) = R_\Omega + \sum_{k=1}^{n} \frac{R_k}{1 + s R_k C_k} + Z_W(s), \quad s = j\omega $$
物理驱动的 AI 需要的特征远超曲线拟合。要从连续体方程（浓溶液理论、Butler-Volmer 动力学）中提取有物理边界的指标：
- 高频欧姆截距 $R_\Omega$：$\omega \to \infty$ 的极限，由电解液电导率和集流体接触电阻主导。
- 电荷转移弧：$-\max(Z_{im})$，对应双电层电容 $C_{dl}$ 与电荷转移电阻 $R_{ct}$，受多孔电极比表面积调制。
- Warburg 扩散尾 $\sigma_W$：$\omega \to 0$ 的渐近相位偏移，描述锂在石墨/NMC 晶格中的嵌入。
注意”线性扰动”这个前提。EIS 的整套数学建立在小信号线性化上——如果扰动幅值过大，系统的非线性会让阻抗的定义本身失效。这一点在后面的合法性检验里会再出现。
二、在 PyBaMM 里生成频域响应
先说一个会卡住很多人的版本问题：EIS 功能曾经在独立包 pybammeis 里，现已并入 PyBaMM 核心。网上的老教程写的是 import pybammeis，新版 PyBaMM 里则是 pybamm.EISSimulation。装错包或抄错导入，报的错和阻抗毫无关系，很容易往错误方向排查。
import numpy as np
import pybamm
# 版本兼容：新版在核心包里，旧版需要单独装 pybammeis
if hasattr(pybamm, "EISSimulation"):
EISSimulation = pybamm.EISSimulation
else:
import pybammeis # pip install pybammeis
EISSimulation = pybammeis.EISSimulation
model = pybamm.lithium_ion.DFN(options={
"surface form": "differential", # 见下方说明，这个不是可选项
"particle shape": "spherical",
"SEI": "solvent-diffusion limited",
})
params = pybamm.ParameterValues("Chen2020")
params["SEI kinetic rate constant [m.s-1]"] = 1e-15
eis = EISSimulation(model, parameter_values=params)
# 对数分布频率向量：1 mHz 到 10 kHz
frequencies = np.logspace(-3, 4, 60)
eis.solve(frequencies)
"surface form": "differential" 必须设置。它让 DAE 求解器正确处理双电层的电容特性；缺了它，频域响应在结构上就是错的——而且不会报错，你会得到一条形状看着还行、物理上没有意义的曲线。
怎么取出原始阻抗数组
官方 README 只演示了 nyquist_plot()，没有记录如何取原始数组，而这个属性名在各版本间变动过。与其抄一个可能对不上的写法，不如让代码自己找：
def extract_impedance(eis):
"""在不同 PyBaMM/pybammeis 版本间稳健地取出复阻抗数组。
硬编码属性名的脚本会在升级后拿到 AttributeError。"""
for attr in ("solution", "impedances", "Z", "impedance"):
z = getattr(eis, attr, None)
if z is None:
continue
z = np.asarray(z).ravel()
if np.iscomplexobj(z):
return z
raise AttributeError(
f"没找到复阻抗数组。该对象的可用属性: "
f"{[a for a in dir(eis) if not a.startswith('_')]}"
)
z = extract_impedance(eis)
nyquist_tensor = np.column_stack((frequencies, z.real, z.imag))
第一次跑的时候，失败分支会把对象上所有可用属性打印出来——比翻文档快，也比猜准。
三、可辨识性：为什么单张 Nyquist 图不够
这是整篇最重要的一节。一条阻抗曲线对应的物理参数组合不是唯一的。
直观理解：中频圆弧的大小由 $R_{ct}$ 决定，而 $R_{ct}$ 同时受交换电流密度、电极比表面积和温度影响。比表面积下降一半、交换电流密度上升一倍，圆弧几乎不动。你在数据里看到的是它们的乘积，不是各自的值。
让神经网络从单张谱去回归”活性物质损失”，在这种简并下它只能学到一个条件均值，并且在训练集分布之外立刻崩掉。打破简并要靠扫描——同一块电芯在多个工作点下的谱，联合起来才能分离这些纠缠的参数：
| 扫描维度 | 分离了什么 | 物理依据 |
|---|---|---|
soc |
热力学状态 vs 动力学退化 | 开路电势随嵌入分数 $\theta$ 变化，而欧姆电阻基本不随 SOC 变 |
temperature_c |
活化能不同的过程 | 各速率常数遵循不同的 Arrhenius 依赖，温度扫描把它们拉开 |
protocol_id |
路径依赖的滞后 | DST、WLTP 等负载曲线留下不同的浓度梯度历史 |
实践建议：把”同一电芯的一组扫描”作为一个训练样本，而不是把每条谱当独立样本。后者不仅浪费了打破简并的信息，还会造成严重的组内泄漏——同一电芯的相邻 SOC 点分落训练和测试，模型靠插值就能作弊。
四、合法性检验：Kramers-Kronig
并不是所有算出来的曲线都是合法的阻抗响应。KK 关系对任何线性、因果、稳定、时不变的系统必然成立——违反了，就说明你的谱违背了这四条前提中的至少一条，通常是线性化在深度退化工作点上失效，或者雅可比矩阵条件数太差。
工程上最实用的做法不是直接算 KK 积分（需要外推到无穷频率，本身就有误差），而是用线性 KK 检验：用一组纯 RC 元件（它们天生满足 KK）去拟合实测谱，看残差。
import numpy as np
def linear_kk_residual(freq, z, n_rc=40):
"""用 n_rc 个时间常数对数分布的 RC 元件拟合阻抗。
RC 网络恒满足 KK，所以拟合不上的部分就是 KK 违规量。
返回相对残差；经验上超过百分之几就该怀疑这条谱。"""
w = 2 * np.pi * np.asarray(freq)
taus = np.logspace(np.log10(1 / w.max()), np.log10(1 / w.min()), n_rc)
# 基函数: 每个 RC 的 1/(1+jw*tau)，外加一个常数项代表 R_ohm
A = np.column_stack([1.0 / (1.0 + 1j * w[:, None] * taus), np.ones((len(w), 1))])
# 复数最小二乘：实部虚部堆成实数方程组一起解，
# 分开拟合会丢掉两者必须共享同一组参数这个约束
A_ri = np.vstack([A.real, A.imag])
z_ri = np.concatenate([z.real, z.imag])
coef, *_ = np.linalg.lstsq(A_ri, z_ri, rcond=None)
z_fit = A @ coef
return float(np.linalg.norm(z - z_fit) / np.linalg.norm(z))
resid = linear_kk_residual(frequencies, z)
if resid > 0.02:
print(f"KK 残差 {resid:.3f} 偏高，这条谱标记为可疑，不进监督训练")
注意实部虚部要堆成一个方程组联合求解。分开拟合会丢掉”两者必须由同一组 RC 参数产生”这个约束——而这个约束恰恰就是 KK 关系的实质。
五、标签 schema：哪些字段不能省
只存 Z_re 和 Z_im 的数据集是不可用的。半年后你拿到一条曲线，无法判断它的变化来自温度、SOC 还是退化机制。
| 字段 | 含义 | 为什么不能省略 |
|---|---|---|
frequency_hz |
频率采样点 | 阻抗张量必须与频率轴严格对齐；换了频率网格的样本不能直接混用 |
z_re, z_im |
复阻抗实部与虚部 | 构成 Nyquist 与 Bode 特征 |
soc, temperature_c |
工作点 | 区分退化效应与热力学状态效应，见第三节 |
model_name, parameter_set |
DFN/SPMe/SPM 与参数来源 | 不同物理模型的谱不可混淆 |
lli, lam_neg, lam_pos, soh |
退化标签 | 监督目标；正负极 LAM 必须分开 |
solver_status |
求解器是否收敛 | 防止失败样本被当作真实物理响应 |
kk_residual |
线性 KK 检验残差 | 量化这条谱的合法性，让下游可以按阈值筛选而不是一刀切 |
sweep_id |
同一电芯一组扫描的分组键 | 做 group split 的依据，也是打破参数简并的单位 |
六、三类必查的异常
每次生成后必须扫一遍，这三类都不会抛异常：
- 实部出现负值：无源系统的阻抗实部必须非负。出现负值一定是数值伪像，直接丢弃而不是取绝对值。
- 低频尾部尖峰：低频段求解最慢也最容易不收敛。检查
solver_status，同时对比相邻频点的差分——真实的 Warburg 尾是平滑的，尖峰是求解失败的痕迹。 - 重复仿真不一致：同一 SOC/温度/参数跑两次，结果应当逐位一致（这是确定性 PDE 求解）。不一致说明有未固定的随机源或求解器状态污染，必须先查清楚再继续生成。
七、拆分与泄漏
最容易犯的错误是把同一条退化轨迹的相邻工作点同时放进训练和测试。模型会记住这条仿真轨迹的局部形状，而不是学到可泛化的电化学。按 sweep_id 或 cell_design_id 分组拆分：
from sklearn.model_selection import GroupShuffleSplit
splitter = GroupShuffleSplit(test_size=0.2, random_state=42)
train_idx, test_idx = next(splitter.split(X, y, groups=sweep_id))
# 断言写进流水线
assert not (set(sweep_id[train_idx]) & set(sweep_id[test_idx])), "分组泄漏"
八、生成管线
cd pybamm-ai-data-lab
python3 -m venv .venv
source .venv/bin/activate
pip install -r requirements.txt
python src/run_all.py --backend pybamm --workers 8 --precision float64--precision float64 不是可以省的性能优化项。EIS 求解涉及病态矩阵求逆，float32 在低频段会累积出可见的误差，表现为 Warburg 尾的抖动——那看起来像物理噪声，其实是精度不够。
最后一条原则：降级 surrogate 的插值结果绝不能混进主训练张量。DAE 求解器的容差失败标记的是刚性物理状态（例如零下温度接近 0% SOC），必须显式标注或 mask，不能静默插补——那等于把”我们不知道”伪装成”我们测到了”。
References
- PyBaMM EIS Simulation DAE Solvers
- pybamm-eis（已归档，功能并入 PyBaMM 核心）
- Physics-based battery model parametrisation from impedance data
- PyBOP: battery model optimisation and parameterisation
本文第三节讲的”单张阻抗谱欠定”，在时域拟合里是同一个问题：PyBaMM 参数拟合与可辨识性 讲怎么判断优化器返回的数字是测出来的还是凑出来的。