ARTICLE DETAIL

资讯详情

深耕网站建设与运营推广的一线实战洞察。

VTST scripts 使用指南:VASP 中 NEB 过渡态计算与势垒分析

VTST scripts 使用指南:VASP 中 NEB 过渡态计算与势垒分析 简介面向VASP使用者的过渡态工具脚本集源自VTST项目用于化学反应过渡态搜索、活化能垒计算与反应路径分析适合从事第一性原理计算和催化机理研究的学生及科研人员。压缩包共152个文件以Perl脚本104个pl和Python脚本19个py为主另含Shell脚本10个sh、Perl模块7个pm及Gnuplot绘图脚本6个gnu等总大小342KB便于快速部署其中pl与pm脚本承担核心计算逻辑py脚本用于数据处理与可视化sh脚本用于任务批处理gnu脚本绘制NEB路径和振动模式图。脚本功能覆盖几何优化、频率分析、NEB最小能量路径搜索、VASP输入输出解析以及能量势垒对比等基本覆盖过渡态计算全流程。目前已有585人学习下载得到不少计算化学用户的认可。借助这套脚本研究者可自动化完成多个中间步骤并结合并行计算与参数调优辅助工具减少人为操作误差显著提升VASP过渡态计算效率与结果可重复性。1. 把vtst scripts_script_当黑匣子用是你吃亏的开始做表面催化或扩散能垒计算的人迟早会在检索框里输入vtst scripts这几个字。VTST 是给 VASP 用的过渡态计算扩展而标题里的vtst scripts_script_指的就是它配套的那一批 Perl 分析脚本生成中间镜像、提取能量、统计受力、画势垒曲线都靠它们。很多新手把它当成黑匣子以为解压就能用结果卡在“补丁没编进去、INCAR 参数漏写、输出看不懂”三件事上。这篇文章按安装、建模、运行、分析、排错、进阶的顺序来写保证你能照着跑出一个像样的能垒曲线也能在算崩的时候知道去查哪里。2. 先把 vtst scripts 装明白编译进 VASP 的前置条件与最小命令2.1 补丁和脚本是两件事别混为一谈VTST 的名称会让人误以为是一个“下载即用”的程序实际上要拆成两层来看。第一层是打进 VASP 源码的补丁代码官方发布包里常见的目录是 vtstcode5、vtstcode6它们会被复制进 VASP 的 src 目录参与重新编译。第二层是独立于 VASP 运行的 Perl 脚本nebmake.pl、nebef.pl、nebarrier.pl 都在这一层。补丁的作用是让 VASP 在电子自洽计算之外再算每个镜像沿路径方向的弹簧力和切线力这是 NEB 算法真正发生的地方脚本只负责计算前准备结构和计算后解析输出不参与势能面的求值。理解这个分层你排错时就能少走一半弯路。如果 INCAR 里写了 ICHAIN0 但 VASP 报 Unknown tag那是补丁没编译进去如果 nebarrier.pl 读出来的能量对不上那是分析脚本的运行环境或目录结构出了问题。我见过不少团队把“脚本用不了”归咎于 VTST 本身最后发现连补丁都没进 VASP白白浪费了一周。所以安装阶段的核心任务不是“把文件放对位置”而是“确认补丁确实编进了可执行文件”。2.2 先确认 VASP 大版本再选 vtstcode 目录不同 VASP 版本对 VTST 的支持方式不完全一样常见做法是VASP 5.4.x 配 vtstcode5VASP 6.x 配 vtstcode6。这个对应关系最好以你手里发布包为准因为官方会随版本维护目录最怕的是从旧教程里复制一个 vtstcode4 过来硬编进新源码。先在你的 VASP 源码目录里确认一下大版本可以用下面的命令# 进入 VASP 源码根目录看 main.F 里的版本标记 grep VASP src/main.F | head -1这段命令的作用只是帮你确认手里的 VASP 是从哪个版本展开的避免补丁目录选错。选型的逻辑是VTST 补丁和 VASP 主程序紧密耦合连 FORTRAN 的模块接口都依赖具体版本所以别指望跨版本通用。如果你同时维护两套 VASP建议把编译好的可执行文件改名区分比如 vasp_std_541_vtst别直接覆盖无补丁版本。分析脚本则不需要严格匹配版本因为读的是 OUTCAR 文本跨版本相对宽容但也要避免用太老的脚本去解析新版输出格式。如果你在编译时遇到未定义符号优先怀疑 VASP 主版本与 vtstcode 目录不对应其次是 makefile 的链接库没有重新生成。很多集群上 VASP 用多个 makefile 入口切换编译器后必须重新链接否则 VTST 的模块已经被编译但主程序仍挂在旧目标文件上。这个坑特别隐蔽因为 make 不会提示任何 VTST 相关错误。2.3 编译补丁的最小命令序列这里的操作步骤以我自己的习惯为例先解压 VTST把对应的 vtstcode 目录里的源文件复制进 VASP 源码再按发布包里的 README 修改 main.F最后重新编译。# 假设 VTST_HOME 是解压后的 VTST 目录VASP_SRC 是 VASP 源码目录 cp ${VTST_HOME}/vtstcode5/*.F ${VASP_SRC}/src/ cp ${VTST_HOME}/vtstcode5/*.H ${VASP_SRC}/src/ # 编辑 src/main.F在电子结构解析循环之前加入 NEB 调用 # 具体插入哪一行要以你手里 README 中的 diff 说明为准 cd ${VASP_SRC} make clean make -j8这段命令的逻辑是把 VTST 的 FORTRAN 模块直接混进 src 目录再通过 VASP 的 makefile 一并编译链接。命令里最需要留神的是 main.F 的修改不要凭记忆加同一个 VASP 版本在不同补丁发布版里的插入位置可能不一样。参数说明make -j8中的 8 是你愿意用于编译的进程数如果机器内存小于 16GB我建议降到 4否则某个 FORTRAN 编译进程吃满内存容易在链接阶段突然退出。分析脚本不需要编译但建议单独放一个目录并加入 PATH避免每次敲路径# 把 vtstscripts 目录加进环境变量 export VTST_SCRIPTS${VTST_HOME}/vtstscripts export PATH${VTST_SCRIPTS}:${PATH} # 检查 Perl 语法能跑到这一步说明脚本文件没损坏 perl -c ${VTST_SCRIPTS}/nebef.plperl -c只做语法检查不实际解析 OUTCAR因此它不报错只能说明脚本能加载环境。真正的验证要到编译验证那一步才有意义。常见做法是把这个 export 写进~/.bashrc但我更建议在每一个 NEB 项目目录里写一个小的env.sh因为你可能同时存在多个 VASP 版本全局 PATH 太宽松反而容易调错脚本。2.4 编译后怎么确认补丁生效跑一次最小 NEB别用 grep 去二进制里找“VTST”字符串那种验证不可靠。最省时间的方法是直接构造一个最小 NEB 算例跑几十秒就知道补丁在不在。下面是一份极简 INCAR用很小的晶胞、两个很接近的结构来测试# 在测试目录里放两个几乎一样的 POSCAR再从父目录提交 cat INCAR EOF SYSTEM vtst_test ENCUT 300 EDIFF 1E-4 IBRION 2 ISIF 0 NSW 2 ICHAIN 0 IMAGES 2 POTIM 0.0 EOF mpirun -np 4 /path/to/vasp_std test.log这段配置的逻辑是ICHAIN0 明确要求 VTST 的 NEB 驱动IMAGES2 让工作量足够小如果补丁没编译进去VASP 会直接忽略 ICHAIN 或者在日志里报 Unknown tag。参数说明POTIM0.0 在这里不是笔误NEB 阶段不依赖 VASP 的离子步长VTST 会用自己的优化器控制步长EDIFF 不需要设得太小测试补丁阶段 1E-4 足够。跑完后看 test.log如果出现镜像数相关的统计输出说明补丁生效可以开始正式建模如果连电子步都没走完就退出先去确认上一步的编译过程。集群调度场景下不需要交互式终端把上面的 mpirun 写进作业脚本即可。注意不要在作业脚本里同时使用调度器环境和 mpirun 的-np抢占同一批核数否则 VTST 的镜像并行会看到多于实际核数的 CPU导致进程分配混乱。3. 让 vtst scripts 替你生成路径从两个 POSCAR 到一条 NEB 链3.1 先处理初末态结构一致性和弛豫口径NEB 的前提是初态和末态已经各自优化到位否则整条路径都会往错误方向漂。准备两个端点结构时有三个一致性要求几乎是铁律晶格常数必须完全一致原子顺序必须完全一致固定原子必须一致。常见做法是从弛豫好的 CONTCAR 里复制而不是手动调整 POSCAR如果初态和末态用了不同尺寸的晶胞nebmake.pl 不会帮你重排只会按坐标一一对应出来的镜像原子坐标会像乱码。取端点文件的命令很简单# 初态结构从第一个弛豫完成的 CONTCAR 复制 cp ../relax_react/CONTCAR POSCAR_init # 末态结构从第二个弛豫完成的 CONTCAR 复制 cp ../relax_prod/CONTCAR POSCAR_final强制固定表面底层的做法也建议在生成镜像之前做而不是生成之后再逐个目录改。你可以在两个端点 POSCAR 里把底两层原子坐标标记为T然后让脚本带着标记一起插值大多数情况下脚本能把T保留下来但你仍然要抽查中间目录因为少数脚本版本会在插值后丢掉选择性动力学标记。这个细节直接决定后面能不能收敛不值得省。3.2 用 nebmake.pl 生成镜像目录并核对编号准备一个干净的目录把两个端点 POSCAR 放进去然后执行nebmake.pl POSCAR_init POSCAR_final 4 ls -Fnebmake.pl 是 vtst scripts 里最常用的生成器参数 4 表示中间镜像数量。它会在当前目录下生成 00 到 05 共六个目录其中 00 是初态05 是终态01、02、03、04 是要跑的中间镜像。注意INCAR 里的 IMAGES 必须和这里的 4 一致不是六个目录都算镜像。如果ls -F看到多了一个backup之类的目录先移走再提交VTST 按目录名遍历多余目录会干扰统计。逻辑上脚本用的是线性插值对每个原子的分数坐标做差值然后乘上晶格生成新的 POSCAR。对扩散这类初末态位移较小的体系线性插值够用对涉及分子转动或化学键重组的大位移过程线性插值容易让中间镜像原子重叠导致能量出现上万 eV 的尖峰。遇到这种情况不要急着改 INCAR先检查 02 或 03 的 POSCAR 里是不是有原子对距离小于 1 埃。常见做法是改用 IDPP 方法生成初猜路径或者临时增加镜像数让每段位移更小。如果你想快速扫描所有中间目录的原子间距可以用 ASE 跑一小段 Pythonpython3 - EOF import ase.io, numpy as np for i in range(1, 5): atoms ase.io.read(f{i:02d}/POSCAR) d atoms.get_all_distances(micTrue) d d[d 0] print(i, round(d.min(), 3)) EOF这段代码会打印每个中间镜像的最小原子间距如果某个值小于 0.8说明插值已经有原子重叠风险。参数说明micTrue是让 ASE 按最小镜像约定计算距离避免周期性边界造成的假近距离读入的是 POSCAR 而不是 CONTCAR因为初猜路径只看插值结构。3.3 INCAR 参数怎么和脚本生成的目录配合把下面的参数块写进父目录的 INCAR并去掉 VASP 常规几何优化的 IBRION/POTIM 组合ICHAIN 0 LCLIMB .TRUE. IOPT 7 LDNEB .TRUE. SPRING -5.0 IMAGES 4这里面每个参数都有自己的任务。ICHAIN0 是总开关用来告诉 VASP 调用 VTST 的 NEB 模块LCLIMB.TRUE. 表示使用 CI-NEB也就是爬坡版本让最高镜像沿势能面爬向鞍点IOPT7 选择 L-BFGS 优化器这是我在多数金属表面体系里的默认选择收敛稳定LDNEB.TRUE. 开启 NEB 切线力投影SPRING-5.0 是相邻镜像之间弹簧常数单位是 eV/埃平方IMAGES4 必须与 nebmake.pl 生成时的数字一致否则目录数量对不上。如果这是一条很陡的路径我通常把 IOPT 改成 3也就是 FIRE。FIRE 对初猜的敏感度更低不容易在第一步就把镜像甩飞缺点是收敛后期慢。更稳妥的流程是先 LCLIMB.FALSE. 跑几百步让路径形状稳定下来看 nebef.dat 的大趋势再把 LCLIMB 改成 .TRUE. 继续跑。注意重新跑之前把每个子目录的 WAVECAR 删掉或转移不要让上一阶段的波函数污染新阶段的电子步。3.4 提交方式在父目录一次跑完不是循环进子目录很多第一次用 VTST 的人会下意识地对 00 到 05 循环运行 VASP这是一个破坏性操作。NEB 的弹簧力必须由主进程统一计算子目录里的独立 VASP 算出来的只是孤立构型的能量完全不是过渡态路径。正确做法是在包含这些子目录的父目录提交一次作业# 父目录包含 00~05 六个子目录 mpirun -np 16 /path/to/vasp_std neb_run.logVTST 会把 16 个进程按镜像和 k 点做两级并行。这里的核数选择有讲究中间镜像为 4 时我习惯用 8 或 12 核让计算进程尽量平均分到 6 个结构上而不是堆在某个镜像里如果你只有一个大节点也可以把 NPAR 显式设为 1避免 k 点并行与镜像并行互相嵌套。运行中判断进度的方法不是盯着父目录的 log而是看每个子目录里的 OUTCAR 电子步# 看某个镜像最新一步的总能 tail -3 03/OUTCAR grep energy 03/OUTCAR | tail -1上面两行的逻辑是tail -3看这个镜像当前有没有在推进grep energy看它最新的能量值。如果03的能量连续多个小时不变而其他镜像还在更新多半是负载不均或进程卡死按第 5 章的 5.5 节处理。4. 算完之后的 vtst scripts 分析能量、受力、势垒曲线和轨迹4.1 用 nebef.pl 提取每个镜像的能量和受力NEB 跑完后最原始的数据分散在 00 到 05 各子目录的 OUTCAR 里。建议别用 grep 去一个个找而是直接信任 vtst scripts 里的 nebef.pl# 在父目录执行自动扫描所有子目录 nebef.pl nebef.dat head -8 nebef.datnebef.pl 的逻辑是按 VTST 约定的目录名依次进入读取每个 OUTCAR 的总能和原子受力然后按镜像顺序汇总输出。它的输出列通常包括镜像编号、能量以及各方向的受力分量这正好是后面画势垒曲线和检查收敛的原料。参数说明命令不加参数直接重定向是最常见用法如果某个子目录的 OUTCAR 不完整脚本不会停下来报错而是少给一行数据所以 head 之后要数一下行数够不够镜像数加表头。很多人的第二反应是把 nebef.dat 里的能量绝对值当成能垒其实不对。nebef.pl 输出的每个镜像的电子总能是一个带有随机参考的绝对值能垒必须由 nebarrier.pl 或自己取相对差值得到。所以当你看到-142.35 eV之类的数字时不要奇怪它为什么没有落在你心里预期的范围绝对值没有意义镜像之间的差才有意义。另一个常见误用是NEB 还在跑就跑 nebef.pl 看中间结果。如果某个子目录的 OUTCAR 正好写到一半脚本会把这一行也读进去造成这个镜像能量异常低或异常高。我的习惯是至少等父目录 log 里出现终止标志或者直接等作业调度系统显示完成再做提取。4.2 用 nebarrier.pl 判断能垒和鞍点位置nebef.dat 只是一堆数字nebarrier.pl 负责把它们变成能垒结论。命令很简单nebarrier.pl这个脚本会读各镜像能量找到最高点和相对最低能量端点之间的差值并指出鞍点位于哪个镜像附近。它的价值在于让你快速判断这条路径是否可信如果过渡态镜像在路径中间且左右两边能量对称下降那是健康的 NEB 结果如果最高镜像落在端点说明初态或末态根本没优化好或者是 IMAGES 太少真实鞍点没有被采样到。需要注意单位和小数位。VTST 的输出通常以 eV 为单位但能垒差可能只有零点零几 eV如果 INCAR 里 EDIFF 设成 1E-3那能垒的精度就会被电子步收敛误差吃掉。所以正式 NEB 收敛阶段我一般把 EDIFF 降到 1E-5 或更低虽然多花时间但 nebarrier.pl 报出来的小数点后两位才有意义。如果你发现能垒值在两次相同设置的计算里相差 0.1 eV先检查电子步收敛标准而不是怀疑脚本算错。如果路径出现两个相近的峰nebarrier.pl 只会给你最高一个。这时候我会把 nebef.dat 导入绘图软件看整体形状两个峰中间夹着一个明显低谷说明真实过渡态不是简单一跳中间很可能有个亚稳中间体。继续往下做之前先确认初态、中间体和末态各自都有独立的结构优化迂回路径跑出来的单一能垒一般不接受。4.3 用 nebmovie.pl 和 gnuplot 把结果变成能放进报告的样子数值之外评审最常看的是曲线和路径动画。曲线部分直接用 nebef.dat 喂给 gnuplot 即可plot nebef.dat using 1:2 with linespoints title NEB path上面这行代码要求 nebef.dat 的第一列是镜像编号、第二列是能量。如果列顺序不同先用 head 看清表头再改 using 参数。注意 gnuplot 默认会把端点 00 和 05 也画进去这没关系势垒曲线本来就该包含端点。轨迹动画方面vtst scripts 里常见的是 nebmovie.pl它会把所有镜像的结构按路径顺序合成一个文件然后用 VESTA 或类似工具播放。做扩散体系时我一般直接看动画里最高镜像附近的原子摆动方向这比纯数字更能发现“路径中间藏了一个不该有的过渡态”。4.4 一个快速检查所有镜像收敛状态的组合命令有时候不想等完整脚本只想知道谁还没收敛。可以用一个循环直接看各子目录最后一步的最大受力for i in 00 01 02 03 04 05; do printf %s $i grep FORCES: max atom $i/OUTCAR | tail -1 done这段命令的逻辑是每个 OUTCAR 里会周期性地输出当前步的最大原子受力tail -1取最后一行就是最新状态。它不会替代 nebef.pl 的完整统计但能帮你秒级定位是哪个镜像卡住如果 03 的受力是其他镜像的十倍单独看 03 的 CONTCAR 和 OUTCAR 会更有效率。参数说明镜像编号按你实际的目录数量改别用0*通配符因为多出来的备份目录会被一起扫进去。5. vtst scripts 使用中最常见的五个翻车现场下面五条都是我在真实项目里反复遇到过的每条按现象、原因、解决给出。你不用同时排查所有项先看现象对号入座通常能省掉半天试错。5.1 编完 VTST 后 INCAR 写 ICHAIN 却报 Unknown tag现象VASP 运行几秒就退出日志里出现Unknown tag ICHAIN。 原因VTST 补丁没有真正编译进可执行程序。最常见的是 src 里拷贝了 vtstcode 的 .F 文件但忘了改 main.F或者改了 main.F 却在 make 时用了旧的目标文件。 解决回到第二章的编译命令先make clean再重新编译。不要把新的可执行文件直接覆盖到系统公共目录建议单独存成vasp_std_vtst然后用 2.4 节的最小 NEB 测试一遍确认这个文件名对应的版本是带补丁的。5.2 nebmake.pl 生成的中间镜像原子重叠第一步能量上万现象第一个 SCC 循环还没结束某个镜像的能量冲到几千甚至上万 eV然后几何优化直接发散。 原因线性插值对化学键重新取向类过程不友好两个原子在中间被插到几乎重合的位置产生的短程排斥力把优化器逼疯。 解决先检查对应目录 POSCAR找到距离小于 1 埃的原子对。然后在初末期结构之间改用 IDPP 路径生成如果当前版本的 vtst scripts 里没有现成工具可以临时增加镜像数比如从 4 改成 6 或 8把每段位移摊薄。更保守的办法是先把初态到末态拆成两段分别做 NEB再接成一条长路径。5.3 力不下降能量曲线像锯齿一样震荡现象跑了两千多步nebef.dat 里的能量曲线在某个镜像附近反复上下跳动受力收敛不到 0.02 eV/埃以下。 原因最常见是 IOPT7 在陡峭势能面入口处步子太大其次是 SPRING 设置过大弹簧力把镜像拖离势能面。 解决把 IOPT 改成 3也就是 FIRE先跑几百步看曲线是否稳定如果仍然震荡把 SPRING 从 -5.0 改到 -3.0并在 INCAR 里给优化器更多步数。注意改 SPRING 后可以保留子目录里的 CONTCAR 作为新初猜但最好把旧 OUTCAR 归档避免新旧数据混淆。5.4 所有 OUTCAR 都在nebef.pl 却输出空文件现象nebef.pl 不报错生成的 nebef.dat 只有表头或者少了几行。 原因父目录里混入了多余编号目录或者某个子目录名被改成了03_bak脚本按固定通配符扫描时漏掉了一个镜像少数情况是某子目录里根本没有 OUTCAR因为作业在进入该镜像的 SCF 前被终止。 解决先用for i in 0*; do echo $i; done列出目录确认只有 00 到 IMAGES1把多余目录移到别处再重跑 nebef.pl。如果重跑还缺单独检查缺失子目录是否有 OUTCAR没有的话补一次该镜像的静态计算或者直接重跑 NEB。别手动改 nebef.dat脚本内部的列对齐你很难补对。5.5 多镜像并行时某个镜像的 OUTCAR 长时间不动现象父日志还活着但 03 的 OUTCAR 大小在十分钟内没有变化其他镜像已更新多步。 原因NPAR 与镜像并行嵌套导致进程组分配不均或者节点之间的 MPI 通讯卡在某个镜像的 k 点循环里。 解决在 INCAR 里显式设置NPAR 1把并行维度收窄到只做镜像并行这是最稳的做法适合中小体系。如果你坚持要用多 k 点并行至少把 NCORE 和 NPAR 的乘积控制在实际可用核数之内避免超过物理核数。这个现象在跨节点并行的 NEB 里尤其常见所以我现在的习惯是 NEB 阶段不用跨节点作业顶多单节点 16 核跑慢一点也比卡死强。6. 想再往前走CI-NEB 收敛后如何用 vtst scripts 验证鞍点6.1 先看鞍点位置稳不稳定NEB 结果里最该信的不是能垒数值而是鞍点镜像编号。如果你在相同设置下跑了两次nebarrier.pl 指出的最高镜像从 03 跳到 04说明路径还没有稳定已有的能垒值只能当作参考。一个快速判断方式是看能量曲线顶部两侧如果最高点两侧各有一点近乎对称的下降说明鞍点区域采样够密如果最高点旁边紧挨着另一镜像能量只差 0.01 eV那就要加密镜像重新跑。6.2 用 nebef.dat 单独查最高镜像的受力CI-NEB 收敛后鞍点镜像的原子受力应该整体接近零切向力尤其应该接近零否则它不满足鞍点的梯度条件。用 awk 从 nebef.dat 里把对应镜像的行抽出来看# 假设最高镜像是 03第二列为能量 awk $13 {print $2, $3} nebef.dat这里$13是镜像编号条件$2和$3是按需要调整的位置列。逻辑很简单把编号 03 的行的能量和受力打印出来判断它是不是真的收敛到鞍点。如果受力里某个分量还大于 0.05 eV/埃说明 CI-NEB 还没跑透不能直接拿去算频率。参数说明awk 的列号以 nebef.pl 实际输出为准先 head 一行再套用不要硬套我的示例。6.3 最后用频率计算做判决脚本验证到头真正的判决还是频率。把最高镜像目录里的 CONTCAR 复制成新的 POSCAR放到独立目录做一次静态弛豫和频率计算如果结果里恰好有一个虚频而且虚频对应的原子振动方向正好指向初态和末态这个过渡态才成立。注意有些体系的频率计算需要修改 IBRION 和内存参数和 NEB 的设置不能直接混用。我自己的习惯是每个 NEB 跑完先把父目录和所有子目录的 OUTCAR 打包归档再开始 CI-NEB 精跑如果不归档精跑翻车后连原始路径都没法恢复等于白等几天。希望帮到你。本文还有配套的精品资源点击获取
返回列表