VASP 入门：从四个文件到可信结果

本章把“一个硅晶胞”从纸上的结构变成可追溯的计算任务。目标不是记住一大串参数，而是知道每一份输入在决定什么、每一个输出能回答什么，以及何时应该停下来检查。建议先读完量子化学基础中的周期性、平面波和自洽场，再按本章顺序操作。

VASP 是需要许可的软件，网页和输入模板不会提供程序或 POTCAR。先确认你属于获得许可的用户组，再向课题组或计算中心取得已安装模块、势文件位置、队列名称与作业提交方法。学校的 Linux 集群通常已经编译好 VASP；不建议刚入门就同时处理编译、MPI 和材料模型三个问题。

在登录节点上可以编辑文件、查看已有输出和做轻量预检查；正式计算应提交到计算节点。CPU 核数不是越多越好：两个原子的练习只适合小规模资源测试。不同编译版本、MPI、CPU/GPU 后端的并行设置并不通用，先采用计算中心已经验证的启动方式。记录准确版本号和编译环境，遇到错误时这些信息比“我用了 VASP”更有用。

开始前的自查：能否在文本编辑器中看到普通文本文件？是否知道自己当前所在目录？是否能查看作业状态？势文件是 PBE 的哪个发布版？如果这些问题仍不清楚，先完成本章“集群操作”再提交任务。

核对依据：官方安装说明 · VASP 许可

把每项计算放入独立目录，会比在同一个文件夹中反复覆盖文件更容易定位问题。一次常见计算至少需要以下四个文件，名称区分大小写。它们互相约束：结构写了什么元素，势文件就必须按同样顺序提供；选择什么泛函，也必须核查所用势与计算设定是否相容。

POSCAR 的缩放因子乘在三个晶格矢量上；使用 Direct 时，三个数是相对于晶格矢量的分数坐标，不是 Å。若使用 Cartesian，坐标也受缩放因子影响。元素行是便于人阅读的标签，实际离子种类还取决于 POTCAR 和原子数分组，所以只改元素名称而不换势文件，并没有真正换掉元素。

POTCAR 从授权势库取得，不手工改写其中参数。含 A、B 两种元素的结构按 A B 排列时，POTCAR 也应先 A 后 B；而 MAGMOM 通常按每个原子排列，LDAUU 等则按每种元素排列。这种“按原子”与“按种类”的区别会在磁性章节再次出现。

INCAR 建议每行一个参数，统一使用英文半角符号。注释用于说明选择理由，不要把中文说明直接接在数值后而没有注释符。不要在同一文件中留下多个互相矛盾的 ENCUT 或 ISPIN；真正被程序读取的值，应在 OUTCAR 中核对。

核对依据：POSCAR · POTCAR · INCAR · OUTCAR

问题：把离子位置固定，求这个周期性硅结构在 PBE 近似下的基态电子密度。这里的“静态”表示原子不移动，不是关闭电子迭代。先从一个非磁性、两个原子的体系熟悉操作，之后再处理过渡金属和表面。

下列 POSCAR 使用金刚石结构的两个原子原胞，晶格参数 5.43 Å 是学习用初始值，不是本教程得到的优化结果。三个原胞矢量并非互相垂直，体积是常规立方晶胞的四分之一。不要把后面的原胞 k 路径直接用到含八个原子的常规晶胞上。

Diamond Si primitive cell; a=5.43 A
5.4300000000
0.0000000000 0.5000000000 0.5000000000
0.5000000000 0.0000000000 0.5000000000
0.5000000000 0.5000000000 0.0000000000
Si
2
Direct
0.0000000000 0.0000000000 0.0000000000
0.2500000000 0.2500000000 0.2500000000

建立 si-scf 目录，保存下面的 KPOINTS。第三行的 Gamma 表示网格以 Γ 为中心；第四行 8 8 8 表示三个倒格矢方向的划分数。这是自动生成的均匀网格，第二行 0 表示由网格规则生成，而不是“计算零个 k 点”。对称性会将网格约化为不可约 k 点。

Si uniform Gamma-centered mesh
0
Gamma
8 8 8
0 0 0

