ALAMODE位移与力常数：完整输入、单位和独立检查
Si64原生位移图样、正负位移、GPU原子力、extract单位转换及ALM二阶拟合。
阅读本篇目录
本篇接着ALAMODE入门，把位移、原子力和拟合文件逐步展开。完整复跑入口仍为资料包根目录的bash run.sh；以下单步命令适合在新的工作副本中理解和重算各阶段。
1 用ALM生成位移图样#
默认复跑后，进入模型目录可看到原生输入和图样：
module load alamode/1.5.0
cd results/model
cat suggest.in
cat si64.pattern_HARMONIC
本例图样完整内容为：
Basis : C
1: 1
1 1 0 0
它指定第1个原子沿笛卡尔x方向位移。图样给出方向，幅度由用户在生成结构时指定，例如0.01 Å。ALM命令为：
alm suggest.in > suggest.log
NORDER=1代表只拟合二阶力常数，即谐波项；NORDER=2才加入三阶项。不要把它理解为只存在“一阶力常数”。
本例suggest.in完整输入，含64个原子位置
&general
PREFIX = si64
MODE = suggest
NAT = 64; NKD = 1
KD = Si
/
&interaction
NORDER = 1
/
&cell
20.5224257265890415
1 0 0
0 1 0
0 0 1
/
&cutoff
Si-Si None
/
&position
1 0.0000000000000000 0.0000000000000000 0.0000000000000000
1 0.1250000000000000 0.1250000000000000 0.1250000000000000
1 0.0000000000000000 0.2500000000000000 0.2500000000000000
1 0.1250000000000000 0.3750000000000000 0.3750000000000000
1 0.2500000000000000 0.0000000000000000 0.2500000000000000
1 0.3750000000000000 0.1250000000000000 0.3750000000000000
1 0.2500000000000000 0.2500000000000000 0.0000000000000000
1 0.3750000000000000 0.3750000000000000 0.1250000000000000
1 0.0000000000000000 0.0000000000000000 0.5000000000000000
1 0.1250000000000000 0.1250000000000000 0.6250000000000000
1 0.0000000000000000 0.2500000000000000 0.7500000000000000
1 0.1250000000000000 0.3750000000000000 0.8750000000000000
1 0.2500000000000000 0.0000000000000000 0.7500000000000000
1 0.3750000000000000 0.1250000000000000 0.8750000000000000
1 0.2500000000000000 0.2500000000000000 0.5000000000000000
1 0.3750000000000000 0.3750000000000000 0.6250000000000000
1 0.0000000000000000 0.5000000000000000 0.0000000000000000
1 0.1250000000000000 0.6250000000000000 0.1250000000000000
1 0.0000000000000000 0.7500000000000000 0.2500000000000000
1 0.1250000000000000 0.8750000000000000 0.3750000000000000
1 0.2500000000000000 0.5000000000000000 0.2500000000000000
1 0.3750000000000000 0.6250000000000000 0.3750000000000000
1 0.2500000000000000 0.7500000000000000 0.0000000000000000
1 0.3750000000000000 0.8750000000000000 0.1250000000000000
1 0.0000000000000000 0.5000000000000000 0.5000000000000000
1 0.1250000000000000 0.6250000000000000 0.6250000000000000
1 0.0000000000000000 0.7500000000000000 0.7500000000000000
1 0.1250000000000000 0.8750000000000000 0.8750000000000000
1 0.2500000000000000 0.5000000000000000 0.7500000000000000
1 0.3750000000000000 0.6250000000000000 0.8750000000000000
1 0.2500000000000000 0.7500000000000000 0.5000000000000000
1 0.3750000000000000 0.8750000000000000 0.6250000000000000
1 0.5000000000000000 0.0000000000000000 0.0000000000000000
1 0.6250000000000000 0.1250000000000000 0.1250000000000000
1 0.5000000000000000 0.2500000000000000 0.2500000000000000
1 0.6250000000000000 0.3750000000000000 0.3750000000000000
1 0.7500000000000000 0.0000000000000000 0.2500000000000000
1 0.8750000000000000 0.1250000000000000 0.3750000000000000
1 0.7500000000000000 0.2500000000000000 0.0000000000000000
1 0.8750000000000000 0.3750000000000000 0.1250000000000000
1 0.5000000000000000 0.0000000000000000 0.5000000000000000
1 0.6250000000000000 0.1250000000000000 0.6250000000000000
1 0.5000000000000000 0.2500000000000000 0.7500000000000000
1 0.6250000000000000 0.3750000000000000 0.8750000000000000
1 0.7500000000000000 0.0000000000000000 0.7500000000000000
1 0.8750000000000000 0.1250000000000000 0.8750000000000000
1 0.7500000000000000 0.2500000000000000 0.5000000000000000
1 0.8750000000000000 0.3750000000000000 0.6250000000000000
1 0.5000000000000000 0.5000000000000000 0.0000000000000000
1 0.6250000000000000 0.6250000000000000 0.1250000000000000
1 0.5000000000000000 0.7500000000000000 0.2500000000000000
1 0.6250000000000000 0.8750000000000000 0.3750000000000000
1 0.7500000000000000 0.5000000000000000 0.2500000000000000
1 0.8750000000000000 0.6250000000000000 0.3750000000000000
1 0.7500000000000000 0.7500000000000000 0.0000000000000000
1 0.8750000000000000 0.8750000000000000 0.1250000000000000
1 0.5000000000000000 0.5000000000000000 0.5000000000000000
1 0.6250000000000000 0.6250000000000000 0.6250000000000000
1 0.5000000000000000 0.7500000000000000 0.7500000000000000
1 0.6250000000000000 0.8750000000000000 0.8750000000000000
1 0.7500000000000000 0.5000000000000000 0.7500000000000000
1 0.8750000000000000 0.6250000000000000 0.8750000000000000
1 0.7500000000000000 0.7500000000000000 0.5000000000000000
1 0.8750000000000000 0.8750000000000000 0.6250000000000000
/
&cell的比例因子使用bohr，&position使用分数坐标。本例谐波小模型采用Si-Si None，包含超胞内可能的二阶相互作用；高阶模型需要结合体系设置有限截断，不能直接照搬None扩大参数规模。
2 用原生displace.py生成结构#
python "$ALAMODE_ROOT/tools/displace.py" --LAMMPS Si64.lammps --mag 0.01 --prefix harm -pf si64.pattern_HARMONIC
生成的harm1.lammps对应0.01 Å位移。本资料包同时使用正、负位移，并检查0.005、0.01、0.02 Å三个幅度。ALAMODE 1.5.0的LAMMPS读取器需要明确的xy xz yz头行；包内立方胞使用三个0倾斜量，保持原几何不变。
3 在GPU上计算静态力#
资料包通过python new-case.py p100和sbatch run.sh建立并提交新力计算。每个位移结构的核心LAMMPS输入为：
package kokkos neigh half
newton on
units metal
atom_style atomic/kk
boundary p p p
read_data structure.lmp
mass 1 28.0855
pair_style sw/kk
pair_coeff * * Si.sw Si
neighbor 0.3 bin
neigh_modify every 1 delay 0 check yes
dump f all custom 1 force.dump id xu yu zu fx fy fz
dump_modify f sort id format float %24.17e
thermo 1
thermo_style custom step atoms pe press fmax
run 0
脚本在分配的GPU资源内执行：
module load lammps/30Sep2026-p100-kokkos
lmp -k on g 1 -sf kk -in force.in
Si.sw为实际使用的Stillinger–Weber势文件。metal单位下，坐标和力分别为Å、eV/Å；xu yu zu保留展开坐标，sort id保证提取时原子顺序一致。
GPU作业结束后，回到资料包根目录执行：
bash run-new.sh
入口把work/p100及可选work/v100的新原子力接入独立new-workflow，重新拟合并生成new-workflow/results和new-workflow/figures；包内原始参考数据与原results保留。只重算P100时，V100仍使用资料包提供的相同结构参考力。若修改材料或结构，应同步准备匹配的ALM结构与位移输入。
4 用extract.py提取位移—力数据#
在results这一层，新建的工作副本可执行：
python "$ALAMODE_ROOT/tools/extract.py" --LAMMPS model/Si64.lammps --unit Rydberg force-cases/p100/h010p/force.dump force-cases/p100/h010n/force.dump > DFSET_raw
head DFSET_raw
每个结构连续64行，每行6个量：
ux uy uz Fx Fy Fz
位移单位为bohr，力单位为Ry/bohr。不能把LAMMPS的eV/Å数值不转换就直接交给ALM。脚本对正负位移和每个原子的力都做了独立单位核对。
5 为什么正负位移要配对#
为了隔离局域谐波项，本例取：
u = (u_plus - u_minus) / 2
F_odd = (F_plus - F_minus) / 2
将(u, F_odd)及其反号各写一份，构成DFSET_central。原始DFSET_raw仍完整保留，另将偶分量(F_plus + F_minus)/2 - F_reference保存为even-force-eVA.csv。
这样既能拟合谐波响应，也能看到原始力中的非谐成分。本例偶分量最大绝对值随位移幅度从0.005、0.01增至0.02 Å，约为6.17e-5、2.47e-4、9.87e-4 eV/Å；它们没有被删除或冒充成零。该观察也不等于已经完成三阶力常数拟合。
6 拟合并输出Hessian#
fit.in与前述结构输入相同，关键变化为：
&general
PREFIX = si64
MODE = optimize
HESSIAN = 1
NAT = 64; NKD = 1
KD = Si
/
&optimize
DFSET = DFSET_central
LMODEL = ols
ICONST = 11
/
结构、&interaction、&cell、&cutoff、&position各段仍需保留，不能只用这个片段代替完整文件。完整文件在results/fit/u010/fit.in。
cd fit/u010
alm fit.in > fit.log
ls si64.fcs si64.xml si64.hessian
ICONST=11在拟合时施加平移不变性约束。si64.xml供ANPHON使用；si64.hessian是完整192×192二阶矩阵，单位Ry/bohr²。0.01 Å数据的训练相对残差约1.27e-6。
7 用没有参与拟合的结构检查#
另取三个64原子随机小位移结构，典型坐标扰动为0.001 Å，直接再用LAMMPS计算力；这些结构没有进入上述谐波拟合。用输出Hessian预测F=-H×u后，三组力RMSE约为4.92e-5、5.69e-5、5.26e-5 eV/Å。
这说明本例在小位移附近的谐波近似一致，不代表大位移、熔化、反应或长时间MD都适用。继续按色散、DOS与热力学使用拟合得到的XML。