LAMMPS自旋NEB操作实例

普通NEB优化原子位置，自旋NEB优化磁矩方向。本例原子始终固定在原点，只让一个单位磁矩在外场和单轴各向异性能量下翻转，实际比较 spin、spin/cg、spin/lbfgs 三种优化器与7个路径图像。

这是一个有解析能量基准的经典单自旋教学模型，不是某种铁、镍或其他材料的DFT翻转能垒。它便于先弄清文件格式、固定模长、路径、力矩和单位，再换成合适的真实磁性模型。

下载完整输入、21个实际图像、角度差分数据与图片，在个人目录执行：

mkdir spin-neb
cd spin-neb
curl -fL https://yuansuan.top/assets/downloads/full-function/lammps-spin-neb-wave82-2c81dd6c57bd.zip -o examples.zip
unzip examples.zip
sbatch submit.slurm

脚本申请7个CPU核、一个P100资源，每个图像1个MPI进程，顺序完成三种优化器和单点检查。实际GNEB使用CPU/MPI，没有GPU力计算或速度对比。模块为 lammps/30Sep2026-p100-gpu（GCC13.3、OpenMPI5.0.10、CUDA12.8）；分析用 ase/3.29.0 的NumPy/Matplotlib。

完成后：

module purge
module load ase/3.29.0
python analyze.py work
python export_spin_path.py work
python plot_results.py work
head -10 work/analysis/spin-neb.csv

已有实算在 reference/，可将后处理参数 work 换为 reference 先看结果。新运行换一个解压目录，不覆盖自己的work。

1 自由度和模型能量#

设置磁矩模长为1，单位方向 s=(sx,sy,sz)，保持 sx²+sy²+sz²=1。各向异性K=0.01 eV，易轴沿z；磁场10 T沿x。模型能量为 E=−b·sx−K·sz²，在当前30Sep2026实现中 b=2 μB H=0.001157802 eV。

两个端点的sx相同，约0.0578901，sz分别为±0.998322962。初始路径故意加入一个方位角弯曲，供GNEB实际优化；不是直接把已经精确的解析路径当作优化结果。

2 自旋端点文件与原子位置文件不同#

自旋文件第一行是记录数，后面每行是 ID 磁矩模长 x y z sx sy sz。本例第7图像文件完整如下：

1
1 1 0 0 0 0.05789010000000043 2.126847770019621e-18 -0.9983229619326554

普通原子位置NEB的末态文件只有 ID x y z，不能直接拿来替代这个格式。这里xyz均为0，因为晶格位置被冻结；变化的是后三个方向分量。

第一图像由输入中的 set atom ... spin 给定，each 模式为其他图像读取各自文件。单独修改spin-1.coords不能代替修改第一图像初态。

3 完整GNEB输入#

variable image uloop 7
units metal
atom_style spin
atom_modify map array
region box block -5 5 -5 5 -5 5 units box
create_box 1 box
create_atoms 1 single 0 0 0 units box
mass 1 1
set atom 1 spin 1 0.05789010000000001 0 0.9983229619326554
pair_style zero 3.0
pair_coeff * *
fix external all precession/spin zeeman 10 1 0 0 anisotropy .01 0 0 1
fix_modify external energy yes
fix band all neb/spin 1.0
compute magnetic all property/atom sp spx spy spz fmx fmy fmz
thermo_modify norm no format float %.16g
min_style spin
min_modify discrete_factor 20
neb/spin 0.0 1.0e-7 4000 4000 20 each spin-${image}.coords
write_dump all custom neb-${image}.dump id type x y z c_magnetic[*] modify sort id format float %.16g
unfix band
run 0
print "ENERGY ${image} $(pe:%.16g)" file energy-${image}.dat screen no
write_dump all custom physical-${image}.dump id type x y z c_magnetic[*] modify sort id format float %.16g

关键设置：

• atom_style spin 保存磁矩模长和方向，atom_modify map array 提供按ID映射。
• precession/spin 同时给出Zeeman外场和各向异性；fix_modify external energy yes 把对应磁能计入总势能，否则不能正确读取本例能垒。
• fix band all neb/spin 1.0 定义自旋路径的图像间约束。
• min_style spin 优化自旋方向，保持原子位置固定；不能换成普通原子位移最小化器来做同一个任务。
• neb/spin ... each spin-${image}.coords 为各图像读入初始方向，先普通GNEB再爬山图像。

程序通过7个MPI分区启动，一个图像一个进程。整包脚本已配好资源；已分配计算资源内的原生关系是 mpirun -np 7 lmp -partition 7x1 -in in.lammps。

4 看方向路径和磁能曲线#

head -12 work/cases/spin/physical-4.dump
cat work/cases/spin/energy-4.dat
head -10 work/analysis/spin-neb.csv

虚线是输入中弯曲的初始方向，实线和箭头来自实际输出。各图像方向模长保持1，所有原子位置仍为0，并没有发生原子迁移。

解析模型能垒为 K(1−b/(2K))²≈0.00887571064 eV，即8.875711 meV。三种优化器的实际结果一致：

迭代计数只说明这个小模型的收敛过程，不是不同硬件或软件的性能排名。

背景为本例解析能量面，节点取实际方向；初始路径的方位角弯曲被消除，路径向最低能量路线靠近。

5 换成真实磁性体系时#

需要合适的自旋交换、各向异性、外场及磁矩参数，并准备正确的多自旋初末态。验证尺寸、路径初猜、图像数和参数来源。本例只含一个局域磁矩，没有验证真实材料、多畴结构或有限温度翻转自由能。

继续：自旋优化器与力矩核对。原生接口见neb/spin、precession/spin。