SYSTEM = Si PBE SCF learning example
GGA = PE
ENCUT = 520
PREC = Accurate
LREAL = .FALSE.
ISPIN = 1
ISTART = 0
ICHARG = 2
ALGO = Normal
EDIFF = 1E-6
NELM = 120
ISMEAR = 0
SIGMA = 0.05
IBRION = -1
NSW = 0
LCHARG = .TRUE.
LWAVE = .TRUE.

POTCAR：只需一个 Si PAW-PBE 数据集，具体版本由你的授权环境决定，下载包不含这个文件。运行结束后，用 OUTCAR 核对离子数为 2、元素与价电子设置正确、最终电子循环达到 EDIFF。不能用“跟别人的总能量小数完全一样”判断成功，因为不同势版本与数值设置会改变总能量基准。

小练习：用可视化软件打开 POSCAR，解释为什么第二个 Si 的分数坐标是 (¼,¼,¼)，再根据晶格矢量把它转换成笛卡尔坐标。读出 OUTCAR 的晶胞体积，与你手算的 a³/4 比较。这比先盯着漂亮的能带图更能排除结构错误。

核对依据：ENCUT · KPOINTS · EDIFF · LREAL · PREC · ALGO

电子循环达到 EDIFF，只说明在当前基组、当前网格和当前近似下完成了求解。它没有告诉你平面波是否够多、布里渊区是否采样充分。真正的收敛测试要比较目标量随着数值设置加密后是否稳定。可以把它理解为给显微镜调分辨率：画面不再明显变化，才能讨论画面中的细节。

1. 先固定几何和所有物理设定。本练习始终使用同一个 si-scf/POSCAR、同一个 POTCAR，不在每个截断能下分别优化结构。
2. 测试 ENCUT。固定一个暂时较密的均匀网格，只改变截断能。可从 max(ENMAX) 以上选连续增大的 5—6 个值；400、450、500、550、600 eV 只在它们全部不低于你的 max(ENMAX) 时才可作为示例序列。
3. 测试 k 网格。使用第一轮已趋于稳定的截断能，比较 4³、6³、8³、10³、12³ 等同一家族 Γ 中心网格；若 12³ 仍不稳定就继续加密。
4. 做交叉检查。在最后选定的 k 网格上再提高一次 ENCUT，确认第一轮结论没有被粗网格掩盖。

每组都从头生成电荷密度并使用独立目录。下面的下载脚本只建立输入，不运行 VASP，也不自动提交几十个作业。它要求使用本页“一行一个参数”的静态 INCAR，会忽略 #/! 后的注释、拒绝重复参数和分号多参数行；NSW>0 或 IBRION 不是 −1 时会停止，防止误把结构优化源作为固定几何测试。k 网格测试也会检查源 ENCUT 不低于势文件最大 ENMAX。目标目录若已存在同样停止，避免覆盖先前结果。脚本只在你自己的计算环境内复制授权 POTCAR，不应把生成目录整体上传到公开仓库。

# 在示例目录的共同上一级运行；先检查实际 ENMAX
python3 prepare_convergence.py encut si-scf conv-encut 400 450 500 550 600
# 先把 si-scf/INCAR 中的 ENCUT 更新为已测试的选定值，再做这一组
python3 prepare_convergence.py kmesh si-scf conv-kmesh 4 6 8 10 12

#!/usr/bin/env python3
"""Prepare independent input directories only. Does not submit jobs or run VASP."""
import argparse
import csv
import math
from pathlib import Path
import re
import shutil

p = argparse.ArgumentParser()
p.add_argument("mode", choices=["encut", "kmesh"])
p.add_argument("source", type=Path)
p.add_argument("destination", type=Path)
p.add_argument("values", nargs="+", type=float)
a = p.parse_args()
source = a.source.resolve()
destination = a.destination.resolve()
for name in ["INCAR", "POSCAR", "KPOINTS", "POTCAR"]:
    if not (source / name).is_file():
        raise SystemExit(f"Missing source input: {name}")
if destination.exists():
    raise SystemExit("Destination exists; choose a NEW directory to avoid overwriting")
incar = (source / "INCAR").read_text(encoding="utf-8")
potcar = (source / "POTCAR").read_text(encoding="utf-8", errors="replace")
limits = re.findall(r"ENMAX\s*=\s*([0-9.]+)", potcar)
if not limits:
    raise SystemExit("ENMAX not found; inspect POTCAR before proceeding")
