VaspCZ：一个提高效率的VASP计算辅助程序

张正德 谈蒙露 任翠兰 怀 平,3,4

1（中国科学院上海应用物理研究所 上海 201800）

2（中国科学院大学 北京 100049）

3（上海科技大学物质科学与技术学院 上海 201210）

4（中国科学院上海高等研究院上海光源科学中心 上海 201204）

ⅤASP（Ⅴienna Ab-initio Simulation Package）是20世纪90年代，由Ⅴienna大学Hafner小组开发的、进行材料电子结构计算和量子力学-分子动力学模拟的商用软件［1-3］。ⅤASP基于密度泛函理论对Kohn-Shan方程进行迭代求解，从头计算求解薛定谔方程，在计算材料的结构参数、状态方程、电子性质、光学性质、磁学性质等各项实践中，预测了很多先于实验的正确结果，取得了很大的成功。目前，该程序包已更新到5.4.4版，在其基础上，世界各地的科学家也开发了诸多补充算法和辅助程序，以扩展其计算能力、提高计算精度或提升计算效率。

2000 年 ，Henkelman 等［4］发 布 了 ⅤTST（Transition State Tools forⅤASP）工具包，提供了寻找鞍点和计算过渡态的源码和脚本。作为ⅤASP的计算能力的补充，其中的CⅠ-NEB方法［5-6］和Dimer方法［7］，已经成为目前计算体相扩散、晶界扩散、催化反应的重要工具之一。Mathew等［8］于2014年在平面波DFT编码ⅤASP中实现了隐式溶剂化模型，描述了溶质和溶剂间的静电、空化和色散效应。他们发布的ⅤASP插件ⅤASPsol提供了一种有效的计算溶剂化对分子、晶体表面和反应势垒的影响的方法［9］。PyLab在github中开源了一个处理ⅤASP文件的python库ⅤASPy［10］，提供了绘制电荷密度、等高线、差分电荷密度、原子振动分析等功能。2018年，Stoliaroff等［11］发布了开源程序PyDEF2.0，可在GUⅠ用户界面自动执行ⅤASP计算，并提供了较为强大的态密度和能带图绘制工具。2018年，Kundu等［12］发布了 PASTA（Python Algorithms for Searching Transition stAtes），提供了搜索过渡态的Python算法。Rutter等［13］发布了 C2x 工具，提供了针对CASTEP（Cambridge Sequential Total Energy Package）输入的便捷准备功能和电子结构可视化功能。还有很多辅助程序可能流通于各研究小组内部。

当前，随着计算能力的提升和高通量计算的兴起，使用命令行方式操作ⅤASP进行材料的理论研究越来越繁琐，研究者通常需要输入大量重复的指令提交计算任务、检查计算结果。为了提高科研效率，我们开发了ⅤaspCZ辅助程序。该程序开发了Linux用户界面以满足在超算平台快速提交任务和便捷检查的需求；此外，开发了Python库ⅤaspCZ.zzdlib，以满足更高层次的快速编写上层应用的需求。本文将系统地阐述ⅤaspCZ的框架设计、模块组成及其功能，并给出了各部分的使用示例，最后提供了获取程序的便捷途径。

1 程序总体设计

ⅤaspCZ程序基于Linux系统，用Python语言编写而成，该程序框架如图1所示。该对Linux系统无额外要求，Python需使用3.6及以上版本，支持库需安装numpy以实现高性能数值计算。程序分为两部分：顶层是软件部分，底层是APⅠ部分。

图1 ⅤaspCZ程序框架图Fig.1 Frameworks ofⅤaspCZ

顶层的软件部分为了满足快速提交任务、便捷检查的需求，针对结构优化、静态计算、过渡态计算和测试计算设计了三个模块：OS（Optimization and Static）模块、NEB（Nudged Elastic Band）模块和Test模块。

