“NI”这个符号一出现搞计算材料的基本都明白是什么意思。镍周期表第28号元素3d过渡金属里的老常客催化、磁性、电池、合金哪哪儿都有它。但真跑到DFT平面波代码里算镍的时候事情就没那么舒坦了——3d电子局域性强和平面波基底天生不对付你要是随便拿一套赝势糊弄过去后面无论是算晶格常数、磁矩还是表面吸附能都会给你埋一堆雷。ATOMPAW生成PAW_PBE赝势这件事本质上是把“原子尺度的全电子量子力学问题”和“晶体尺度的平面波计算”做一个干净利落的接口。PAW方法Projector Augmented Wave投影缀加波是现在平面波DFT代码里最主流的方案之一。而ATOMPAW就是那套能从原子全电子计算出发手搓出PAW赝势数据的开源工具链。它不依赖哪个商业代码的闭源赝势库生成过程透明参数自己可控最重要的是你可以按需生成PAW_PBE赝势——也就是配PBE交换关联泛函的那一套。这篇文章我就用Ni当靶子完整盘一遍从电子构型分析、原子数据准备、截断半径选择到ATOMPAW实际运行、赝势质量检查和排错的全过程。适合正在折腾VASP、Quantum ESPRESSO、ABINIT这类平面波代码又对官方赝势库的“黑盒参数”不放心的人。哪怕你之前没碰过ATOMPAW只要知道一点DFT基础按着后面的思路走一遍也能搞出靠谱的Ni赝势来。1. 整体设计与核心思路为什么拿Ni开刀最合适1.1 镍的电子结构难点在哪镍的基态电子排布写出来是 [Ar] 3d8 4s2但固体镍里实际上3d和4s轨道之间的电荷转移很敏感经常出现3d9 4s1的“有效”构型。也就是说你在原子阶段选什么参考构型、初态构型都会直接影响赝势在固体环境中的行为。更深一层的问题是3d轨道在径向尺度上比4s要紧凑得多空间上又和3s、3p半芯态部分重叠。平面波展开要求所有价电子波函数在截断能范围内足够平滑可3d轨道那几根尖峰必须在PAW框架里靠原子轨道投影来恢复。这就是为什么对Ni来说传统的模守恒赝势往往需要很高的截断能容易超过600-800 eV而PAW方法可以在这个问题上松一大口气。另外Ni是有磁性的。铁磁态下Ni的磁矩大概是每原子0.6 μB左右这个磁矩对电子结构极其敏感。你用不同的赝势生成方案哪怕只是截断半径差0.1 bohr、或者芯态/价态划分差一个3p轨道最后的磁矩都可能差出0.05-0.1 μB。对很多表面催化计算来说这个差异直接决定吸附能能不能比得准实验。1.2 PAW相比模守恒与超软赝势的机制优势先快速梳理一下PAW的底层逻辑因为ATOMPAW的所有参数设置都是围绕这套逻辑展开的。模守恒赝势的思路是直接把原子核和芯态一起抹掉用一个平滑的赝势替代让价电子波函数在外面和真实波函数一模一样。但在处理3d这种强局域轨道时为了达到“一模一样”截断能往往被推得很高计算成本爆炸。超软赝势USPP走的是另一条路允许赝波函数在芯区不守恒用多一套补偿电荷来修正代价是引入额外的投影算符和电子密度修正项。PAW则可以理解成“把全电子波函数找回来的超软赝势”——它保留了一个原子轨道基组展开的变换关系真正把全电子波函数显式重建出来。对Ni这种体系PAW的好处特别明显。一个原因是3d态的化学环境变化很大从金属到氧化物到表面吸附体系d带的位置、占据都在变PAW通过原子投影通道的灵活扩展能更好保持“可转移性”。另一个原因是PAW计算出的磁矩和全电子结果更贴合超软赝势在磁矩上经常出现莫名其妙的偏离PAW就好得多。用一句话给不熟悉的人解释的话模守恒赝势是“换了一个假的势让它产生和真的一样的波函数”超软是“换了个假的波函数再用补偿电荷找补”PAW是“波函数暂时换成好算的算出结果后把真实的芯区细节装回去”。1.3 为什么固定选PBE泛函项目的名字里带着PAW_PBE那泛函就是Perdew-Burke-Ernzerhof的GGA泛函。这几乎是现在凝聚态DFT计算的默认选项没有太多争议空间。PBE对金属键、表面吸附、大部分氧化物结构都能给出一个和实验比较接近的结果。对于Ni无论是晶格常数、弹性模量还是磁矩PBE通常给得都比LDA更稳。LDA普遍高估结合能把晶格常数压得偏小Ni这种磁性金属更是重灾区。PBE虽然也略微高估晶格常数但误差范围在1%以内的情况非常多见文献可比性也强。还有一个非常实际的原因VASP自带的PAW赝势里所谓“PAW_PBE”其实就是用类似方法配合PBE生成的Quantum ESPRESSO的pslibrary里也是PBE占大头。所以你手搓一套PAW_PBE赝势出来可以直接用别人的文献参数对比验证不至于因为泛函不一致而鸡同鸭讲。ATOMPAW本身也内置了PBE的处理路径你只要在模版里选对泛函生成的径向数据和交换关联势都是配套好的不用自己额外去查交换关联的数值表。2. 生成Ni赝势前必须想清楚的四个物理参数2.1 参考态与初态怎么定PAW赝势生成的第一步是给原子作全电子计算。这里需要两套电子构型参考态配置reference configuration被用来构建PAW基函数和投影算符的那个原子态相当于“锚点”。初态配置initial configuration生成电子密度初始猜测时用的通常选一个容易收敛的离子态。对Ni常见的选择是参考态3d8 4s2初态3d9 4s1。为什么这样挑因为参考态要尽量接近晶体里的实际电子分布。固体Ni是d带接近填满、s带几乎退化的状态3d8 4s2在各个体系间是比较好折中的一个选择。初态选3d9 4s1则是因为全电子计算里这个配置比较容易自洽收敛特别在使用对数网格快速估算能级时3d9 4s1比3d8 4s2的更少出现轨道翻转。有些人会把参考态写成3d10理由是Ni在化合物里经常出现d10闭壳层状态。但这样做会让4s态在原子层面的行为偏软到金属Ni里容易出现微弱的不对劲。我自己的经验是除非你专门针对Ni(II)配合物做高精度计算否则3d8 4s2作为参考态是最省心的。2.2 芯态和半芯态的取舍在PAW生成里“哪些电子算芯态、哪些算价态”这个问题比很多人想象的更重要。Ni的原子序数是281s到2p一共10个电子这些离3d太远肯定算芯态没悬念。真正要掰扯的是3s、3p这8个电子。3s2 3p6的能量深度大约在-80到-100 eV量级不算深空间上和3d有相当大的重叠。如果把这8个电子直接扔进芯态计算会明显变软平面波截断能要求也低但后果是镍原子在压力下的行为、晶格动力学、以及一些磁性性质会出问题。因为3p轨道和3d轨道的径向重叠在压缩环境下会被“压”出不可忽略的电子云重排效应你要是把它们全当芯等于假设这套电子云永远僵硬不动。很多标题里带“PAW_PBE Ni”的赝势在VASP库里有几种版本其中常区分“3p作为价态”和“3p作为半芯态”的版本根本原因就在这。ATOMPAW里处理半芯态不算复杂就是把3s、3p也一并列入价态集合但要留意它们的径向波函数在芯区还有幅度这会加大径向网格密度要求也让最终的赝势文件大一点。配3p半芯态的Ni赝势在截断能相同的情况下会比纯3d4s价态版本稍微慢一些测量下来大概多5%-10%的计算耗时但对磁性体系的准确度收益完全值得这个代价。2.3 截断半径rc怎么选截断半径大概是整个ATOMPAW输入参数里最微妙的一个。它决定了全电子波函数在多大半径内被投影波函数替代也决定了赝势的“软硬”程度。rc取得小比如1.4 bohr全电子波函数被保留的范围大精度高更“硬”但平面波截断能要求很高。rc取得大比如2.2 bohr赝势很软平面波收敛快但会牺牲一部分对局域化学环境的敏感度表面吸附、磁性分裂这类计算容易出偏差。Ni比较常用的rc范围在1.6-2.0 bohr之间。比如VASP标准PAW_PBE里的Ni对应的等效截断能在250-300 eV附近这差不多就是rc取到1.7-1.8 bohr的感觉。如果你做高精度表面催化我会建议rc往1.6再靠一靠配合包含3p半芯态会把吸附能计算里的基组重叠误差压得更低。另外径向网格也不是越多越好。ATOMPAW默认使用对数网格起点一般在10^-5 bohr量级终点要能覆盖到rc以外足够远。对Ni网格点取300左右通常够了如果你发现生成的径向波函数在近核区出现非物理振荡优先怀疑的是网格点数不足而不是rc错了。2.4 投影通道与角动量通道的设置PAW投影算符的数量直接决定了精确度上限。每个角动量通道l至少需要一个投影函数但对于Ni的d通道强烈建议至少配两个通道。原因在3d散射势的复杂性。d轨道的形状、能量位置在金属和绝缘体中变化很大一个投影函数只能覆盖一定的能量范围两个或三个投影函数才能把散射行为更好地插值到不同化学环境。具体到ATOMPAW的参数表格里通常会用“每通道投影数目”和“局域轨道数目”这类选项Ni的设置可以参考s通道1个或2个投影p通道1-2个d通道2-3个。有人会贪多把所有通道都加满投影函数。确实更精确但PAW投影矩阵会变大计算成本上升而且有些代码对过多投影通道支持得不太好。一个实践经验是先跑一遍原子自洽看看每个投影通道的占据数是否都在合理范围内0-1之间如果某个通道的占据数突然接近2或者出现怪异的负值那就说明投影通道设置过度了。3. 实战ATOMPAW生成Ni_PBE赝势的完整流程3.1 环境准备与编译要点ATOMPAW的编译不算难主要依赖Fortran编译器和LAPACK/BLAS库。常见组合是gfortran OpenBLASIntel机器上也可以用ifort MKL后者在大数组运算时略快几秒但日常用真感觉不出来。典型的配置流程就是在源码目录里运行configure指定FC和LAPACK路径make。之后会生成可执行文件atompaw以及后处理工具genpaw用于把中间数据转换为特定代码格式。如果你只需要VASP格式的POTCAR编译时最好带上VASP相关的支持开关目标代码是Quantum ESPRESSO就选qe的输出格式ABINIT自然就是abinit格式。这些输出格式在配置阶段明确指定后面生成转换时会顺很多。实际建议是全都编上反正源码里就那几条开关多花不了两分钟。3.2 输入文件长什么样ATOMPAW的输入文件在不同版本之间略有差异但核心内容块是一致的。下面我放一个针对Ni的参考模板里面的字段顺序建议以你本机版本的官方文档为准——我写这个是让你看懂每块在干什么而不是让你无脑照抄。# ATOMPAW 输入示例Ni参考配置 3d8 4s2 # 字段说明原子名 Z 处理状态数 截断半径 Ni 28 4 1.8 # 每个态一行主量子数 n角动量 l占据数 f能量初始猜测 3 2 8 0.0 4 1 2 0.0 3 1 6 0.0 3 0 2 0.0 # 投影与局域通道设定 # 此处略——视版本不同以交互输入或后续标志段完成第一行Ni是标识第二行28是原子序数第三行的4代表你要显式处理四个价态含半芯态千万别把那10个芯态也算进去。第四行1.8是截断半径单位bohr。接下来的四行分别是3d8、4s2、3p6、3s2。很多版本里每个状态还需要初始能量填0.0通常表示让程序自动在参考态附近找初始能级。你要是遇到收敛困难再手动给一个大约的轨道能量。至于投影通道数量的给定方式每个版本处理得不一样有的放在交互提问里有的放在文件后段。我习惯先用一个保守默认s通道1个投影、p通道1个、d通道2个。后续如果做d带中心类的高精度计算再把d通道加到3个重新生成一遍对比一次便知差别。3.3 运行ATOMPAW并理解输出运行方式通常就一条命令atompaw atompaw.in或者按你版本支持的参数方式atompaw -i atompaw.in -o atompaw.out程序启动后会先做全电子原子自洽然后在此基础上构造PAW数据集。整个过程从几秒到几十秒不等如果你在HPC集群的登录节点上跑也基本无压力。运行结束后你会看到一批输出文件。核心的是原子波函数数据文件里面包含径向网格、各态的对数导数、投影函数矩阵、补偿电荷参数、交换关联势数据等。这些文件有些是给后续转换用的中间文件有些是最终可以放进代码的赝势文件。别急着打开二进制文件看——没意义。你应该看的是atompaw输出日志里那几行摘要总能、各态占据数、投影函数的范数以及径向函数的平滑度指标。只要这些指标没报警文件就是可用的。3.4 转换成VASP/QE/ABINIT格式同一套原子数据和投影参数要配合不同DFT代码还需要做格式转换。目标代码是VASP用genpaw或对应脚本生成POTCAR。VASP要求POTCAR里有完整的原子数据头、径向网格、投影函数以及PBE交换关联势信息。生成后你可以用grep VRHFIN看一眼头信息确认元素、版本号这些元数据没乱。目标代码是Quantum ESPRESSO转为UPF格式。QE的UPF文件本质是XML里面包含同样的数据集你可直接在pseudopotentials目录里引用。目标代码是ABINIT通常生成pspcode的PAW格式文件再配合abinit的输入文件使用。转换工具的具体命令每个版本有区别但不外乎输入中间文件名输出文件名。转换完成后强烈建议做一个最小计算测试一个Ni原子或一个Ni晶胞都行把截断能设到400 eV左右看看总能量能否在几个电子伏的合理范围内收敛。3.5 输入参数走向验证一个可落地的操作顺序实际工程里我建议你按这个顺序操作不要跳步先用参考态3d8 4s2跑一遍全电子原子计算确认总能和文献参考一致Ni原子的PBE总能大约在-2700到-2750 eV量级不同轨道数设置略有差异。用默认rc1.8、默认投影通道生成第一版赝势。用这版赝势算一个Ni的fcc晶胞对比实验晶格常数约3.52 Å和磁矩约0.6 μB/原子。如果偏差大先怀疑截断能不够再怀疑投影通道或rc。记录不同截断能下的总能量变化画一个收敛曲线确认你用于生产计算的截断能落在平台区。如果对计算结果精度仍不满意回到ATOMAPAW调整rc或投影通道数量重新生成。这套流程的优点是把每一次参数改动的效果都暴露在可观察的物理量上而不是在“赝势文件本身好不好看”这个层面打转。4. 赝势质量怎么检验才算过关4.1 原子级检验对数导数与全电子对比ATOMPAW生成的日志里会包含各角动量通道的对数导数曲线。这是赝势质量的第一道关卡在rc附近和几个典型能量点赝势的对数导数应该和全电子波函数的结果重合得非常好。如果曲线在某个能量点出现明显分叉说明这个通道的投影函数不足以描述散射需要增加投影或调整rc。另一个直观检验是把PAW计算出来的原子总能和真正的全电子计算对比。ATOMPAW在生成过程中本身也是在作全电子计算的所以它内部就有一组基准数据。后处理工具里一般也带原子测试功能能直接输出差值。一般来说1 meV/原子以内的差异是合格的。4.2 晶体级检验晶格常数、磁矩、体弹模量原子级检验通过后真正有说服力的是晶体结果。拿fcc Ni来说实验晶格常数约3.52 Å铁磁状态下每个原子的磁矩约0.6 μB。用PBE泛函配合靠谱的PAW赝势计算结果通常在3.50-3.54 Å之间磁矩在0.6-0.65 μB之间。假如算出来晶格常数偏到3.48以下或3.56以上说明赝势的径向行为有偏优先检查是否漏了3p半芯态、或者rc给得太软。体弹模量也可以作为辅助指标。Ni的实验体弹模量大约在186-190 GPa范围。PBE对这个量的误差普遍比晶格常数大一些但200 GPa以上就要警惕了。很多情况下体弹模量偏大是“太硬”的赝势表现这时可以稍微加大rc或者减少投影通道试试。4.3 截断能收敛曲线这个检验不值得跳过。用你生成的Ni赝势分别以300、350、400、450、500、600 eV的截断能计算fcc Ni晶胞的总能量记得固定晶格常数然后看能量差。一个合格的PAW_PBE Ni赝势在400 eV附近应当已经进入能量平台区500 eV到600 eV之间每原子能量变化应该小于0.01 eV。如果你的赝势在600 eV下还在明显漂移那多半是rc取得太小或者径向网格精度不够。很多新手的误区是“反正PAW硬截断能设很高就好了”可高截断能会拖慢计算又不带来准确度。一张收敛曲线能告诉你的不是“多少eV才能算准”而是“在多少eV之后增加截断能纯属浪费机时”。4.4 转移性检验从金属到表面再到氧化物真正高质量的PaW_PBE赝势不应该只在一个体系里表现好。Ni催化研究经常横跨fcc金属、表面、纳米颗粒、氧化物界面你的赝势最好在这些场景里都保持稳定。转移性检验最实用的做法是跑两个额外体系Ni(111)表面看表面总能是否明显低于同原子数的随机构型算一下表面能和文献值对比NiO或Ni与CO的吸附模型看吸附能是否在合理的1-2 eV范围内波动。表面上吸附体系的误差如果超过0.1 eV影响到催化趋势判断时再去回头调赝势参数。通常这时候你会发现在金属表面计算“看起来不错”的赝势在吸附体系里会暴露3d投影通道不足的问题。5. 常见错误、报错排查与实操避坑5.1 电子构型翻转 / 自洽不收敛Ni原子计算里最常见的坑就是3d和4s轨道能量接近自洽过程中电子占据容易在3d和4s之间跳来跳去。这会导致全电子原子计算的“配置翻转”最终生成赝势里的orbital occupation和你的预期不一致。出现这种情况时建议在输入文件里给每个态都填上明确的初始能量不要用0.0自动猜测。比如3d填0.3、4s填0.1之类的粗糙参考值至少让程序从接近目标解的轨道出发。另一个办法是把初态配置改成3d9 4s1这个构型自洽收敛明显更稳。5.2 投影通道数到底该不该加有一个容易被人误解的倾向投影函数越多越准。实际经验是对Ni来说s通道一个、p通道一个、d通道两个已经能覆盖绝大部分常规DFT计算。只有在非常精细的电子谱、强关联修正DFTU或者自旋极化明显的体系中d通道加到3个才能看出微弱好处。加投影通道的代价不只是计算时间还可能引入数值近似问题。投影矩阵如果病态在DFT代码里会导致电子步难以收敛、甚至出现负电子密度。我遇到输出文件里出现“ghost state”或者能带里出现明显伪带的多数都是投影通道设置过头或者rc与投影函数不匹配导致的。5.3 磁性计算的初始磁矩用新生成的Ni赝势算磁性体系时第一性原理的初始磁矩设置也很关键。对fcc Ni在DFT输入里设置初始磁矩约0.7 μB/原子通常能很快收敛到物理态。如果你设置初始磁矩0可能会掉进顺磁解磁矩收敛成0然后你又开始怀疑赝势有问题。其实赝势没问题是初始磁矩卡在了错误的势阱里。此外自旋极化计算里要对s和d通道分开看磁矩贡献。Ni的磁矩主要来自3d通道4s的磁矩贡献很小如果生成的赝势在4s上也显示出很强的自旋分裂那可能是投影通道的s通道设置有问题。5.4 与官方赝势对比时的预期差很多人手搓完一套PAW_PBE赝势立刻拿它跟VASP自带或者pslibrary的Ni赝势对比发现总能量差了几十甚至上百eV马上觉得自己生成错了。这个判断是错的。不同赝势的截止势定义、芯态处理、补偿电荷形式不一样总能量的绝对零点约定也不同。真正应该对比的不是总能量而是能量差量——晶格常数、磁矩、吸附能、能带宽度这些相对量。我的习惯是先把两种赝势算出来的物理量放一张表里如果晶格常数差在0.01 Å以内、磁矩差在0.02 μB以内这个赝势就是合格的。绝对总能差多少根本不用管。5.5 后处理转换失败的常见原因genpaw或转换脚本报错十有八九不是脚本坏了而是你安装ATOMPAW时没有把目标代码的支持模块编译进去。也就是说你在configure时没有启用VASP/QE/ABINIT的输出接口后处理工具自然不认你的目标格式。处理办法是回到configure阶段重新检查开关选项。如果确实忘记了也不用重新生成原子数据把中间文件存好重新编译一个有对应接口的ATOMPAW再转换就行。原子数据本身不依赖这些接口不需要重算。5.6 一个容易被忽略的硬件相关性ATOMPAW在部分老版本里对BLAS库的线程数有“怪癖”——你用了多线程OpenBLAS反而会因为在很小的原子计算中频繁创建线程开销大而变慢甚至偶发非确定性结果。跑单原子级别计算时我通常建议把线程数设置为1比如在运行前设置环境变量export OMP_NUM_THREADS1这样不仅稳定速度也可能反而更快。等到做大规模晶体计算、真正需要多节点并行时再把MPI/OpenMP拉满。几个真正影响Ni赝势成败的细节最后再聊几个我实际操作中踩过、后来觉得非常值得分享的细节。第一个是径向网格的“尾端”。很多人看rc小于2.0 bohr就以为径向网格到2.2 bohr足够实际上PAW数据集里有一部分形状函数和补偿电荷需要延伸到rc之外一段距离通常要到2.5-3.0 bohr。截断得太早会在实空间网格化时产生电荷振荡表现为总能出现莫名其妙的微小跳跃。第二个是参考态能量填0.0并不总是安全。虽然ATOMPAW设计了自动猜测机制但Ni这种3d-4s竞争强烈的元素我建议至少有一次尝试用手动给定能量。哪怕只是加一个非常粗糙的估计值自由度变小自洽路径反而更稳定。第三个是验证时不要只看一个体系。哪怕你的fcc Ni算得很漂亮如果打算后续做表面和吸附最好从一开始就多跑一个表面体系检验。因为表面体系对d投影通道的要求通常更高等大批量生产计算做完才发现赝势在表面问题上不可靠那才叫真难受。生成PAW_PBE赝势这件事最大的成就感不是“我跑通了”而是你能清楚说出你这套赝势在物理上做了哪些取舍参考态是什么、半芯态含不含、rc为什么取这个值、投影通道配到几级。这些决定恰恰是商业赝势库文档里经常语焉不详的地方。你自己生成一遍手里才真正有数。 
长河君 · 2026/10/14 7:29:55 · 阅读 1.2k 用ATOMPAW生成Ni的PAW_PBE赝势：从电子构型到质量检验全攻略
“NI”这个符号一出现搞计算材料的基本都明白是什么意思。镍周期表第28号元素3d过渡金属里的老常客催化、磁性、电池、合金哪哪儿都有它。但真跑到DFT平面波代码里算镍的时候事情就没那么舒坦了——3d电子局域性强和平面波基底天生不对付你要是随便拿一套赝势糊弄过去后面无论是算晶格常数、磁矩还是表面吸附能都会给你埋一堆雷。ATOMPAW生成PAW_PBE赝势这件事本质上是把“原子尺度的全电子量子力学问题”和“晶体尺度的平面波计算”做一个干净利落的接口。PAW方法Projector Augmented Wave投影缀加波是现在平面波DFT代码里最主流的方案之一。而ATOMPAW就是那套能从原子全电子计算出发手搓出PAW赝势数据的开源工具链。它不依赖哪个商业代码的闭源赝势库生成过程透明参数自己可控最重要的是你可以按需生成PAW_PBE赝势——也就是配PBE交换关联泛函的那一套。这篇文章我就用Ni当靶子完整盘一遍从电子构型分析、原子数据准备、截断半径选择到ATOMPAW实际运行、赝势质量检查和排错的全过程。适合正在折腾VASP、Quantum ESPRESSO、ABINIT这类平面波代码又对官方赝势库的“黑盒参数”不放心的人。哪怕你之前没碰过ATOMPAW只要知道一点DFT基础按着后面的思路走一遍也能搞出靠谱的Ni赝势来。1. 整体设计与核心思路为什么拿Ni开刀最合适1.1 镍的电子结构难点在哪镍的基态电子排布写出来是 [Ar] 3d8 4s2但固体镍里实际上3d和4s轨道之间的电荷转移很敏感经常出现3d9 4s1的“有效”构型。也就是说你在原子阶段选什么参考构型、初态构型都会直接影响赝势在固体环境中的行为。更深一层的问题是3d轨道在径向尺度上比4s要紧凑得多空间上又和3s、3p半芯态部分重叠。平面波展开要求所有价电子波函数在截断能范围内足够平滑可3d轨道那几根尖峰必须在PAW框架里靠原子轨道投影来恢复。这就是为什么对Ni来说传统的模守恒赝势往往需要很高的截断能容易超过600-800 eV而PAW方法可以在这个问题上松一大口气。另外Ni是有磁性的。铁磁态下Ni的磁矩大概是每原子0.6 μB左右这个磁矩对电子结构极其敏感。你用不同的赝势生成方案哪怕只是截断半径差0.1 bohr、或者芯态/价态划分差一个3p轨道最后的磁矩都可能差出0.05-0.1 μB。对很多表面催化计算来说这个差异直接决定吸附能能不能比得准实验。1.2 PAW相比模守恒与超软赝势的机制优势先快速梳理一下PAW的底层逻辑因为ATOMPAW的所有参数设置都是围绕这套逻辑展开的。模守恒赝势的思路是直接把原子核和芯态一起抹掉用一个平滑的赝势替代让价电子波函数在外面和真实波函数一模一样。但在处理3d这种强局域轨道时为了达到“一模一样”截断能往往被推得很高计算成本爆炸。超软赝势USPP走的是另一条路允许赝波函数在芯区不守恒用多一套补偿电荷来修正代价是引入额外的投影算符和电子密度修正项。PAW则可以理解成“把全电子波函数找回来的超软赝势”——它保留了一个原子轨道基组展开的变换关系真正把全电子波函数显式重建出来。对Ni这种体系PAW的好处特别明显。一个原因是3d态的化学环境变化很大从金属到氧化物到表面吸附体系d带的位置、占据都在变PAW通过原子投影通道的灵活扩展能更好保持“可转移性”。另一个原因是PAW计算出的磁矩和全电子结果更贴合超软赝势在磁矩上经常出现莫名其妙的偏离PAW就好得多。用一句话给不熟悉的人解释的话模守恒赝势是“换了一个假的势让它产生和真的一样的波函数”超软是“换了个假的波函数再用补偿电荷找补”PAW是“波函数暂时换成好算的算出结果后把真实的芯区细节装回去”。1.3 为什么固定选PBE泛函项目的名字里带着PAW_PBE那泛函就是Perdew-Burke-Ernzerhof的GGA泛函。这几乎是现在凝聚态DFT计算的默认选项没有太多争议空间。PBE对金属键、表面吸附、大部分氧化物结构都能给出一个和实验比较接近的结果。对于Ni无论是晶格常数、弹性模量还是磁矩PBE通常给得都比LDA更稳。LDA普遍高估结合能把晶格常数压得偏小Ni这种磁性金属更是重灾区。PBE虽然也略微高估晶格常数但误差范围在1%以内的情况非常多见文献可比性也强。还有一个非常实际的原因VASP自带的PAW赝势里所谓“PAW_PBE”其实就是用类似方法配合PBE生成的Quantum ESPRESSO的pslibrary里也是PBE占大头。所以你手搓一套PAW_PBE赝势出来可以直接用别人的文献参数对比验证不至于因为泛函不一致而鸡同鸭讲。ATOMPAW本身也内置了PBE的处理路径你只要在模版里选对泛函生成的径向数据和交换关联势都是配套好的不用自己额外去查交换关联的数值表。2. 生成Ni赝势前必须想清楚的四个物理参数2.1 参考态与初态怎么定PAW赝势生成的第一步是给原子作全电子计算。这里需要两套电子构型参考态配置reference configuration被用来构建PAW基函数和投影算符的那个原子态相当于“锚点”。初态配置initial configuration生成电子密度初始猜测时用的通常选一个容易收敛的离子态。对Ni常见的选择是参考态3d8 4s2初态3d9 4s1。为什么这样挑因为参考态要尽量接近晶体里的实际电子分布。固体Ni是d带接近填满、s带几乎退化的状态3d8 4s2在各个体系间是比较好折中的一个选择。初态选3d9 4s1则是因为全电子计算里这个配置比较容易自洽收敛特别在使用对数网格快速估算能级时3d9 4s1比3d8 4s2的更少出现轨道翻转。有些人会把参考态写成3d10理由是Ni在化合物里经常出现d10闭壳层状态。但这样做会让4s态在原子层面的行为偏软到金属Ni里容易出现微弱的不对劲。我自己的经验是除非你专门针对Ni(II)配合物做高精度计算否则3d8 4s2作为参考态是最省心的。2.2 芯态和半芯态的取舍在PAW生成里“哪些电子算芯态、哪些算价态”这个问题比很多人想象的更重要。Ni的原子序数是281s到2p一共10个电子这些离3d太远肯定算芯态没悬念。真正要掰扯的是3s、3p这8个电子。3s2 3p6的能量深度大约在-80到-100 eV量级不算深空间上和3d有相当大的重叠。如果把这8个电子直接扔进芯态计算会明显变软平面波截断能要求也低但后果是镍原子在压力下的行为、晶格动力学、以及一些磁性性质会出问题。因为3p轨道和3d轨道的径向重叠在压缩环境下会被“压”出不可忽略的电子云重排效应你要是把它们全当芯等于假设这套电子云永远僵硬不动。很多标题里带“PAW_PBE Ni”的赝势在VASP库里有几种版本其中常区分“3p作为价态”和“3p作为半芯态”的版本根本原因就在这。ATOMPAW里处理半芯态不算复杂就是把3s、3p也一并列入价态集合但要留意它们的径向波函数在芯区还有幅度这会加大径向网格密度要求也让最终的赝势文件大一点。配3p半芯态的Ni赝势在截断能相同的情况下会比纯3d4s价态版本稍微慢一些测量下来大概多5%-10%的计算耗时但对磁性体系的准确度收益完全值得这个代价。2.3 截断半径rc怎么选截断半径大概是整个ATOMPAW输入参数里最微妙的一个。它决定了全电子波函数在多大半径内被投影波函数替代也决定了赝势的“软硬”程度。rc取得小比如1.4 bohr全电子波函数被保留的范围大精度高更“硬”但平面波截断能要求很高。rc取得大比如2.2 bohr赝势很软平面波收敛快但会牺牲一部分对局域化学环境的敏感度表面吸附、磁性分裂这类计算容易出偏差。Ni比较常用的rc范围在1.6-2.0 bohr之间。比如VASP标准PAW_PBE里的Ni对应的等效截断能在250-300 eV附近这差不多就是rc取到1.7-1.8 bohr的感觉。如果你做高精度表面催化我会建议rc往1.6再靠一靠配合包含3p半芯态会把吸附能计算里的基组重叠误差压得更低。另外径向网格也不是越多越好。ATOMPAW默认使用对数网格起点一般在10^-5 bohr量级终点要能覆盖到rc以外足够远。对Ni网格点取300左右通常够了如果你发现生成的径向波函数在近核区出现非物理振荡优先怀疑的是网格点数不足而不是rc错了。2.4 投影通道与角动量通道的设置PAW投影算符的数量直接决定了精确度上限。每个角动量通道l至少需要一个投影函数但对于Ni的d通道强烈建议至少配两个通道。原因在3d散射势的复杂性。d轨道的形状、能量位置在金属和绝缘体中变化很大一个投影函数只能覆盖一定的能量范围两个或三个投影函数才能把散射行为更好地插值到不同化学环境。具体到ATOMPAW的参数表格里通常会用“每通道投影数目”和“局域轨道数目”这类选项Ni的设置可以参考s通道1个或2个投影p通道1-2个d通道2-3个。有人会贪多把所有通道都加满投影函数。确实更精确但PAW投影矩阵会变大计算成本上升而且有些代码对过多投影通道支持得不太好。一个实践经验是先跑一遍原子自洽看看每个投影通道的占据数是否都在合理范围内0-1之间如果某个通道的占据数突然接近2或者出现怪异的负值那就说明投影通道设置过度了。3. 实战ATOMPAW生成Ni_PBE赝势的完整流程3.1 环境准备与编译要点ATOMPAW的编译不算难主要依赖Fortran编译器和LAPACK/BLAS库。常见组合是gfortran OpenBLASIntel机器上也可以用ifort MKL后者在大数组运算时略快几秒但日常用真感觉不出来。典型的配置流程就是在源码目录里运行configure指定FC和LAPACK路径make。之后会生成可执行文件atompaw以及后处理工具genpaw用于把中间数据转换为特定代码格式。如果你只需要VASP格式的POTCAR编译时最好带上VASP相关的支持开关目标代码是Quantum ESPRESSO就选qe的输出格式ABINIT自然就是abinit格式。这些输出格式在配置阶段明确指定后面生成转换时会顺很多。实际建议是全都编上反正源码里就那几条开关多花不了两分钟。3.2 输入文件长什么样ATOMPAW的输入文件在不同版本之间略有差异但核心内容块是一致的。下面我放一个针对Ni的参考模板里面的字段顺序建议以你本机版本的官方文档为准——我写这个是让你看懂每块在干什么而不是让你无脑照抄。# ATOMPAW 输入示例Ni参考配置 3d8 4s2 # 字段说明原子名 Z 处理状态数 截断半径 Ni 28 4 1.8 # 每个态一行主量子数 n角动量 l占据数 f能量初始猜测 3 2 8 0.0 4 1 2 0.0 3 1 6 0.0 3 0 2 0.0 # 投影与局域通道设定 # 此处略——视版本不同以交互输入或后续标志段完成第一行Ni是标识第二行28是原子序数第三行的4代表你要显式处理四个价态含半芯态千万别把那10个芯态也算进去。第四行1.8是截断半径单位bohr。接下来的四行分别是3d8、4s2、3p6、3s2。很多版本里每个状态还需要初始能量填0.0通常表示让程序自动在参考态附近找初始能级。你要是遇到收敛困难再手动给一个大约的轨道能量。至于投影通道数量的给定方式每个版本处理得不一样有的放在交互提问里有的放在文件后段。我习惯先用一个保守默认s通道1个投影、p通道1个、d通道2个。后续如果做d带中心类的高精度计算再把d通道加到3个重新生成一遍对比一次便知差别。3.3 运行ATOMPAW并理解输出运行方式通常就一条命令atompaw atompaw.in或者按你版本支持的参数方式atompaw -i atompaw.in -o atompaw.out程序启动后会先做全电子原子自洽然后在此基础上构造PAW数据集。整个过程从几秒到几十秒不等如果你在HPC集群的登录节点上跑也基本无压力。运行结束后你会看到一批输出文件。核心的是原子波函数数据文件里面包含径向网格、各态的对数导数、投影函数矩阵、补偿电荷参数、交换关联势数据等。这些文件有些是给后续转换用的中间文件有些是最终可以放进代码的赝势文件。别急着打开二进制文件看——没意义。你应该看的是atompaw输出日志里那几行摘要总能、各态占据数、投影函数的范数以及径向函数的平滑度指标。只要这些指标没报警文件就是可用的。3.4 转换成VASP/QE/ABINIT格式同一套原子数据和投影参数要配合不同DFT代码还需要做格式转换。目标代码是VASP用genpaw或对应脚本生成POTCAR。VASP要求POTCAR里有完整的原子数据头、径向网格、投影函数以及PBE交换关联势信息。生成后你可以用grep VRHFIN看一眼头信息确认元素、版本号这些元数据没乱。目标代码是Quantum ESPRESSO转为UPF格式。QE的UPF文件本质是XML里面包含同样的数据集你可直接在pseudopotentials目录里引用。目标代码是ABINIT通常生成pspcode的PAW格式文件再配合abinit的输入文件使用。转换工具的具体命令每个版本有区别但不外乎输入中间文件名输出文件名。转换完成后强烈建议做一个最小计算测试一个Ni原子或一个Ni晶胞都行把截断能设到400 eV左右看看总能量能否在几个电子伏的合理范围内收敛。3.5 输入参数走向验证一个可落地的操作顺序实际工程里我建议你按这个顺序操作不要跳步先用参考态3d8 4s2跑一遍全电子原子计算确认总能和文献参考一致Ni原子的PBE总能大约在-2700到-2750 eV量级不同轨道数设置略有差异。用默认rc1.8、默认投影通道生成第一版赝势。用这版赝势算一个Ni的fcc晶胞对比实验晶格常数约3.52 Å和磁矩约0.6 μB/原子。如果偏差大先怀疑截断能不够再怀疑投影通道或rc。记录不同截断能下的总能量变化画一个收敛曲线确认你用于生产计算的截断能落在平台区。如果对计算结果精度仍不满意回到ATOMAPAW调整rc或投影通道数量重新生成。这套流程的优点是把每一次参数改动的效果都暴露在可观察的物理量上而不是在“赝势文件本身好不好看”这个层面打转。4. 赝势质量怎么检验才算过关4.1 原子级检验对数导数与全电子对比ATOMPAW生成的日志里会包含各角动量通道的对数导数曲线。这是赝势质量的第一道关卡在rc附近和几个典型能量点赝势的对数导数应该和全电子波函数的结果重合得非常好。如果曲线在某个能量点出现明显分叉说明这个通道的投影函数不足以描述散射需要增加投影或调整rc。另一个直观检验是把PAW计算出来的原子总能和真正的全电子计算对比。ATOMPAW在生成过程中本身也是在作全电子计算的所以它内部就有一组基准数据。后处理工具里一般也带原子测试功能能直接输出差值。一般来说1 meV/原子以内的差异是合格的。4.2 晶体级检验晶格常数、磁矩、体弹模量原子级检验通过后真正有说服力的是晶体结果。拿fcc Ni来说实验晶格常数约3.52 Å铁磁状态下每个原子的磁矩约0.6 μB。用PBE泛函配合靠谱的PAW赝势计算结果通常在3.50-3.54 Å之间磁矩在0.6-0.65 μB之间。假如算出来晶格常数偏到3.48以下或3.56以上说明赝势的径向行为有偏优先检查是否漏了3p半芯态、或者rc给得太软。体弹模量也可以作为辅助指标。Ni的实验体弹模量大约在186-190 GPa范围。PBE对这个量的误差普遍比晶格常数大一些但200 GPa以上就要警惕了。很多情况下体弹模量偏大是“太硬”的赝势表现这时可以稍微加大rc或者减少投影通道试试。4.3 截断能收敛曲线这个检验不值得跳过。用你生成的Ni赝势分别以300、350、400、450、500、600 eV的截断能计算fcc Ni晶胞的总能量记得固定晶格常数然后看能量差。一个合格的PAW_PBE Ni赝势在400 eV附近应当已经进入能量平台区500 eV到600 eV之间每原子能量变化应该小于0.01 eV。如果你的赝势在600 eV下还在明显漂移那多半是rc取得太小或者径向网格精度不够。很多新手的误区是“反正PAW硬截断能设很高就好了”可高截断能会拖慢计算又不带来准确度。一张收敛曲线能告诉你的不是“多少eV才能算准”而是“在多少eV之后增加截断能纯属浪费机时”。4.4 转移性检验从金属到表面再到氧化物真正高质量的PaW_PBE赝势不应该只在一个体系里表现好。Ni催化研究经常横跨fcc金属、表面、纳米颗粒、氧化物界面你的赝势最好在这些场景里都保持稳定。转移性检验最实用的做法是跑两个额外体系Ni(111)表面看表面总能是否明显低于同原子数的随机构型算一下表面能和文献值对比NiO或Ni与CO的吸附模型看吸附能是否在合理的1-2 eV范围内波动。表面上吸附体系的误差如果超过0.1 eV影响到催化趋势判断时再去回头调赝势参数。通常这时候你会发现在金属表面计算“看起来不错”的赝势在吸附体系里会暴露3d投影通道不足的问题。5. 常见错误、报错排查与实操避坑5.1 电子构型翻转 / 自洽不收敛Ni原子计算里最常见的坑就是3d和4s轨道能量接近自洽过程中电子占据容易在3d和4s之间跳来跳去。这会导致全电子原子计算的“配置翻转”最终生成赝势里的orbital occupation和你的预期不一致。出现这种情况时建议在输入文件里给每个态都填上明确的初始能量不要用0.0自动猜测。比如3d填0.3、4s填0.1之类的粗糙参考值至少让程序从接近目标解的轨道出发。另一个办法是把初态配置改成3d9 4s1这个构型自洽收敛明显更稳。5.2 投影通道数到底该不该加有一个容易被人误解的倾向投影函数越多越准。实际经验是对Ni来说s通道一个、p通道一个、d通道两个已经能覆盖绝大部分常规DFT计算。只有在非常精细的电子谱、强关联修正DFTU或者自旋极化明显的体系中d通道加到3个才能看出微弱好处。加投影通道的代价不只是计算时间还可能引入数值近似问题。投影矩阵如果病态在DFT代码里会导致电子步难以收敛、甚至出现负电子密度。我遇到输出文件里出现“ghost state”或者能带里出现明显伪带的多数都是投影通道设置过头或者rc与投影函数不匹配导致的。5.3 磁性计算的初始磁矩用新生成的Ni赝势算磁性体系时第一性原理的初始磁矩设置也很关键。对fcc Ni在DFT输入里设置初始磁矩约0.7 μB/原子通常能很快收敛到物理态。如果你设置初始磁矩0可能会掉进顺磁解磁矩收敛成0然后你又开始怀疑赝势有问题。其实赝势没问题是初始磁矩卡在了错误的势阱里。此外自旋极化计算里要对s和d通道分开看磁矩贡献。Ni的磁矩主要来自3d通道4s的磁矩贡献很小如果生成的赝势在4s上也显示出很强的自旋分裂那可能是投影通道的s通道设置有问题。5.4 与官方赝势对比时的预期差很多人手搓完一套PAW_PBE赝势立刻拿它跟VASP自带或者pslibrary的Ni赝势对比发现总能量差了几十甚至上百eV马上觉得自己生成错了。这个判断是错的。不同赝势的截止势定义、芯态处理、补偿电荷形式不一样总能量的绝对零点约定也不同。真正应该对比的不是总能量而是能量差量——晶格常数、磁矩、吸附能、能带宽度这些相对量。我的习惯是先把两种赝势算出来的物理量放一张表里如果晶格常数差在0.01 Å以内、磁矩差在0.02 μB以内这个赝势就是合格的。绝对总能差多少根本不用管。5.5 后处理转换失败的常见原因genpaw或转换脚本报错十有八九不是脚本坏了而是你安装ATOMPAW时没有把目标代码的支持模块编译进去。也就是说你在configure时没有启用VASP/QE/ABINIT的输出接口后处理工具自然不认你的目标格式。处理办法是回到configure阶段重新检查开关选项。如果确实忘记了也不用重新生成原子数据把中间文件存好重新编译一个有对应接口的ATOMPAW再转换就行。原子数据本身不依赖这些接口不需要重算。5.6 一个容易被忽略的硬件相关性ATOMPAW在部分老版本里对BLAS库的线程数有“怪癖”——你用了多线程OpenBLAS反而会因为在很小的原子计算中频繁创建线程开销大而变慢甚至偶发非确定性结果。跑单原子级别计算时我通常建议把线程数设置为1比如在运行前设置环境变量export OMP_NUM_THREADS1这样不仅稳定速度也可能反而更快。等到做大规模晶体计算、真正需要多节点并行时再把MPI/OpenMP拉满。几个真正影响Ni赝势成败的细节最后再聊几个我实际操作中踩过、后来觉得非常值得分享的细节。第一个是径向网格的“尾端”。很多人看rc小于2.0 bohr就以为径向网格到2.2 bohr足够实际上PAW数据集里有一部分形状函数和补偿电荷需要延伸到rc之外一段距离通常要到2.5-3.0 bohr。截断得太早会在实空间网格化时产生电荷振荡表现为总能出现莫名其妙的微小跳跃。第二个是参考态能量填0.0并不总是安全。虽然ATOMPAW设计了自动猜测机制但Ni这种3d-4s竞争强烈的元素我建议至少有一次尝试用手动给定能量。哪怕只是加一个非常粗糙的估计值自由度变小自洽路径反而更稳定。第三个是验证时不要只看一个体系。哪怕你的fcc Ni算得很漂亮如果打算后续做表面和吸附最好从一开始就多跑一个表面体系检验。因为表面体系对d投影通道的要求通常更高等大批量生产计算做完才发现赝势在表面问题上不可靠那才叫真难受。生成PAW_PBE赝势这件事最大的成就感不是“我跑通了”而是你能清楚说出你这套赝势在物理上做了哪些取舍参考态是什么、半芯态含不含、rc为什么取这个值、投影通道配到几级。这些决定恰恰是商业赝势库文档里经常语焉不详的地方。你自己生成一遍手里才真正有数。