max_enmax = max(map(float, limits))
# Restrict this educational generator to an unambiguous static source INCAR.
settings = {}
for number, raw in enumerate(incar.splitlines(), start=1):
    clean = re.split(r"[!#]", raw, maxsplit=1)[0].strip()
    if not clean:
        continue
    if ";" in clean or clean.count("=") != 1:
        raise SystemExit(f"Use one unambiguous tag per line (INCAR line {number})")
    key, value = (part.strip() for part in clean.split("=", 1))
    key = key.upper()
    if key in settings:
        raise SystemExit(f"Duplicate INCAR tag: {key}")
    settings[key] = value
try:
    nsw = int(settings.get("NSW", "0"))
    ibrion = int(settings.get("IBRION", "-1"))
except ValueError:
    raise SystemExit("NSW and IBRION must be plain integers")
if nsw > 0 or ibrion != -1:
    raise SystemExit("Source must be STATIC: NSW <= 0 and IBRION = -1; do not use a relaxation INCAR")
if a.mode == "kmesh":
    try:
        cutoff = float(settings["ENCUT"].replace("D", "E").replace("d", "e"))
    except (KeyError, ValueError):
        raise SystemExit("kmesh source must set an explicit numeric ENCUT")
    if not math.isfinite(cutoff) or cutoff < max_enmax:
        raise SystemExit(f"Source ENCUT must be >= max ENMAX ({max_enmax} eV)")
values = a.values
if any(not math.isfinite(x) or x <= 0 for x in values):
    raise SystemExit("Values must be positive finite numbers")
if len(set(values)) != len(values):
    raise SystemExit("Duplicate values")
if a.mode == "encut" and min(values) < max_enmax:
    raise SystemExit(f"All ENCUT values must be >= max ENMAX ({max_enmax} eV)")
if a.mode == "kmesh" and any(x != int(x) for x in values):
    raise SystemExit("k mesh values must be positive integers")
# This generator expects the distributed one-tag-per-line INCAR, not an arbitrary INCAR.
base = []
for line in incar.splitlines():
    key = line.split("=", 1)[0].strip().upper()
    if key not in {"ISTART", "ICHARG", "LWAVE", "LCHARG"} and not (a.mode == "encut" and key == "ENCUT"):
        base.append(line)
base += ["ISTART = 0", "ICHARG = 2", "LWAVE = .FALSE.", "LCHARG = .FALSE."]
destination.mkdir(parents=True)
rows = []
for value in values:
    label = f"{value:g}"
    run = destination / (("e" if a.mode == "encut" else "k") + label)
    run.mkdir()
    for name in ["POSCAR", "POTCAR"]:
        shutil.copy2(source / name, run / name)
    text = "\n".join(base) + "\n"
    if a.mode == "encut":
        text += f"ENCUT = {value:g}\n"
        shutil.copy2(source / "KPOINTS", run / "KPOINTS")
    else:
        n = int(value)
        (run / "KPOINTS").write_text(f"Gamma mesh convergence\n0\nGamma\n{n} {n} {n}\n0 0 0\n", encoding="utf-8")
    (run / "INCAR").write_text(text, encoding="utf-8")
    rows.append([run.name, a.mode, label, "prepared_not_run"])
with (destination / "manifest.csv").open("w", newline="", encoding="utf-8") as f:
    writer = csv.writer(f)
    writer.writerow(["directory", "variable", "value", "status"])
    writer.writerows(rows)
print(f"Prepared {len(rows)} independent folders; submit manually after review")

用统一的能量定义比较，比如同一种展宽下的 energy(sigma→0)，并同时记录自由能和熵项；不要把一组 TOTEN 与另一组 E0 混在一张表里。对这个两原子原胞，单胞能量差除以 2 才得到 eV/atom。针对相稳定性、吸附能或势垒，最终还要检验相应能量差本身，而不只是孤立的总能量。

课堂上可暂取“与更精细两组相比小于 1 meV/atom”作为能量收敛练习目标，但这不是对所有课题的统一标准。若你要讨论 0.5 meV/atom 的磁各向异性，这样的门槛显然不够。至少要求最后几个点形成稳定趋势，不能恰好两点接近就宣称完成。

核对依据：截断能与可比性 · k 点网格规则 · 体积与数值收敛

