LAMMPS神经网络势的力与压力分析

在神经网络势实操的完整包目录中，先从参考结果生成表格：

bash analyze.sh

head -n 12 tables/energies-forces.csv

energy_peratom_eV 是整盒能量除以原子数。pressure_bar 是LAMMPS原生压力，FD_pressure_bar 则由两次额外固定结构能量求导得到。两者的来源不同，适合相互检查。

1 铜模型：从描述符重新计算网络能量#

head -n 12 reference/potentials/Cu.nn.mliap.model

铜网络有30个输入、两个各2节点的tanh隐藏层和1个线性输出，共71个可训练参数；文件还保存两组30维描述符归一化数组。分析脚本从原生 sna/atom 输出读取B，再按 (B−scale0)/scale1 归一化，逐层计算 tanh(Wx+b)，最终获得每个原子的能量。

这项前向计算使用NumPy，未调用LAMMPS网络求值器；它复用了原生描述符，没有独立重写SNA的球谐展开。硅SO3模型不提供这一项独立前向对照，CSV相应字段留空，不能把空值当作零误差。

铜16个参考条件的整盒重建能量最大差约1.09×10⁻¹² eV；图的左侧每个点是一个真实条件。索引依次为P100/V100，再为0.98/1.00/1.02/扰动结构，每个结构下为普通/Kokkos路径。

2 原子力与坐标能量导数#

输入对原子1的x、17的y、43的z坐标分别作±10⁻⁵ Å位移，固定其它坐标和晶胞，用 −[E(x+δ)−E(x−δ)]/(2δ) 与原生力分量比较。先完成单坐标±δ，再恢复该坐标，随后处理下一个坐标。

所有32个参考条件的最大差约1.14×10⁻⁸ eV/Å，低于本例10⁻⁵ eV/Å的绝对容差。这里核对的是预选的三个分量，不是所有坐标或所有可能构型。理想高对称结构的力可能接近零，扰动结构更有助于观察非零受力。

3 压力与能量—体积导数#

对当前结构施加相同的三个方向微应变ε=10⁻⁵，坐标随盒子仿射缩放，得到E₊、V₊和E₋、V₋。以 −(E₊−E₋)/(V₊−V₋) 得到近似压力，再从eV/ų换算成bar。它对应固定分数坐标的各向同性变形，没有在每个体积重新优化原子。

点来自实际能量与压力，虚线是相等关系。最大压力差约0.003215 bar，本例绝对容差0.01 bar。图用GPa显示，1 GPa=10000 bar。

这些零步输入没有速度和动能，原生压力对应构型维里项；有温度的MD还包含动能压力，不能直接用本页固定结构导数去对照完整的有限温度压力。求有限温度自由能、热膨胀或EOS需要相应采样和统计，不能从这组微应变数值直接推出。

4 找到自己的原始文件#

reference/条件名/base.dump 包含基准坐标、原子力和逐原子能量；铜文件还包含30个SNA分量。properties.dat 保存±位移和±应变的实际总能量、体积。自己的结果保存在对应的 new-work/ 路径。

重新绘制包内参考图：

bash plot.sh

同模型的普通/Kokkos及P100/V100结果可以对照；不同模型、不同结构、不同能量零点不能混合。输入和算法的详细选项见LAMMPS ML-IAP和压力定义。