VASPBERRY入门与MoS₂波函数准备
从原子高度优化到完整SOC二维网格，准备实际可用于Berry曲率、Chern和Z₂分析的VASP输入。
阅读本篇目录
VASPBERRY读取VASP波函数进行Berry曲率与拓扑后处理。本系列使用鸢算科技集群实际计算的三原子1H-MoS₂，分为三个步骤：本篇准备结构与波函数，第二篇画曲率，第三篇计算Chern数和Z₂。
下载输入、数值表、短命令脚本和五张图的绘图程序。资料包不包含POTCAR和WAVECAR；新计算使用您有权使用的势文件。
1 先看模型和流程#
采用固定面内晶格常数约3.171634 Å、c方向15 Å的单层模型。先固定晶胞和Mo位置，只优化两个S的z坐标；随后所有SOC网格沿用同一优化结构。这里的原子高度优化不改变晶格常数，也不优化真空厚度。
真实结构优化进行了3次DFT力评估，自由原子的最大残余力约0.001458 eV/Å，低于0.02 eV/Å。SOC分别使用6×6×1、12×12×1和18×18×1完整网格，实际计算了36、144和324个k点。
如果先了解结果，把下载文件保存为 vaspberry.zip，运行：
python3 -m zipfile -e vaspberry.zip .
cd vaspberry-MoS2-wave281
bash run.sh
ls analysis/figures
这一步只读取已有数值并重画五张图。完整的新计算从下面的结构优化开始。
2 优化原子高度#
在资料包根目录执行：
bash prepare.sh
依次输入 1（结构优化）、1（P100）或 2（V100 16G），再按照提示执行：
cd work/relax
sbatch run.sh
任务使用VASP 6.6.1的常规 vasp_std，GPU版本由NVHPC 25.3编译。初始POSCAR完整如下，F F T 表示只允许z坐标变化：
1H-MoS2
1.00000000000000
3.1716343032172727 0.0000000000000000 0.0000000000000000
1.5858171516086363 2.7467158781003151 0.0000000000000001
0.0000000000000000 0.0000000000000000 15.0000000000000000
Mo S
1 2
Selective dynamics
Direct
0.0000000000000000 0.0000000000000000 0.5000000000000000 F F F
0.3333333333333333 0.3333333333333333 0.6064573703372080 F F T
0.3333333333333333 0.3333333333333333 0.3935426296627920 F F T
完整优化INCAR：
SYSTEM=MoS2-height-preparation
ENCUT=450
PREC=Accurate
EDIFF=1E-6
EDIFFG=-0.02
IBRION=2
ISIF=2
NSW=80
POTIM=0.5
ISPIN=1
ISYM=-1
NBANDS=24
KPAR=1
ISMEAR=0
SIGMA=0.02
ALGO=Normal
NELM=160
LREAL=.FALSE.
LASPH=.TRUE.
LWAVE=.FALSE.
LCHARG=.FALSE.
对应KPOINTS：
MoS2 atomic height preparation
0
Gamma
6 6 1
0 0 0
EDIFF=1E-6 控制电子自洽，负的 EDIFFG=-0.02 指定力阈值。ISIF=2 固定晶胞；POSCAR选择性约束进一步限制哪些原子和方向可动。
计算结束先检查：
grep 'reached required accuracy' OUTCAR
tail OSZICAR
ls -lh CONTCAR
确认优化收敛，保留CONTCAR。不能仅以Slurm任务结束判断结构优化成功。
3 准备完整SOC网格#
回到资料包根目录，再运行：
cd ../..
bash prepare.sh
选择 2（SOC）、网格 6、GPU类型。程序从您完成的 work/relax/CONTCAR 复制结构，再生成独立SOC目录：
cd work/soc-mesh6
sbatch run.sh
SOC必须使用 vasp_ncl，本例的实际任务命令为：
module load vasp/6.6.1-p100
mpirun --bind-to none -np 1 vasp_ncl > vasp.log 2>&1
V100 16G改为 vasp/6.6.1-v100，资源申请也同时改为对应队列和GPU类型；准备脚本会一起处理。
完整SOC INCAR：
SYSTEM=MoS2-full-SOC-mesh6
ENCUT=450
PREC=Accurate
EDIFF=1E-6
IBRION=-1
NSW=0
ISYM=-1
NBANDS=48
KPAR=1
LSORBIT=.TRUE.
SAXIS=0 0 1
MAGMOM=0 0 0 0 0 0 0 0 0
GGA_COMPAT=.FALSE.
ISMEAR=0
SIGMA=0.02
ALGO=Normal
NELM=160
LREAL=.FALSE.
LASPH=.TRUE.
LWAVE=.TRUE.
LCHARG=.FALSE.
LORBIT=11
ISYM=-1 关闭对称性约化，KPAR=1、1个MPI进程便于保持完整波函数记录；LWAVE=.TRUE. 保留后处理输入。MAGMOM 为零的非磁模型不施加外场；这也是后续讨论时间反演对称性的物理前提之一。
KPOINTS采用明确列出的完整二维网格，而不是高对称路径。资料包已经包含全部点；下面的等价生成命令展示6×6网格如何写出：
python3 - <<'KPOINTS'
n = 6
with open('KPOINTS', 'w') as f:
f.write(f'Full2DexplicitGammaSOCmesh\n{n*n}\nReciprocal\n')
for j in range(n):
for i in range(n):
f.write(f'{i/n:.16f} {j/n:.16f} 0 {1/(n*n):.16f}\n')
KPOINTS
要做12、18网格，回到根目录运行准备脚本，选择相应网格。所有网格用同一优化结构，输出分别保留，不能覆盖前一档WAVECAR。
4 确认选取了整个占据子空间#
grep -E 'NELECT|NBANDS|NKPTS' OUTCAR
ls -lh WAVECAR
本例使用PBE Mo、S势，Mo价电子数6、每个S为6，总电子数18。在SOC自旋器表示下，实际输出中前18条带占据为1、其后为0，因此选择 --bands 1:18。不能把普通非自旋简并计算的“电子数除以2”直接搬过来。
实际采样的最小直接能隙分别为1.6371、1.6388、1.6388 eV。换势、掺杂、磁态或电子数以后，必须重新查看占据和带隙，重新选择子空间。
接下来按Berry曲率实操处理WAVECAR，再看Chern与Z₂。本教程使用已安装的VASPBERRY 1.6.3，命令接口可对照对应官方版本。