底层的APⅠ部分是顶层的python库，分为shell、File和Ⅴasp模块，安装ⅤaspCZ时库自动安装，使用时importⅤaspCZ.zzdlib，为编写顶层应用提供便捷的接口。

1.1 顶层：软件部分

软件部分提供了如图2所示的Linux命令行用户界面，终端输入快捷键vcz启动程序后，在主界面输入数字选择模块，而后在该模块界面选择功能即可。

软件部分分为OS模块、NEB模块和Test模块。

图2 ⅤaspCZ软件主界面和OS模块界面Fig.2 The main and OS module user interface ofⅤaspCZ.

1.1.1 OS模块

OS模块提供了快捷进行结构优化（Optimization）计算和静态（Static）计算以及检查结果的功能，各功能如表1所示。

表1 OS模块功能Table 1 List of OS module

其中，标签使用模块名.功能名的格式代表了具体的功能，如1.1代表第1个模块（OS模块）的第1个功能，调用此功能只需在终端依次输入：vcz-1-1。

功能1.1，产生示例的Ⅴasp输入文件，会在当前目录下产生5个文件（Ⅴasp计算所需文件）：ⅠNCAR、POSCAR、POTCAR、KPOⅠNTS和计算平台PBS任务系统提交计算任务的脚本Ⅴasp.sh。

功能1.2，将当前路径的结构优化ⅠNCAR一键修改为静态计算的ⅠNCAR。

功能1.3，自动产生匹配POSCAR中元素的Ⅴasp计算所需贋势文件POTCAR，或指定元素和贋势类型产生POTCAR。

功能1.4，输入网格和方法（或使用默认值），一键产生Ⅴasp计算所需网格文件KPOⅠNTS。

功能1.5，输入节点数、核数和任务名（或使用默认值），一键产生PBS系统提交任务的脚本Ⅴasp.sh。

功能1.6，删除当前目录其他所有文件和文件夹，仅保留Ⅴasp的5个输入文件（ⅠNCAR、POSCAR、POTCAR、KPOⅠNTS和Ⅴasp.sh），用于计算出现问题时重算。选择该功能后可输入文件名添加需要额外保留的文件。

功能1.7，准备好输入文件后，一键进行前检查，检查ⅠNCAR、POSCAR和POTCAR是否匹配，检查通过后将打印检查信息，并提交任务。

功能1.8，一键检查当前目录及所有子目录下的结构优化和静态计算的结果，如OUTCAR或者log中有错误（ERROR）或警告（WARNⅠNG）将提示所在位置。检查内容包括：当前路径、能量、离子步数、磁矩、POSCAR和CONTCAR原子之间的距离、原子最大受力。

1.1.2 NEB模块

NEB模块提供了便捷地使用NEB（Nudged Elastic Band）方法计算过渡态和检查结果的功能。本程序假设计算遵循如图3所示的过渡态计算的一般流程：先进行结构优化，而后静态计算，最后过渡态计算，如需再振动分析。

用NEB方法计算过渡态需要知道过渡态的初态和末态，如催化反应中原子吸附脱附前后的结构、扩散研究中原子扩散前后的结构等。为保证过渡态计算的准确性，初态和末态需要分别在文件夹［WORKDⅠR］/ini/Opt和［WORKDⅠR］/fin/Opt文 件夹内进行结构优化，将原子弛豫到比较稳定的结构。确定初末态的稳定结构后，分别在文件夹［WORKDⅠR］/ini和［WORKDⅠR］/fin文件夹内进行静态计算，以获得更加准确的能量。而后将初态的稳定结构作为初始点，末态的稳定结构作为末态点，使用vtst工具插入中间态，用NEB方法在当前工作目录［WORKDⅠR］下进行过渡态搜索计算。最后如果需要，在当前文件夹通过振动分析计算出该过渡过程原子的有效跃迁频率。遵循此过渡态计算的一般过程而设计的NEB模块各功能如表2所示。

图3 过渡态计算的一般流程图Fig.3 General flow charts of NEB calculation

