VASP + Phonopy 算完声子以后，很多人会看到一串 ZPE、Fvib(T)、Svib(T)、Cv(T)，然后不知道怎样放进反应自由能图。关键是先把问题拆开：VASP 总能给 ΔE，Phonopy 或频率计算给的是振动热校正，电化学模型再给电位、pH、气相和溶剂项。Phonopy 输出不是反应台阶本身，它只提供某个结构的振动贡献；反应台阶需要把反应前后各物种的校正按同一参照相减。
对吸附与催化来说，真正要交付的是每一步 ΔG，而不是单个结构的热力学表。比如 OER 的 *OH、*O、*OOH，HER 的 *H，CO2RR 的 *COOH 和 *CO，都要把吸附态、干净表面、气相或 CHE 参照放在同一套能量口径下处理。只把某个吸附态的 ZPE 加到总能上，不会自动得到反应自由能。
VASP 总能只给电子能量差
结构优化或静态单点得到的 EDFT 是电子基态能量，通常对应 0 K 下不含振动熵的能量。吸附能和反应能可以先用 EDFT 相减得到 ΔE，但实验讨论的电催化台阶通常需要Gibbs 自由能 ΔG。两者之间的差别，主要来自零点能、熵、温度、溶剂、电位和 pH 等项。
ΔG = ΔEDFT + ΔZPE – TΔS + ΔGsolv + ΔGU + ΔGpH
这个式子不是要求每篇工作都把所有项算到极致，而是提醒每一项的来源要清楚。ΔEDFT 来自同一泛函和同一参照的总能差；ΔZPE 和 ΔS 来自频率或声子热力学；ΔGU 和 ΔGpH 来自电化学口径；溶剂项则来自显式水、隐式溶剂、经验校正或文献口径，不能和真空总能混成一个不说明来源的数字。
图1 亚稳 IrO2 结构 OER 活性计算、标度关系和火山图示意。DOI：10.1021/jacsau.2c00709
反应台阶要按反应式相减
若某一步写成 A → B，那么热校正也要写成 GB – GA。对于吸附态，A 和 B 可能是同一个 slab 上的不同中间体；对于气相参照，可能涉及 H2、H2O、O2、CO2、CO 或 N2。不能把每个结构的热校正单独列完就结束，还要按反应式把它们组合成 ΔZPE 和 ΔS。热校正的单位对象是物种，反应自由能的单位对象是台阶。
ZPE 来自所有正频率振动模式
零点能可以理解为量子振动在 0 K 下仍保留的能量。对一个结构，频率越高，ZPE 越大。吸附态中的 O-H、C-H、N-H 伸缩振动常常对 ZPE 贡献明显，所以 *OH、*OOH、*H 的热校正不能随手忽略。
ZPE = 1/2 Σi hνi
如果 Phonopy 给的是整个超胞的声子热力学量，就要确认它对应的是哪一个结构。做吸附自由能时，常见处理是对吸附物及其附近可动原子做频率或局域振动校正，而不是把巨大 slab 的全部低频模式粗暴塞进每个台阶。对于低频表面模式，数值噪声和约束方式会明显影响熵项，需要说明处理口径。
图2 掺杂 MoS2 的声子态密度示意，可用于理解振动模式如何进入热校正。DOI：10.1021/acsomega.3c02623
振动自由能和熵随温度变化
Phonopy 的热力学输出里，Fvib(T) 与 Svib(T) 是温度相关量。常用写法可以把振动自由能写成零点能加上温度项，熵则来自自由能对温度的导数。实际做反应台阶时，通常在 298.15 K 或目标反应温度下取数。
Fvib(T) = ZPE + kBT Σi ln[1 – exp(-hνi/kBT)]
Svib(T) = -∂Fvib(T)/∂T
需要特别区分固体吸附态和气相分子。气相 H2、H2O、CO2、CO 等通常具有平动、转动和振动熵，不能只用一个小盒子里的分子频率随便代替。很多电催化工作采用 CHE 口径处理 H+ + e–，并对气相分子使用热化学表或统一校正值。
以 OER 台阶为例
OER 常见四步可写成 * → *OH → *O → *OOH → *。每一步都要把对应吸附态的 EDFT、ZPE、S 和参照分子组合起来。若采用 CHE，H+ + e– 常用 1/2 H2 代替，电位项再按转移电子数修正。
ΔGi(U,pH) = ΔGi(0,0) – n e U + m kBT ln10 × pH
这里 n 是该步转移电子数，m 是与质子数相关的 pH 系数。若每一步都是一对质子电子转移，电位项斜率相同；若某一步是化学吸附、偶联或脱附，电位项就不一定直接降低它的自由能。把台阶拆清楚后，热校正和电位校正才不会串位。
图3 Mo 氧化物催化过程自由能图、表面构型和 PDOS 对照示意。DOI：10.1038/s41467-024-53047-1
以 HER 台阶为例
HER 中最常见的是 * + 1/2H2 → *H，目标量是 ΔGH*。这一步的 ΔE 来自 E*H – E* – 1/2EH2；热校正则来自 *H 与参照 H2 的 ZPE 和熵差。若在碱性 HER 中讨论水解离，还要把 H2O、*OH 或界面水结构列入对应反应式。
ΔGH* = E*H – E* – 1/2EH2 + ΔZPE – TΔS
图4 分子催化剂 HER 自由能路径示意，展示不同步骤怎样组成反应自由能图。DOI：10.1021/acsomega.5c10226
这个式子看起来简单，但最容易出错的是参照不一致。例如 *H 的频率来自一个 slab，干净 * 的能量来自另一套 k 点或不同磁性态，H2 的能量又来自不同 ENCUT，这样得到的 ΔG 就不是一个干净台阶。热校正只是在总能差基础上补项，不能补救参照混乱。
先确定每个反应物种，再分别取热校正
比较清楚的顺序是：先完成每个中间体结构优化；再对最终结构做同口径静态单点；随后为需要热校正的物种生成位移结构，用 VASP 计算力常数，再由 Phonopy 得到频率、ZPE 和热力学函数。最后建立一张反应台阶工作表，把每个物种的 EDFT、ZPE、TS、溶剂项、电位项按反应式相减。
phonopy –vasp -d –dim=”2 2 1″
# 对生成的位移结构分别计算力
phonopy –fc vasprun.xml
phonopy -t –tmin=0 –tmax=800 –tstep=10
这段片段对应有限位移法的基本流程。dim 的选择要根据体系尺寸和相互作用范围决定；slab 体系通常不希望真空方向重复过多。若只关心吸附物局域振动，也可以用 VASP 频率计算处理可动原子组，但要说明约束范围。
图5 HER 催化循环与弛豫结构示意，说明自由能台阶必须绑定具体中间体构型。DOI：10.1021/acsomega.5c10226
虚频和低频要单独处理
若目标是最低能吸附态，明显虚频通常说明结构没有真正停在局部极小点，应该重新优化或检查位移、力收敛和约束。对于很小的低频模式，尤其是表面吸附物的受限平动和转动，数值会对熵项很敏感。此时可以采用统一截断、类谐振子处理或只报告 ZPE 主导校正，但口径必须一致。
写自由能图时，不建议把每个 Phonopy 输出截图堆上去。更有用的是说明：哪些物种做了振动校正，气相分子使用何种热化学口径，电位和 pH 项怎样加入，最终每个台阶 ΔG 是怎样由这些项相加得到的。这样才能看清曲线每一个拐点从哪里来。
还有一个容易忽略的边界：若反应前后 slab 原子数、吸附覆盖度或磁性态发生改变，热校正不能只按名称相减。比如 *OH 覆盖下的表面和干净表面若已经发生明显重构，ΔG 不只包含 *OH 的振动贡献，还包含工作表面本身改变带来的能量差。此时应先定义反应起点，再决定哪个结构的热力学量进入台阶。
实际整理时，可以为每个中间体保留一行：结构名称、EDFT、ZPE、TS、Gcorr、溶剂项、电位项、pH 项和最终 G。随后每个台阶只做相邻两行的差。这样能避免把校正项手动加错位置，也能清楚检查某个异常台阶究竟来自总能差、ZPE、熵还是电位修正。
要点一：Phonopy 给的是单个结构的振动热力学量，反应台阶要按反应式相减。
要点二：ΔG 通常由 ΔEDFT、ΔZPE、-TΔS、溶剂项、电位项和 pH 项共同组成。
要点三：吸附态、干净表面和气相参照必须使用同一能量口径，否则热校正无法修正参照错误。
要点四：气相分子不能简单当作 slab 里的局域振动处理，平动、转动和热化学参照要单独说明。
要点五：自由能图的核心是每一步 ΔG，而不是单独展示某个结构的声子热力学表。
