ORCA使用D4：开关对照与结果提取

给ORCA输入增加 D4即可启用相应色散修正。容易出错的地方是把输出中的色散项又加到最终总能量上。本篇用四个固定模型的PBE与PBE-D4实际单点结果说明输入、提取和检查方法。

下载输入、真实参考结果与绘图脚本。解压后进入资料包所在目录，先画已有结果：

bash view.sh

图片和逐条件数值分别写入 analysis/figures、analysis/bookkeeping.csv。脚本使用集群的 ASE3.29.0、NumPy、Matplotlib 环境，画参考数据不提交计算。结果目录已存在时，指定一个新目录：

bash view.sh --output another-analysis

本组在鸢算科技集群实际计算了中性水、相同几何的带电水，以及氧原子间距2.9 Å、4.2 Å的两个水分子对。水分子对按同一朝向侧向摆放，属于固定结构的能量记账示例，没有做氢键构型优化。带电水总电荷为+1、多重度为2；其余为中性单重态。不同主程序使用不同基组与电子芯层处理，只比较各自开关D4的差，不能直接相减两种软件的绝对DFT能量。

1 建立目录并提交#

本例使用ORCA6.1.1官方预编译程序及OpenMPI4.1.8运行环境，模块为 orca/6.1.1，PBE/def2-SVP、VeryTightSCF、DEFGRID3。官方二进制所用编译器不根据MPI编译器推测。主程序输出的内置D4版本为3.4.0；独立对照程序为4.2.0，本篇相同结构及参数的数值逐项比较，不能据此认为全部版本和参数都无差别。

python3 prepare.py orca
cd cases/orca-water-plain
sbatch job.slurm
cd ../orca-water-d4
sbatch job.slurm

每个模板申请1张对应GPU、1个CPU核心，核时用户使用模板里的核时分区；套餐用户按账户规则选择分区。ORCA在这里走CPU程序路径，申请GPU资源不表示ORCA本例有GPU加速。CP2K采用对应GPU构建。不要在登录节点直接运行完整DFT计算。

准备与运行分开：prepare.py按所选软件建立输入目录，随后只提交所需的目录。若已有 cases，请重新解压到新目录，避免覆盖已有结果。只有选择VASP时需要您已有合法授权的O、H PAW势，它们不在下载包中；选择CP2K或ORCA不会读取VASP势文件。

2 查看完整ORCA输入#

中性水、D4开启：

! PBE def2-SVP EnGrad VeryTightSCF DEFGRID3 D4
%pal nprocs 1 end
%maxcore 1000
* xyzfile 0 1 structure.xyz

EnGrad同时计算单点能量与解析梯度；nprocs 1与本例申请的一个CPU核一致。将第一行的 D4删除便是plain对照。带电水输入为：

! PBE def2-SVP EnGrad VeryTightSCF DEFGRID3 D4 UKS
%pal nprocs 1 end
%maxcore 1000
* xyzfile 1 2 structure.xyz

xyzfile后面依次是总电荷、多重度、坐标文件。改变电荷后应同步修改独立D4的 --charge，不能沿用中性体系设置。

ORCA主程序由完整路径直接启动一次，不能用 mpirun多次启动驱动程序。下载包的Slurm脚本加载模块后调用的核心命令为：

"$ORCA_ROOT/orca" orca.inp > orca.out

3 提取能量并检查收敛#

grep 'SCF CONVERGED AFTER' orca.out
grep 'ORCA TERMINATED NORMALLY' orca.out
grep 'Dispersion correction' orca.out
grep 'FINAL SINGLE POINT ENERGY' orca.out

中性水的最终能量为 −76.272750392776 Hartree，色散项约为 −0.000196920 Hartree。最终能量已经包含色散项。与同一设置的plain结果相减，差值为 −0.000196920342 Hartree。色散打印行的小数位少于最终能量，因此可以看到末位舍入差。

总电荷+1的水在相同几何下，D4修正约为 −0.000098867458 Hartree。它说明EEQ电荷模型会影响修正项，不是水分子的电离能。

4 从engrad读出梯度，再换算为力#

在任一已结束的ORCA案例目录执行：

python3 - <<'PY'
import numpy as np
lines=open('orca.engrad').read().splitlines()
v=[float(x) for x in lines if not x.startswith('#') and len(x.split())==1]
n=int(v[0])
energy=v[1]
gradient=np.array(v[2:2+3*n]).reshape(n,3)
forces=-gradient*27.211386245988/0.529177210903
print('Energy / Hartree:',energy)
np.savetxt('forces-eV-A.txt',forces,header='Fx Fy Fz / eV Angstrom^-1')
PY

orca.engrad保存的是总能量梯度，力要加负号。其原始梯度单位为Hartree/Bohr，不能直接当作eV/Å。比较D4修正时，取D4开启和关闭的两组总力之差，而不是把总力直接与独立D4的力相等。

5 批量比较自己的结果#

回到根目录：

python3 d4.py
bash view.sh --cases cases --d4 d4-results --output my-analysis --partial

已完成的模型需要成对plain/d4结果。analysis/bookkeeping.csv是下载包整组参考表；自己的表在指定输出目录。本例四种模型的D4梯度传递与独立工具吻合。

def2-SVP属于小基组，水分子对还可能有基组叠加误差。上述结果演示启用、提取和避免重复修正，没有证明结合能、优化几何或电离能的收敛。

参考：ORCA色散修正手册、ORCA官方色散教程。