LAMMPS局域超动力学操作

本组使用LAMMPS30Sep2026的原生REPLICA功能：四个粒子中，1号粒子沿x跨越解析双阱，另外三个局域谐振子作为陪伴自由度。所有粒子还在各自的y/z势阱内运动。使用LJ约化单位，图中X是占位标签，不代表某种真实元素；能垒1不是1 eV。

模型为1号粒子的 U=(x²−1)²+5(y+3)²+5z²，其他粒子为 U=5(x+1)²+5[y−(2id−5)]²+5z²。初态x=−1，另一稳态x=+1，解析能垒为1。这里用来学习输入、事件和时钟；材料扩散率还需要适合材料的势函数、事件定义、去相关和统计收敛。

这个解析外场模型仅演示原生命令、偏置表达式、事件文件和方法时钟。邻近对、qfactor和Vmax的选择没有证明满足真实过渡态处偏置为零等动力学条件；以下方法时钟不能作为本模型已验证的真实跃迁时间或材料扩散率。实际材料要检查这些条件、事件去相关和参数收敛。

下载原生输入、实际输出、CSV和图片，解压到个人目录：

mkdir hyper-demo
cd hyper-demo
curl -fL https://yuansuan.top/assets/downloads/full-function/lammps-hyper-wave84-73ba3ce03eab.zip -o examples.zip
unzip examples.zip
sbatch submit.slurm

提交脚本会依次运行各个小算例，然后自动提取数据、导出结构和画图。使用 lammps/30Sep2026-p100-gpu 基础模块，原生GCC13.3、OpenMPI5.0.10、CUDA12.8构建；后处理用 ase/3.29.0。这些命令走CPU/MPI路径，并不因模块名字含gpu而自动使用GPU力计算。脚本申请一个CPU核和一个P100配额；三个方法在单进程内依次运行。

完成后，在解压目录运行：

module purge
module load ase/3.29.0
python analyze.py work
python export_events.py work
python plot_results.py work

已有实算结果位于 reference/，想先读图时可把后处理参数work换为reference。新计算请重新解压到新目录，保留自己的既有结果。

1 从全局模式换到局域模式#

局域模式可以在不同局部区域选择多个偏置对，并调整每个邻近对的Cij系数。本例仍是四粒子操作模型，用于读懂选项和输出，不用于证明大体系的空间独立性或动力学加速效果。

2 完整局域输入#

variable replica equal 1
units lj
atom_style atomic
atom_modify map array sort 0 0.0
boundary f f f
region box block -5 5 -8 8 -5 5 units box
create_box 1 box
mass 1 1
create_atoms 1 single -1 -3 0 units box
create_atoms 1 single -1 -1 0 units box
create_atoms 1 single -1 1 0 units box
create_atoms 1 single -1 3 0 units box
pair_style zero 2.5
pair_coeff * *
neighbor 0.3 bin
neigh_modify every 1 delay 0 check yes
variable dy atom y-(2*id-5)
variable active atom id==1
variable U atom v_active*(x*x-1)*(x*x-1)+(1-v_active)*5*(x+1)*(x+1)+5*v_dy*v_dy+5*z*z
variable Fx atom -v_active*4*x*(x*x-1)-(1-v_active)*10*(x+1)
variable Fy atom -10*v_dy
variable Fz atom -10*z
fix field all addforce v_Fx v_Fy v_Fz energy v_U
fix_modify field energy yes
compute physicalT all temp
compute_modify physicalT extra/dof 0
thermo 1000
thermo_style custom step c_physicalT pe ke
thermo_modify norm no format float %.16g
displace_atoms all random .01 .01 .01 19475 units box
velocity all create 0.3 32613 mom no rot no dist gaussian temp physicalT
fix integrate all nve
variable bathseed equal 12341+17*v_replica
fix bath all langevin 0.3 0.3 0.2 ${bathseed} zero no
compute event all event/displace 1.0
timestep .005
min_style cg
comm_modify cutoff 4.0
dump events all custom 1000000 events.dump id type x y z
dump_modify events sort id format float %.16g
fix boost all hyper/local 2.3 .10 .10 .3 1.5 10.0 1.3 bound .1 reset 0 check/bias 20 error
variable g0 equal f_boost
variable b0 atom v_g0
variable g1 equal f_boost[1]
variable b1 atom v_g1
variable g2 equal f_boost[2]
variable b2 atom v_g2
variable g3 equal f_boost[3]
variable b3 atom v_g3
variable g4 equal f_boost[5]
variable b4 atom v_g4
variable g5 equal f_boost[6]
variable b5 atom v_g5
variable g6 equal f_boost[7]
variable b6 atom v_g6
variable g7 equal f_boost[23]
variable b7 atom v_g7
variable g8 equal f_boost[24]
variable b8 atom v_g8
variable g9 equal f_boost[4]
variable b9 atom v_g9
dump trajectory all custom 1 trajectory.dump id type x y z vx vy vz v_b0 v_b1 v_b2 v_b3 v_b4 v_b5 v_b6 v_b7 v_b8 v_b9
dump_modify trajectory sort id format float %.16g
hyper 30000 20 boost event min 0.0 1.0e-8 500 1000 dump events
print "FINAL_HYPER $(f_boost[23]:%.16g)" file hyper-time.dat screen no
undump trajectory
unfix boost
unfix bath
run 0
write_dump all custom final.dump id type x y z vx vy vz fx fy fz modify sort id format float %.16g
print "COMPLETE $(step:%.0f) $(pe:%.16g)" file complete.dat screen no