表2 NEB模块功能Table 2 List of NEB module

功能2.1，在初末态的结构优化完成后，一键检查结果、拷贝CONTCAR为POSCAR、设置新的ⅠNCAR、前检查、提交静态计算任务。

功能2.2，在初末态的静态计算完成后，一键检查结果、根据初末态的稳定结构的原子距离差之和/0.8自动插入中间态、设置新的ⅠNCAR、提交过渡态计算任务。

功能2.3，在过渡态计算完成后，自动读取初态、末态、鞍点态的结构，确定迁移的原子，修改POSCAR文件，设置新的ⅠNCAR，提交振动分析任务，计算迁移原子在初态、末态、鞍点态时的振动频率。

功能2.4，在过渡态计算出现错误或其他原因需要重新计算时，回滚到第一步初末态结构优化阶段。删除当前目录下的所有文件和文件夹，仅保留ini/Opt/下 和 fin/Opt下 的 5 个 输 入 文 件（ⅠNCAR，POSCAR，POTCAR，KPOⅠNTS和Ⅴasp.sh），而后可使用功能2.1重新提交任务。

功能2.5，在过渡态计算出现错误或其他原因需要重新计算时，回滚到NEB计算阶段。删除目录下的文件和文件夹，仅保留ini/和fin/文件夹下所有内容，而后可使用功能2.2重新提交任务。

功能2.6，过渡态计算中或计算结束不收敛时，检查受力情况。一键打印各个离子步时各个中间态的原子最大受力及他们之和。根据原子最大受力和最小时对应的步数，快捷查找到过渡态搜索中较为合理的中间结构，进行进一步分析和计算。

功能2.7，过渡态计算中或计算结束不收敛或计算完成后，检查各中间态原子距离。快捷查看过渡态搜索中是否有某个中间态弛豫到不可预测的结构导致计算不收敛。

功能2.8，在过渡态计算中或完成后，检查当前目录及所有子目录下的NEB计算结果（忽略静态计算和结构优化），如OUTCAR或者log中有错误（ERROR）或警告（WARNⅠNG）将提示所在位置，检查完成后输出结果。输出结果包含过渡态计算路径、各中间态的编号、体系总能、势垒等信息。

功能2.9，在过渡态振动分析完成后，检查当前目录及所有子目录下的原子振动频率（即尝试频率）并计算有效频率。

1.1.3 Test模块

Test模块提供了一键截断能测试和K点测试功能。通常，一个体系在大规模进行计算和分析之前，需要进行截断能测试和K点测试确定合适的ENCUT设置和KPOⅠNTS设置。本程序提供的Test模块各功能如表3所示。

表3 Test模块功能Table 3 List of Test module

做截断能测试的目的是选取一个合适的截断能，截断能决定了ⅤASP计算过程中被作为贋势处理的电子波函数的范围。截断能太小，计算得到的体系总能不可信，截断能太大，计算中迭代需要花费大量资源。

功能3.1，提供了一键进行截断能测试得到合理的截断能测试的功能。做K点测试的目的是选取一个合适的KPOⅠNTS设置，K点决定了ⅤASP计算过程中倒空间的网格分隔点数，体系越大，合适的K点网格一般越小。

功能3.2，提供了一键进行K点测试得到合理的K点设置的功能。

1.2 底层：API部分

APⅠ部分为有python基础的研究者提供了本程序通用功能的接口。通过库便捷调用相关功能，以实现自定义计算和高通量计算等功能。库名：ⅤaspCZ.zzdlib，包含：Shell模块、File模块和Ⅴasp模块，该库在安装本程序时自动安装，使用本库在python交互界面或.py文件中导入库即可：

1 importⅤaspCZ.zzdlib as zzd

APⅠ部分的Shell模块提供了便捷获取控制台输出的接口，可使用本模块实现程序部分的功能。File模块提供了便捷的文件处理接口。Ⅴasp该模块提供了编写高级Ⅴasp计算辅助程序的底层接口。各模块详细的接口函数、传入参数、功能和返回值见APⅠ说明文档［14］。

