VASP与critic2电子密度拓扑
用NaCl实算寻找临界点和键径，理解密度、拉普拉斯与拓扑连接的含义。
阅读本篇目录
Bader电荷回答每个空间盆中分到了多少电子；QTAIM临界点分析进一步观察电子密度在空间中的极值、鞍点及连接。本例沿用NaCl的VASP电荷密度，用critic2提取临界点和键径。
1. 准备一致的电荷密度#
采用2原子NaCl原胞，常规立方晶格参数5.70 Å，PBE PAW64、ENCUT=550 eV、EDIFF=1E-6，FFT电荷网格120×120×120。静态计算设置LAECHG=.TRUE.，得到同一结构、同一网格的AECCAR0和AECCAR2。
VASP使用6.6.1 P100版，NVHPC 25.3、CUDA 12.8、MKL 2023.2及配套MPI构建。后处理使用critic2/1.4模块，实际程序1.4.155，由GCC 13.3构建。
2. 选参考密度并寻找临界点#
在包含CHGCAR、AECCAR0和AECCAR2的目录建立topology.cri：
crystal CHGCAR
load AECCAR0
load AECCAR2
load as "$1+$2"
reference 3
auto
这里CHGCAR用于读取结构，临界点搜索的参考场是AECCAR0+AECCAR2；不要把价电子密度的极值与重建全电子密度的极值混用。附件post.sh在Slurm任务中执行critic2 topology.cri > topology.cro，后处理本身使用CPU。
AUTO会从多个起点寻找密度梯度接近零的位置。具体搜索和增补起点的方法见critic2临界点手册。
3. 阅读本次临界点#
括号中的(3,−3)、(3,−1)、(3,+1)、(3,+3)分别描述核、键、环和笼类临界点的Hessian特征。不要把“非等价点数”与考虑对称倍数后的晶胞总数混淆。
| 类型 | 非等价点数 | 周期原胞中总数 |
|---|---|---|
| 核 | 2 | 2 |
| 键类鞍点 | 2 | 12 |
| 环类鞍点 | 1 | 12 |
| 笼 | 1 | 2 |
本例周期体系的Morse和为2−12+12−2=0。这是搜索结果的整体检查之一，不足以单独证明每个局部数值已收敛。
| 键径连接 | 密度 / e·bohr⁻³ | 拉普拉斯 / e·bohr⁻⁵ |
|---|---|---|
| Na–Cl | 0.01123088 | +0.05395408 |
| Cl–Cl | 0.00438811 | +0.01106495 |
Na–Cl键径长度约5.3859 bohr，即2.8501 Å。Cl–Cl键径来自这个晶体电子密度的拓扑连接；存在一条键径不等同于存在一根可赋予独立解离能的化学键，不能只按点名“bond”判断键强度。
4. 换体系时重点检查什么#
先检查结构、密度文件和参考场，再检查临界点分类、对称倍数、连接端点与梯度。拉普拉斯涉及二阶导数，比密度值更敏感，应增加FFT网格、比较插值与搜索参数；尤其不能用有限网格上的核附近密度直接推导精密核性质。
这套120³网格用于展示流程，尚未进行多网格收敛比较。对于弱相互作用和真空区域，还需辨别低密度数值临界点。不要用单个拉普拉斯正负号概括全部成键性质，可结合Bader与COHP及ELF一起分析。
输入、提取脚本与结果#
先按NaCl目录输入完成VASP计算，再提交post.sh，最后执行extract_qtaim.py提取自己的topology.cro。下载包保留参考数据和图，不包含授权势文件或大体积电荷密度。