核心两行：

fix boost all hyper/local 2.3 .10 .10 .3 1.5 10.0 1.3 bound .1 reset 0 check/bias 20 error
hyper 30000 20 boost event min 0.0 1.0e-8 500 1000 dump events

前四个参数是cutbond、qfactor、Vmax和温度；接着是局域选择距离Dcut=1.5、boostostat时间参数alpha=10、目标方法时钟倍数Btarget=1.3。真实材料要用适合相互作用和事件尺度的Dcut，而不是照搬1.5。

comm_modify cutoff 4.0 显式提供邻域通信范围；分布式计算时要覆盖Dcut与事件位移等所需范围。bound .1 限制Cij在0.9–1.1；check/bias 20 error 定期检查偏置选择的重叠。

3 reset 0不是永远固定Vmax#

本次 reset 0 在事件导致邻近对重建时，可以用当时Cij平均值更新Vmax，并重新调整系数。输入Vmax=0.10，实际最终Vmax约0.113222090。因此分析瞬时偏置时，应读原生当前Vmax列，不能一直代入0.10。

Step0的部分系数统计尚未累计，图中从实际动力学步骤开始展示；不要把初始化统计0误认为Cij越过下界。

4 运行和读数据#

在已分配的计算资源内单独执行：

module purge
module load lammps/30Sep2026-p100-gpu
mkdir manual-local
cp inputs/local/in.lammps manual-local/
cd manual-local
lmp -in in.lammps
cd ..

整包结果：

head -6 work/analysis/local-history.csv
cat work/cases/local/hyper-time.dat
head -10 work/analysis/local-events.csv

目标方法时钟按 Btarget×dt×步数 报告，因此150×1.3=195。瞬时邻域boost受适应过程和偏置参数影响，本次平均并不等于1.3；不能把它与全局模式的逐步积分公式混为一谈。这个短例不能证明已达到目标动力学精度。

5 看事件与回跳#

python export_events.py work
head -8 work/analysis/local-events.xyz
python plot_results.py work

局域XYZ含初态和3个实际淬火事件，共4帧。原始轨迹中可以看到短时间回跳：这正是为什么不能仅数事件就宣称获得可靠稀有事件速率。

图中细线来自原生连续动力学坐标，红点来自淬火稳态；三条曲线的差异只表示这三个有限条件的实际输出，不是效果排名。若用于材料研究，应检查过渡态处偏置、局域独立性、回跳与时间参数，再比较收敛后的统计量。

全局模式与偏置积分；原生参数见 LAMMPS hyper/local。