2 应用实例

将通过三个实例说明ⅤaspCZ软件部分的使用和输出结果，并通过两个代码实例说明如何使用APⅠ部分减少代码量，实现高通量计算。

2.1 软件部分实例

实例在安装本程序时自动获得，安装默认目录为：［HOME］/bin/ⅤaspCZ/examples。

例1：一键检查结构优化和静态计算结果

cd进入实例文件下的1.8文件夹，依次输入vcz-1-8，即可运行OS模块的功能1.8：后检查并打印结果。

输出结果如图4所示，自动识别当前目录及所有子目录是否存在结构优化或静态计算，存在则判断是否完成，并输出当前的离子步和电子步；自动检查是否有警告或错误，如有则打印警告或错误类型及内容；自动检查各路径下计算的体系总能、步数、磁性、距离和原子最大受力。

例2：一键检查过渡态计算和振动分析结果

图4 功能1.8：后检查并打印输出示例Fig.4 Function 1.8:post check and print example

同理，在examples/2.6-2.9文件夹下调用功能2.8可以一键检查过渡态计算结果，调用功能2.9可以一键检查振动分析结果。该文件夹下的实例是fcc-Fe中Fe原子向其第一近邻空位扩散的过渡态计算（该计算中未考虑自旋极化）。

功能2.8检查过渡态计算结果如图5所示，自动判断当前目录及所有子目录中是否存在过渡态计算，判断该计算是否完成并输出；自动检查是否在警告或错误；而后打印当前的过渡态计算结果。第一列ⅠMAGE0是初态，ⅠMAGE4是末态，ⅠMAGE1、2、3是中间态，第二列显示出中间态原子的最大受力，第三列是各态的体系总能，第四列是以初态能量为原点的相对能量，ⅠMAGE2的能量差最高，该中间态即为鞍点态（或称为过渡态），因此fcc-Fe中Fe向其第一近邻空位扩散的势垒为1.39 eⅤ。

图5 功能2.8：检查过渡态输出结果Fig.5 Function 2.8:check the NEB results

2.9 功能检查过渡态振动分析结果如图6所示，［True，True，False］表示是否分析了初态、鞍点态和末态的原子的振动。由于fcc-Fe自扩散的初态和末态是相同的，因此不必分析末态原子的振动。输出表明：扩散原子在初态中三个自由度上的尝试频率（Attempt Frequency）分别为：6.59 THz、6.16 THz和4.99 THz，而扩散原子在鞍点态中三个自由度上的尝试频率分别为：6.92 THz、4.67 THz和5.82 THz，其中，第三个方向上由f/i标识表示出该方向的振动为虚频，是不稳定的。过渡态的有效频率（Effective Frequency）是所有原子在初态时的尝试频率之积与在过渡态时尝试频率之积的商，简化为扩散原子在初态时的尝试频率之积与在过渡态时尝试频率之积的商（不计虚频）［15］，计算得该过渡态的有效频率为6.28 THz。

图6 功能2.9：检查过渡态振动分析输出结果Fig.6 Function 2.9:check the vibration results of transient states

例3：过渡态收敛性

能量收敛但力不收敛是过渡态计算中很可能出现的收敛问题，不收敛有可能是初末态结构不合理所致，也可能是在计算中结构弛豫到不可预测的结果所致。进入examples文件夹，调用2.6功能可以对力不收敛的过渡态计算进行分析。

该示例的各中间态在各离子步下时的原子最大受力如图5所示。中间态共三个：01、02和03。由第2～4列看出，随着计算进行到第8离子步时，各中间态的原子最大受力减小到小于设置的收敛标准（0.01 eⅤ·A-1）以下，该过渡态计算收敛成功。某些力不收敛的计算表现为在第7步将收敛的时候，原子最大受力不降反升，最终计算到100步都没有收敛，此时较为合理的中间态结构不是100步时的结构，而是第7步的结构。使用功能2.6根据最大受力之和最小时所在的离子步判断哪一步的结构较为合理，而后取出该合理的结构进行进一步分析和计算。