SCF 在固定结构上调整电子；结构优化在多个 SCF 之间移动原子，必要时也改变晶胞。能量对原子位移的负导数是力，能量对应变的响应决定应力。优化想要的是：电子先足够收敛，让力可信，再用这些力更新结构。过松的电子收敛会把数值噪声传给离子优化器。

在全新 si-relax 目录中放入上节 Si 的 POSCAR、相同 POTCAR 与已收敛的均匀 KPOINTS，再使用下面 INCAR。520 eV 仅在满足你的体积优化精度要求时才能保留。先核查是否至少达到约 1.3×max(ENMAX) 这一常见的提高截断能起点，再通过残余应力与更高截断能测试决定实际数值。

SYSTEM = Si PBE cell relaxation learning example
GGA = PE
ENCUT = 520
PREC = Accurate
LREAL = .FALSE.
ISPIN = 1
ISTART = 0
ICHARG = 2
ALGO = Normal
EDIFF = 1E-7
NELM = 120
ISMEAR = 0
SIGMA = 0.05
IBRION = 2
ISIF = 3
NSW = 100
POTIM = 0.5
EDIFFG = -0.01
LCHARG = .FALSE.
LWAVE = .FALSE.

理想硅结构中两个原子的分数坐标可能从一开始就具有零净力，但初始晶胞体积仍未必最佳。因此本例不能只看原子坐标有没有变，还应比较晶格、体积和最终应力。对于固定部分原子的表面，优化条件只适用于允许运动的自由度，固定原子保留的力仍值得检查和报告。

Pulay 应力为什么出现？有限平面波基组会随晶胞改变，其不完备误差可能看起来像额外应力。因此“外压接近零”在很低截断能下未必意味着可靠的平衡体积。使用较高截断能、足够精细网格，并从优化后的结构重新开始一轮，是检验晶格稳定性的实用方法；不要随意加一个 PSTRESS 数字把压力“调成零”。

完成条件：看到离子优化停止条件确实满足，最后几个电子循环也都收敛；检查最终受优化自由度上的力、晶格自由度相关应力和结构合理性。查看 CONTCAR 是否出现不合理短键或意外坍塌。最后还需更精细的静态计算复核。几何优化找到的是一个局部驻点，是否为稳定极小值、是否存在更低能结构，还涉及声子和结构搜索。

核对依据：ISIF 自由度 · 离子优化算法 · EDIFFG 的符号 · Pulay 应力与体积

结构优化和最终性质计算通常分开做。前者的重点是找到结构，后者在不移动离子的条件下，用统一、较严格的数值设置生成可比能量、电荷与电子态。特别是在最后一次离子更新后，不应未经检查就把任何旧 CHGCAR 当成对应最终结构的已收敛密度。

在新的 si-static 目录中，把已经核验收敛的 si-relax/CONTCAR 复制为 POSCAR，使用相同 POTCAR。下面 12³ 网格和 520 eV 同样是示例，实际选定值必须来自你的收敛实验。为了减少初学者的重启混淆，这里选择从头做电子自洽；它可能比恰当复用波函数慢，但路径清楚。

SYSTEM = Si final static learning example
GGA = PE
ENCUT = 520
PREC = Accurate
LREAL = .FALSE.
ISPIN = 1
ISTART = 0
ICHARG = 2
ALGO = Normal
EDIFF = 1E-7
NELM = 120
ISMEAR = 0
SIGMA = 0.05
IBRION = -1
NSW = 0
LCHARG = .TRUE.
LWAVE = .TRUE.

Si final uniform Gamma-centered mesh
0
Gamma
12 12 12
0 0 0

这份 INCAR 的 IBRION=-1、NSW=0 让结构固定；与第一项任务相比，电子阈值更紧，网格更密。它保存的 CHGCAR 可以作为下一章 PBE 能带计算的输入。必须同时记录它对应的 POSCAR，而不是只孤立保存 CHGCAR。

电子 EDIFF 比较的是相邻迭代的能量变化，不是“DFT 与实验误差小于 10⁻⁷ eV”。即使数值收敛得极好，泛函近似、温度、缺陷浓度和模型边界条件仍可能带来更大的偏差。把数值误差与模型误差分别写清楚，是学习计算化学的重要一步。

