DeePKS 案例教程：Rh(111) 表面 CO 吸附能计算
作者：梁馨元，邮箱：2201111875@stu.pku.edu.cn
审阅：陈默涵，邮箱：mohanchen@pku.edu.cn
最后更新时间：2026/01/06
前言：密度泛函理论中的“鱼和熊掌”
在计算催化和表面科学中，CO 分子在过渡金属（如 Rh, Pt, Cu）表面的吸附行为是一个经典的测试案例。长期以来，采用密度泛函理论计算的研究者们发现传统的广义梯度近似泛函（如 PBE）在预测 CO 分子在这些体系表面吸附位点时存在明显的位点偏好误差（Site Preference Error），即著名的“CO Puzzle”：PBE 往往预测 CO 倾向于吸附在多配位的中心位（如 FCC 或 HCP 位点），而实验观测结果却显示 CO 优先吸附在单配位的顶位（Top 位点）。
虽然高阶的杂化泛函（如 HSE06）能够纠正这一错误，但其巨大的计算开销使得对复杂表面体系的结构优化和分子动力学模拟变得难以实现。
本教程来自 2025 年由 Xinyuan Liang et al 发表在 AI for Science 1, 025007 期刊上的一篇论文”Investigating CO adsorption on Cu(111) and Rh(111) surfaces using machine learning exchange-correlation functionals"。该研究提出了一种基于 DeePKS 框架的解决方案：利用少量 HSE 高精度电子结构数据训练一个神经网络修正项，并将其加在低阶泛函（PBE）上。这样，我们就能以 PBE 的计算速度，获得媲美 HSE 的预测精度。
为了帮助快速上手，本教程抽取了论文中的核心流程，以 Rh(111)表面吸附 CO 分子为例，演示如何结合 ABACUS 进行表面吸附能计算，包括 PBE 计算、HSE 高精度修正以及 DeePKS 模型的训练与测试。测试所需文件可从这个 Github 链接下载，本教程所有的计算都由 ABACUS 软件的 3.8.0 版本完成。
一、采用 PBE 泛函进行性质测试
在进行复杂的模型训练之前，我们必须首先通过 PBE 泛函建立起物理体系的基准。这一阶段的所有计算都将进行原子的结构优化（Relax），为后续的高精度修正提供稳定的构型。
相关文件存放于：01_adsorption_energy/01_PBE/
1. 吸附能的定义与计算
吸附能是衡量分子在表面结合牢固程度的关键指标。在本教程中，我们定义吸附能如下：
公式含义说明：
- ：吸附能。负值表示吸附过程放热，数值越负说明结合越牢固。
- ：CO 分子吸附在 Rh(111) 表面后的体系总能量。
- ：单独 Rh(111) 表面能量。
- ：真空中孤立 CO 分子的能量。
- ：该吸附体系中包含 CO 分子的数量
当吸附能为负值时，表明 CO 分子吸附于表面会放出能量。吸附能的数值越负，表明放出的能量越多，吸附越稳定。
2. CO 分子
首先，我们需要计算 CO 分子在真空中的能量。
- 路径：
00_CO/ - 操作：执行
bash run.sh $np $abacus。其中 $np 表示计算使用的进程数，$abacus 表示计算使用的 ABACUS 软件可执行文件路径，下同。 - 物理细节：我们将 CO 分子放置在一个较大的真空晶胞中，以防止分子与其周期性镜像产生相互作用。通过 Relax（结构弛豫）计算，我们获得其自由状态下的平衡键长和总能量。
3. 体相 Rh 的晶格常数
由于不同软件和泛函对平衡晶格常数的预测存在微小差异，直接使用实验值可能导致体系内部产生人为应力，故需要使用 ABACUS 软件进行结构优化计算，获得平衡状态下的关键物理参数。
- 路径：
01_cell_lenth_relax/ - 操作：执行
bash workflow_cell_lenth_relax.sh $abacus 关键步骤：
- 初猜构型：利用
STRU_initial开始，其初始参数通常取自实验数据或常用数据库。 - 布里渊区的 k点收敛性测试：这是保证计算精度的第一步。在
01_KPT下，我们需要逐渐加密 k 点网格，直至体系单原子的能量波动收敛在 1 meV 以内。 - 晶胞优化（Cell Relax）：利用收敛的 k 点执行晶胞优化，从而获得最符合当前计算设置的 Rh 晶格常数。这对于后续切割出的表面层间距精度至关重要。
- 初猜构型：利用
4. 表面体系建模
- 路径：
02_slab/ - 自动化构建：运行
bash workflow_slab.sh $abacus 关键步骤：
- 表面构建：
00_make_stru文件夹中，脚本会调用make_stru.py（基于 ASE 库）切割出 Rh(111) 表面。我们构建了一个包含 5 层原子的 Slab 模型，并进行 2×4 重复单元的扩胞。然后调用relax_line.sh，固定底部的 3 层原子以模拟体相环境，仅允许最上方的 2 层原子在计算中自由移动。 - 表面优化：
01_relax文件夹中，对上述构建出来的表面进行了结构优化，获取 Rh(111)表面能量 - 吸附构型构建：
02_make_slab_and_CO文件夹中，通过make_slab_and_CO.py在表面原子的特定几何位置放置 CO 分子。- Top 位：CO 位于表面原子的正上方。
- FCC 位：CO 位于表面三个原子形成的空隙上方，对应第三层原子的位置。 脚本会自动根据层间距估计位置，并为后续的吸附能计算准备好初始结构。
- 表面构建：
5. 吸附构型建模
- 路径：
03_top/和04_fcc/ - 操作： 将上一步中构建的吸附构型初始结构复制到对应文件夹中，并执行
bash run.sh $np $abacus进行结构优化计算。
二、 杂化泛函验证
PBE 泛函往往会预测出错误的 CO 基态吸附位点，如认为 FCC 比 Top 位点稳定，但实验证明 Top 位点应该是更稳定的吸附位点。我们使用 HSE06 杂化泛函来提供高精度的参考。
路径：01_adsorption_energy/02_HSE/
命令：bash workflow_run.sh $np $abacus
- 策略：考虑到 HSE06 计算昂贵，直接做 Relax 的成本极高，我们采取 “PBE 优化构型 + HSE SCF 单点能计算” 的方式进行吸附能计算。
workflow_run.sh会自动提取 PBE 优化后的构型进行 HSE 计算 - 收敛控制：我们在
INPUT中开启了exx_real_number以加速计算，并同时设置了scf_thr和scf_ene_thr来保证电子步的准确性 分析结论：在
01_PBE与02_HSE下均可以执行get_energy.sh提取能量数据，并进行吸附能计算。PBE 与 HSE 的参考能量及吸附能结果总结于下列表格中，其中能量单位均为 eV，且所有计算方法采用的构型均为 PBE 进行结构优化后的结果结果表明，PBE 计算的吸附能差（Top-FCC）为正值，表明它认为 FCC 构型比 Top 构型更稳定。而 HSE 计算的吸附能差为正值，表明它认为 Top 构型比 FCC 构型更稳定。虽然 HSE 的结果与实验结果一致，但是它的计算开销明显偏大，在 56 个 CPU 核（Intel(R) Xeon(R) Gold 6348 CPU at 2.60GHz）上计算吸附结构的单点 SCF 计算，PBE 大致需要 250 s，而 HSE 需要 $1×10^{4}$s，所需时间是 PBE 的 40 倍。这证明了引入 DeePKS 进行精度修正并保持高计算效率的必要性。
三、 DeePKS 模型训练
这一步是本教程的核心——通过机器学习模型，学习从 PBE 到 HSE 的差别。
- 路径：
02_DeePKS_train/ 数据构成：在
systems/目录下，我们准备了 110 个典型构型（含孤立分子、纯表面及多种吸附状态），分布在group.00到group.010中。丰富的样本涵盖了分子与表面在不同距离、不同位点下的相互作用环境，可以根据需要进行使用。各个文件夹中的构型及数量如下：group.00：CO 分子，共 10 帧group.01-02：Rh(111)表面，共 20 帧group.03-010：CO 吸附于 Rh(111)表面的吸附构型，共 80 帧。其中group.010是测试集。
训练配置：
- 进入
iter/目录，所有的训练权重和迭代逻辑均在此定义。 - 环境配置：请在
scf_abacus.yaml中更新 ABACUS 路径，并在machines.yaml中配置计算资源（核心数、节点信息等）。详情可参考 DeePKS-kit Documentation。
- 进入
- 启动：执行
bash run.sh。随着训练的进行，可以在训练日志 iter.*/01.train/log.train 中看到能量和受力的 Loss（损失）值稳步下降，这意味着模型正在逐渐掌握 HSE 的物理规律。比如 DeePKS-Rh 模型的iter.00/01.train/log.train文件内容截取部分如下，trn_loss_energy/force为训练数据集中的能量/受力损失结果，tst_loss_energy/force为测试数据集中的能量/受力损失结果，由该文件可看出，这些损失结果在逐步下降直到达到收敛。
# using seed: 296828622
# load 10 systems with fields {'gvx', 'eig', 'lb_e', 'lb_f'}
# load 1 systems with fields {'gvx', 'eig', 'lb_e', 'lb_f'}
# working on device: cpu
#epoch trn_loss_energy trn_loss_force tst_loss_energy tst_loss_force
0 1.0527e-04 1.6697e-06 1.2355e-04 1.6879e-06
100 7.4401e-07 1.1174e-06 3.4644e-07 1.0646e-06
200 8.6170e-07 1.2734e-06 3.4947e-07 1.0515e-06
300 4.2330e-07 9.8216e-07 3.9029e-07 1.0438e-06
400 7.6001e-07 1.2293e-06 1.4376e-07 1.0125e-06
500 8.2814e-07 1.0442e-06 4.1805e-07 1.0112e-06
600 5.9581e-07 8.9611e-07 7.1089e-07 1.0025e-06
700 7.0637e-07 1.0686e-06 1.7279e-07 9.7932e-07
800 7.2727e-07 1.0796e-06 5.7555e-07 9.9391e-07
900 4.6372e-07 9.1234e-07 1.7266e-07 9.7904e-07
1000 7.4336e-07 9.8193e-07 4.5155e-07 9.6691e-07
1100 8.4172e-07 1.0442e-06 4.8143e-07 9.5933e-07
1200 6.4534e-07 9.9026e-07 1.6678e-07 9.6040e-07
1300 7.0138e-07 7.9872e-07 2.4255e-07 9.5105e-07
1400 6.2933e-07 9.0397e-07 2.0593e-07 9.5360e-07
1500 4.3967e-07 1.0044e-06 1.5530e-07 9.4350e-07
1600 4.3169e-07 9.3393e-07 1.7085e-07 9.4739e-07
1700 8.0506e-07 9.2013e-07 2.2469e-07 9.3172e-07
1800 5.0249e-07 8.4920e-07 2.2132e-07 9.3199e-07
1900 5.9592e-07 8.4586e-07 1.7829e-07 9.2661e-07
2000 3.9274e-07 7.7909e-07 3.6763e-07 9.2786e-07
2100 4.3962e-07 8.2965e-07 1.5695e-07 9.1888e-07
2200 6.2852e-07 9.8707e-07 2.6523e-07 9.1798e-07
2300 5.2529e-07 9.4738e-07 2.4293e-07 9.1359e-07
2400 6.9039e-07 8.7784e-07 3.5536e-07 9.0674e-07
2500 4.2287e-07 9.1985e-07 1.6493e-07 9.1632e-07
2600 5.9064e-07 8.8513e-07 2.5823e-07 9.0907e-07
2700 4.2411e-07 8.5994e-07 2.7212e-07 8.9692e-07
2800 6.6338e-07 9.0960e-07 3.3543e-07 8.8923e-07
2900 6.2811e-07 9.4254e-07 3.0133e-07 8.8903e-07
3000 5.4608e-07 9.2532e-07 2.9315e-07 8.9335e-07
3100 6.4306e-07 9.4999e-07 1.9330e-07 8.7804e-07
3200 5.0587e-07 9.2644e-07 4.5180e-07 8.8116e-07
3300 6.0241e-07 8.0011e-07 1.5113e-07 8.7715e-07
3400 4.3890e-07 8.7079e-07 1.9734e-07 8.7758e-07
3500 3.8168e-07 7.8445e-07 1.4563e-07 8.7392e-07
3600 4.0038e-07 8.3490e-07 3.3003e-07 8.6537e-07
3700 5.5578e-07 7.0195e-07 3.3352e-07 8.6164e-07
3800 3.5600e-07 7.4717e-07 2.5941e-07 8.5808e-07
3900 5.0819e-07 7.5260e-07 1.4789e-07 8.5646e-07
4000 3.6074e-07 6.7551e-07 1.5708e-07 8.5468e-07
4100 3.3839e-07 8.7797e-07 1.6385e-07 8.5325e-07
4200 3.6129e-07 8.0437e-07 1.6003e-07 8.5508e-07
4300 2.8686e-07 6.3710e-07 1.8270e-07 8.4811e-07
4400 5.6940e-07 7.2762e-07 1.9755e-07 8.4645e-07
4500 4.2648e-07 8.4686e-07 1.5987e-07 8.4921e-07
4600 3.2731e-07 7.7073e-07 1.9087e-07 8.4239e-07
4700 3.8769e-07 7.9189e-07 1.7775e-07 8.3565e-07
4800 4.8347e-07 7.7812e-07 3.7766e-07 8.3013e-07
4900 4.3515e-07 7.8276e-07 1.8542e-07 8.3715e-07
5000 3.8764e-07 8.2905e-07 1.7881e-07 8.2561e-07
四、 DeePKS 效果测试
最后，我们将训练好的模型应用于实际计算，检验其是否真正达到了 HSE 的精度。
1. 模型部署
我们将训练产出的 model.pth 转化为 ABACUS 可直接读取的 model.ptg。
- 路径：
01_adsorption_energy/03_DeePKS/ - 转化：运行
python transform.py。为了方便使用，该文件夹下已内置了一个训练好的模型供参考，该模型即为参考文章中的 DeePKS-Rh 模型，以model.ptg格式提供。
2. 性能验证
运行 bash workflow_run.sh $np $abacus。此时仍然是基于 PBE 优化后的结构直接进行 SCF 单点能计算，不过由于 DeePKS 模型的计算效率较高，如果有需要的话，也可以快速地进行结构优化计算。
- 最终表现：再次计算吸附能以及吸附能差值后，可以观察到如下结果，参考结果也总结在下列表格当中。
- 位点偏好修正：DeePKS 成功预测出 Top 位比 FCC 位更稳定，与实验和 HSE 趋势完全一致。
- 精度对标：吸附能的绝对数值非常接近 HSE 的高精度结果，小于 1 meV/atom。
- 效率优势：完成整个计算所需的时间仅比普通 PBE 略多，远远快于 HSE 泛函。上述说明过，进行吸附结构的单点 SCF 计算，PBE/HSE 大致需要 250/$1×10^{4}$ s，而DeePKS-Rh模型只需要$1.2×10^{3}$s，比 HSE 快了将近 8 倍。