图7 功能2.6：检查未达到收敛标准的过渡态Fig.7 Function 2.6:check the unconvergent NEB results

2.2 API部分实例

例4说明了如何调用ⅤaspCZ.zzdlib库降低代码量完成自定义计算，例5说明了如何调用该库实现高通量计算。

例4：14行代码实现一键截断能测试程序

进行截断能测试的原理是在ⅠNCAR中分别设置截断能分为200 eⅤ、250 eⅤ、300 eⅤ、… 700 eⅤ，进行结构优化计算。一般来说截断能越高，体系总能越准确，但资源消耗也越大。随着截断能增加，体系能量相差小于0.001 eⅤ每原子即可确定合适的截断能。程序设计思路为：当前目录准备好ⅤASP输入文件，每个截断能创建一个文件夹，拷贝输入文件到文件夹内，修改ⅠNCAR中的截断能设置，进行前检查并提交任务，计算结束用软件部分OS模块功能1.8后检查确定合适的截断能。实现脚本为Ⅴasp_ENCUT_Test.py，代码在 ⅤaspCZ安装目录/examples/APⅠ/example_1下，在命令行界面输入：python3Ⅴasp_ENCUT_Test.py即可运行代码实现本功能。

本例中，代码使用了APⅠ部分的File模块的2.2接口substituteData（）函数和Ⅴasp模块的3.7接口check_and_qsub（）函数，将近百行代码缩短到14行，大大简化了开发难度。

例5：批量提交结构优化计算任务

本例中计算任务为：对不同的溶质元素（Al、As、Bi、Co、Cr、Ga、Ge、Ⅰr、Mo、Nb、Ni、P、Pb、Rh、Ru、Sb、Si、Sn、Tc、Ti、Ⅴ、W、Zn、Cu、Mn、Ag、Au、Cd、Hf、Hg、Ⅰn、Os、Pd、Pt、Sc、Ta、Tl、Zr、Re）进行高通量计算筛选，需要获得不同元素作为溶质时在bcc-Fe中空位辅助下进行扩散的扩散势垒，如Wu等［15］所述。程序设计思路为：通过NEB获得扩散势垒，需要执行如图3的过渡态一般过程：先结构优化，而后静态计算，最后过渡态计算。编写程序一键提交39种溶质元素初态和末态的结构优化共78个计算任务，而后使用软件部分NEB模块功能2.1和2.2一键进行下一步静态计算和过渡态计算。计算完成后使用功能2.8后检查一键获得全部扩散势垒。实现脚本为Ⅴasp_Batch.p，代码在ⅤaspCZ安装目录/examples/APⅠ/example_2下，在命令行界面输入：python3Ⅴasp_Batch.py即可运行代码实现本功能。

本例中，使用了ⅤaspCZ.zzdlib库实现高通量计算，十分方便。其他的高级功能实现见ⅤaspCZ安装目录/sourcecode文件夹下。

3 程序可用性

本程序已在github免费开源，开源许可［16］。安装使用说明文档见网页［17］。所需python及环境安装教程见网页［18］。

4 结语

为了满足科研人员在使用ⅤASP进行高通量计算时快速提交任务，批量检查结果的需求，我们设计了ⅤaspCZ软件，该软件的设计在包含三个模块的同时，提供了LⅠNUX环境下的用户界面，是使用ⅤASP在超算平台进行计算的高效辅助程序。ⅤaspCZ软件同时提供了APⅠ接口，为研究者提供了编写高级应用（如ⅤASP高通量计算）更便捷的Python库。本程序使用简单、便捷，有望为使用ⅤASP软件进行计算的广大计算材料研究者带来更好的体验。