LAMMPS混合规则与hybrid势实操

给不同原子类型设置势函数时，必须明确“异类原子之间用什么参数”和“多个势如何组合”。本篇用同一个4粒子几何，实际比较三种LJ混合规则、显式交叉参数、hybrid和hybrid/overlay。先走通输入、提交、写表和核对，再把规则迁移到自己的力场。

1. 下载、预览与提交#

下载六组完整输入、原生数值和绘图脚本，在新目录解压：

python3 -m zipfile -e lammps-mixing-hybrid-wave304-*.zip .
cd lammps-mixing-hybrid-wave304
bash preview.sh
sbatch run.slurm
squeue -u "$USER"

preview/重画两张参考图；作业内run.sh顺序运行6份真实LAMMPS输入，输出原生坐标/力、总能和每条件241点交叉势表，再自动分析、画图。native-results/是这次新计算，results/是完整NPZ/CSV/JSON，plots/是新图；reference-results/是附件参考数值。复跑使用新的解压目录。

LAMMPS30Sep2026采用GNU13.3/OpenMPI5.0.10构建，本例加载lammps/30Sep2026-v100-gpu却调用普通CPU样式。模板为V10016分区1GPU1CPU，GPU只随模板分配，不能将此例当成GPU加速测试。后处理用pymatgen/2026.9.24中的NumPy与Matplotlib。

2. 同一个几何，明确两种类型#

4粒子类型依次为1、1、2、2。坐标为(0,0,0)、(1.8,0.1,0)、(3.7,−0.2,0.1)、(5.6,0.2,−0.1)。没有周期边界或时间积分，所有数值采用LJ约化单位；不是特定材料的Å、eV或ps。

自相互作用参数为：1–1的ε=1、σ=1；2–2的ε=0.5、σ=1.4。截断均为8，不作能量平移。关键是接下来怎样定义1–2相互作用。

左图为实际x–y投影，各点标出真实z坐标；中图画的是全部4个原子的Fx，完整三分量保存在NPZ中。连接点只帮助比较，不是额外计算了连续轨迹。

3. arithmetic完整输入#

cat inputs/arithmetic/in.lammps

units lj
atom_style atomic
atom_modify map array sort 0 0.0
boundary f f f
region box block -2 8 -3 3 -3 3 units box
create_box 2 box
create_atoms 1 single 0 0 0 units box
create_atoms 1 single 1.8 .1 0 units box
create_atoms 2 single 3.7 -.2 .1 units box
create_atoms 2 single 5.6 .2 -.1 units box
mass * 1
pair_style lj/cut 8
pair_coeff 1 1 1 1
pair_coeff 2 2 .5 1.4
pair_modify mix arithmetic
neighbor .3 bin
neigh_modify every 1 delay 0 check yes
thermo_style custom step atoms pe
thermo_modify norm no format float %.16g
run 0
write_dump all custom atoms.dump id type x y z fx fy fz modify sort id format float %.16g
print "TOTAL_ENERGY $(pe:%.16g)" file energy.dat screen no
pair_write 1 2 241 r .9 5 cross.table CrossPair

只给1 1与2 2，由pair_modify mix arithmetic生成交叉项。run 0完成固定结构的力和能量计算；随后pair_write输出真实1–2势表，atoms.dump保留所有原子力。

4. 三种混合规则有什么区别#

对本例LJ势，参数定义为：

注意arithmetic并不是把ε也取算术平均。对应输入只改变一行：

pair_modify mix geometric

或：

pair_modify mix sixthpower

这三份实际输入已经分别保存在inputs/geometric、arithmetic和sixthpower中。它们给出不同的势能和力，原因是定义了不同的交叉参数，不能据此评价哪条混合规则对任意真实材料都更准确。

5. 显式交叉参数会覆盖自动生成#

第四个条件直接指定算术规则对应的交叉系数：

pair_style lj/cut 8
pair_coeff 1 1 1 1
pair_coeff 2 2 .5 1.4
pair_coeff 1 2 0.7071067811865476 1.2

实际结果中，它与arithmetic条件的总能和全部原子力相同。这是检查输入是否按预期生效的实用方法：把你认为的混合结果显式写出，再用同一几何比较。

cat results/DATA.json
head results/arithmetic-cross.csv
head results/explicit-arithmetic-cross.csv

6. hybrid：为不同类型对分配不同势#

本例hybrid的关键部分：

pair_style hybrid lj/cut 8 morse 8
pair_coeff 1 1 lj/cut 1 1
pair_coeff 2 2 lj/cut .5 1.4
pair_coeff 1 2 morse .7 1.6 1.9

1–1和2–2用LJ，1–2用Morse。Morse三个参数依次为D=0.7、α=1.6、r₀=1.9。对这个输入，每种类型对只分给一个子势。

不要假设把两个势名字并列后，软件就会自动得到合理的跨类型相互作用。仅在匹配的子势和其混合条件允许时才能生成交叉参数；跨子势、复杂多体模型以及所需物理耦合必须明确设计。

7. overlay：同一类型对叠加贡献#

本例overlay的关键部分：

pair_style hybrid/overlay lj/cut 8 morse 8
pair_coeff * * lj/cut .8 1.1
pair_coeff 1 2 morse .7 1.6 1.9

所有类型对先有ε=0.8、σ=1.1的LJ；只有1–2再加一份Morse。因此1–2能量是两项之和，1–1与2–2没有Morse贡献。本例特意把系数写清楚，不依赖多子势情况下的隐式混合。

overlay叠加能量不等于它天然更符合真实材料；若两个势描述了同一物理贡献，叠加会重复计算。多体势还需要核对环境与元素映射，不能只用本例两体势的直觉拼接。

图中保留全部0.9–5的241个原生点，使用能容纳正负值的symlog纵轴，既展示短程排斥，也保留势阱；曲线没有被人为归一成同一深度。显式算术条件与arithmetic重合。

8. 本次实际结果与核对#

独立程序对全部6个粒子对求和，并核对4×3个力分量及每条件241个势能/径向力点。总能和力按1E−10归一误差检查，势表按1E−8检查。这里核对的是输入、单位与实现；不是某种材料的参数精度或统计收敛。

head native-results/hybrid/cross.table
cat native-results/hybrid/energy.dat
head native-results/hybrid/atoms.dump

势表最后一列是沿粒子间距的径向力−dE/dr，不是某个原子的Fx。逐原子力还需要乘方向，并把所有邻居贡献求和。表中r=5仅是本次写表范围终点，实际势截断为8。

9. 单独分析自己的新结果#

module purge
module load pymatgen/2026.9.24
python analyze.py --data native-results --output my-results
python plot.py --data my-results --output my-plots

脚本专门对应本例4粒子、6种设置和已有参数；换体系时同步修改解析参照和检查条件。参考文件不会自动变成新计算，始终从自己的native-results/提取。

继续看键、角与二面角势能表、长程静电和薄层slab修正。规则定义见LAMMPS pair_modify与hybrid官方说明。