总能量的绝对值通常很大且依赖势的参考基准，研究更常比较受控条件下的能量差。不同元素组成的结构不能直接按总能量大小排名稳定性，需要形成能、化学势或相图等合适定义。不同原子数比较时，也应先明确是每原子、每化学式还是整个超胞。

交付给下一步的数据：确认收敛的 POSCAR、四输入的哈希与参数记录、CHGCAR、必要的 WAVECAR、OUTCAR、OSZICAR、结构化输出及作业日志。不要只留下一个截图或手抄能量；后续若发现某个设置有误，完整记录能让你判断哪些结果需要重做。

核对依据：OSZICAR 字段 · CHGCAR 的内容与重启 · 电子收敛条件

以下命令供 Linux 集群终端使用，不是 Windows PowerShell 命令。先用 pwd 确认位置，再看目录内容。所有正式计算使用新目录；复制前确认目标没有需要保留的文件。不要在不清楚路径时执行递归删除，也不要把“重新跑一遍”理解为直接抹掉上次失败证据。

pwd
ls -lh
head -n 12 POSCAR
cat INCAR
cat KPOINTS
grep -E 'TITEL|ZVAL|ENMAX' POTCAR
sha256sum INCAR POSCAR KPOINTS POTCAR
# 以下从共同上一级目录执行；仅从零建立新目录时执行 mkdir
# 配套包若已有 si-static，先确认未混入旧结果，并略过下面这一行
mkdir si-static
# 核验 si-relax 确实收敛后，再复制；-i 避免静默覆盖
cp -i si-relax/CONTCAR si-static/POSCAR
cp -i si-relax/POTCAR si-static/POTCAR

下面为 Slurm 待配置模板，文件名以 .template 结尾。必须把 YOUR_PARTITION 与 YOUR_VASP_MODULE 替换为计算中心给出的真实值，核对核数、时限、内存要求和启动器后，再另存为 submit.slurm。有些机器采用 mpirun 或特定 srun 插件，应使用管理员给出的方式，不能把示例当作集群通用配置。

配套包中的模板位于共同父目录。完成配置后，将它保存为本次计算目录中的 submit.slurm，例如 si-scf/submit.slurm。提交前先进入保存四个输入文件的实际目录；在父目录执行 sbatch si-scf/submit.slurm 并不会自动把提交目录变成 si-scf。后续任务也分别保存自己的脚本，并从对应目录提交。

#!/bin/bash
#SBATCH --job-name=si-scf
#SBATCH --partition=YOUR_PARTITION
#SBATCH --nodes=1
#SBATCH --ntasks=8
#SBATCH --cpus-per-task=1
#SBATCH --time=00:30:00
#SBATCH --output=slurm-%j.out
#SBATCH --error=slurm-%j.err

# Save the configured script as submit.slurm INSIDE the calculation directory.
# Submit after cd into that directory, where all four input files are present.
# For SOC/noncollinear tasks, replace BOTH vasp_std commands below with vasp_ncl.
set -euo pipefail
: "${SLURM_SUBMIT_DIR:?Submit through sbatch after editing placeholders}"
cd "$SLURM_SUBMIT_DIR"
module load YOUR_VASP_MODULE
export OMP_NUM_THREADS=1
for f in INCAR POSCAR KPOINTS POTCAR; do
  test -s "$f" || { printf 'Missing or empty input: %s\n' "$f" >&2; exit 1; }
done
command -v vasp_std >/dev/null
sha256sum INCAR POSCAR KPOINTS POTCAR > inputs.sha256
srun vasp_std > vasp.out 2>&1

#SBATCH 行是 Slurm 请求资源的指令，不是 VASP 参数。这里的 8 个 MPI 进程和 1 个线程只演示 CPU 作业；它不是性能推荐，更不适用于所有 GPU 环境。set -euo pipefail 让常见脚本错误及失败命令更早暴露，输入文件检查防止在错误目录误启动。输入哈希记录在作业目录里，不包含势文件正文。

# 已完成配置，并保存为 si-scf/submit.slurm；以下从共同父目录开始
# 如果当前已经在 si-scf 中，略过 cd 这一行
cd si-scf
sbatch submit.slurm
squeue -u "$USER"
# 替换为 sbatch 返回的数字作业号
sacct -j YOUR_JOB_ID --format=JobID,State,Elapsed,ExitCode
# 在实际计算目录观察文本，不要边跑边改 INCAR
tail -n 30 OSZICAR
tail -n 40 vasp.out

等待排队时，作业没有产生 OUTCAR 是正常的；运行后突然停止，则检查超时、内存限制和节点错误。取消作业时只取消你确认的具体作业号，避免批量操作影响其他任务。输入文件在计算运行中不应修改；需要改参数时，新建任务目录，并记下修改原因。

提交后要建立三个对应关系：目录对应哪个结构；作业号对应哪次运行；结果对应哪套实际参数。建议在研究笔记中记录“今天想检验什么假设”，而不只是“今天跑了 12 个作业”。模板中的时限 30 分钟不是软件对该任务耗时的保证，实际耗时取决于硬件与编译环境。

核对依据：Slurm sbatch 官方文档 · VASP 环境与安装

最容易浪费算力的做法，是看到失败就把 NELM 和 NSW 都加到很大。先判断失败发生在哪一层：程序能否启动、结构是否合理、电子是否收敛、离子是否收敛、目标性质是否随数值设置稳定。只有找对层次，修改参数才有意义。

若能量持续震荡，可从头初始化排除不兼容 WAVECAR/CHGCAR，确认 ISMEAR 适合体系，再按照官方电子最小化说明逐项调整 ALGO 或混合参数。一次只改变一个有理由的选项，并保留对照。增加 NELM 只适用于确实在缓慢趋近收敛的情形；毫无下降趋势时延长运行往往只是重复失败。

原子数变大、表面真空变厚、体系从绝缘体变金属之后，原来有效的网格与混合参数可能需要重新检验。不要根据某个报错搜索到的一行参数，未经理解就复制进所有任务。遇到程序崩溃时，应保存报错上下文、软件版本、最小输入和运行环境，敏感势文件只按许可渠道共享。

科学问题也需要独立排查。例如“算出的晶格比实验大”可能涉及泛函、温度与零点振动；“带隙偏小”可能来自半局域泛函本身，不能靠加大 SIGMA 修正；“缺陷形成能负数”可能是参考化学势或计量关系出错，不能简单通过提高 ENCUT 解释。

# 计算记录：请逐项填写
科学问题与待检验假设：
结构来源 / 原胞或超胞 / 原子数 / 电荷状态：
VASP 版本 / 编译环境 / 可执行程序：
POTCAR 发布版 / TITEL / ZVAL / ENMAX / SHA-256：
泛函 / U 与 J / SOC / 自旋初态：
ENCUT / PREC / k 网格 / 展宽：
电子与离子停止阈值 / 固定自由度：
收敛测试目标量 / 阈值 / 对照目录：
作业号 / 日期 / 资源 / 退出状态：
电子是否收敛 / 最大力 / 最终应力：
使用的能量定义 / 单位 / 归一化方式：
后处理工具版本 / 脚本 / 图的数据来源：
已知局限与下一步验证：

本章完成练习：画出 ENCUT 与 k 网格收敛图；比较初始结构和优化结构；用最终静态 OUTCAR 写一段真实的“计算方法”；说明你的证据能支持什么、还不能支持什么。暂时没有集群权限时，可以先完成文件阅读、结构可视化和实验设计，但不要填写假想能量充当运行结果。

下一步进入电子结构：DOS、能带、磁性与电荷分析；需要按任务快速定位参数时，使用参数与排错索引。对于催化、声子、分子动力学等任务，请先确保已经理解本章的收敛和记录方法。

核对依据：ALGO 与电子迭代 · NELM 的含义 · 电子最小化方法

先读文件说明与对应章节，再补齐模型、POTCAR 与前置计算。后缀为 .template 的文件包含待改参数，不能原样提交。

正文为独立整理的学习说明。标签含义和版本差异以与你的程序版本对应的官方文档为准。

1. VASP 许可
2. 安装与可执行程序
3. POSCAR
4. POTCAR
5. INCAR
6. KPOINTS
7. ENCUT
8. PREC
9. LREAL
10. EDIFF
11. ALGO
12. NELM
13. ISIF
14. IBRION
15. EDIFFG
16. Volume relaxation
17. OUTCAR
18. OSZICAR
19. CHGCAR
20. Slurm sbatch
21. Electronic minimization