Critic2 随机配对位移与力常数回归
实际DFT计算正负位移，比较全数据拟合、留出数据预测和指定正则强度。
阅读本篇目录
随机位移路线让超胞内各原子沿不同方向移动，再用实际受力拟合力常数。本例生成2组随机位移及其负位移，实际运行4个Si16 QE单点，并比较哪些数据参与拟合、哪些只用来检查预测。
1. 准备与完整重算#
unzip critic2-dft-forces-wave146-*.zip -d critic2-dft-forces
cd critic2-dft-forces
bash run.sh list
Critic2模块为 critic2/1.4（内部1.4.155，GCC/GFortran13.3、OpenMP、LAPACK/BLAS）；QE模块为 qe/7.6-8173m（InteloneAPI2024.2、MPI+OpenMP、MKL）。对照与绘图采用 phono3py/4.8.2 独立Python3.11.16环境内的Phonopy4.8.1与Symfc1.7.3。完整DFT阶段需要已分配至少4个CPU：单MPI进程×4OpenMP线程；Critic2后处理仅使用其中1核、1线程。普通计算按账户资源使用Slurm，包内提供 run.slurm 模板。
在新解压目录运行：
bash run-all.sh
该包的随机分支实际执行4个PBE/ONCV单点；另2个是参考结构与单有限位移例，合计6个SCF。每个结构使用相同的40/160Ry截止能、2³电子k网格和1d-10Ry电子收敛，避免把不同电子设置混进同一次拟合。
2. 生成可复现的正负位移#
bash run.sh 3对应的完整原生输入：
crystal ../../data/unit.vasp
units bohr
vibrations calculator qe
vibrations memory 64 MB
vibrations create_displacements 2 2 2 distance .02 random 2 seed 260630 plusminus template rnd-*.scf.in dataset random_disp.txt
RANDOM 2生成两组随机方向；每个原子的位移长度均为0.02bohr。SEED 260630固定本例随机种子；PLUSMINUS再生成两组严格相反的位移，最后得到4个位移结构，而不是2个。
实际配对顺序为：
| 正位移 | 对应负位移 |
|---|---|
| rnd-001 | rnd-003 |
| rnd-002 | rnd-004 |
不是把相邻的001/002视为一对。原始位移向量和顺序在 random_disp.txt中。
head -n 30 reference/work/random-generate/random_disp.txt
grep 'convergence has been achieved' reference/qe-work/rnd-001/out.out
grep 'convergence has been achieved' reference/qe-work/rnd-003/out.out
Critic2先对每份总力扣除整体漂移，再构造配对力：
F_pair(u) = [ F(+u) - F(-u) ] / 2
这会抵消与位移无关的参考力和偶次位移力项，其中包括三阶势能产生的二次力项；线性力及更高的奇次力仍保留，不能理解成消除了所有非谐性。
3. 两对数据全部拟合#
只处理附带真实力，在新解压目录运行：
bash run.sh 4
cat work/random-all/points.dat
完整输入：
crystal ../../data/unit.vasp
units bohr
vibrations calculator qe
vibrations memory 64 MB
vibrations read_forces dataset random_disp.txt NORIDGE
vibrations write_fc2 FORCE_CONSTANTS full
vibrations calcq file points.dat
q 0 0 0
q .5 .5 0
q .5 .5 .5
q .13 .27 .31
ENDCALCQ
本例对称性把FC2归并为 8个独立参数。两对配对数据提供96个笛卡尔力方程；NORIDGE明确选择无正则最小二乘。参数数目由结构、超胞和对称性决定，换体系后不再固定为8。
4. 只用第一对拟合，预测第二对#
bash run.sh 5
grep 'Training force residual' work/random-one-pair/out.cro
grep 'Held-out data' work/random-one-pair/out.cro
crystal ../../data/unit.vasp
units bohr
vibrations calculator qe
vibrations memory 64 MB
vibrations read_forces dataset random_disp.txt PAIRS 1 NORIDGE
vibrations write_fc2 FORCE_CONSTANTS full
vibrations calcq file points.dat
q 0 0 0
q .5 .5 0
q .5 .5 .5
q .13 .27 .31
ENDCALCQ
PAIRS 1表示只拟合第1组正负配对数据，48个方程；第2组不参与拟合，用来检查预测。不能把训练残差小当作对其他结构同样准确。
5. 指定正则强度#
bash run.sh 6
grep 'Training force residual' work/random-ridge/out.cro
grep 'Held-out data' work/random-ridge/out.cro
crystal ../../data/unit.vasp
units bohr
vibrations calculator qe
vibrations memory 64 MB
vibrations read_forces dataset random_disp.txt PAIRS 1 RIDGE 0.0001
vibrations write_fc2 FORCE_CONSTANTS full
vibrations calcq file points.dat
q 0 0 0
q .5 .5 0
q .5 .5 .5
q .13 .27 .31
ENDCALCQ
RIDGE 0.0001为本例的显式无量纲正则设置。它不是所有材料都适用的最优参数，需结合未拟合数据的误差选择。
默认自动交叉验证要求至少10组用于拟合的数据。使用PLUSMINUS时，每个正负配对算1组，因此至少需要20个位移结构输出。本例只有2对，明确指定 NORIDGE或 RIDGE，没有用这4个SCF宣称已经完成自动交叉验证。
本次实际力残差，经全量FC2独立回算并换成eV/Å：
| 拟合方式 | 训练力RMSE / eV·Å⁻¹ | 未拟合的第2对RMSE / eV·Å⁻¹ |
|---|---|---|
| 2对全部拟合，NORIDGE | 0.00022411 | —（没有留出数据） |
| 只拟合第1对，NORIDGE | 0.00023174 | 0.00022288 |
| 只拟合第1对，RIDGE=0.0001 | 0.00023315 | 0.00022444 |
两对全部用于拟合时，没有留出数据，表中的缺项不是0。单对拟合的第2对结果是实际预测，未把它用于参数拟合。
6. 看输出、画图和替换自己的数据#
bash analyze.sh
ls figures-reference
完整重算则查看 analysis/random-one-pair-forces.csv和 figures/。CSV列出真实配对DFT力、由原生全量FC2回算的力及两者差值。独立Symfc拟合也使用相同的真实力和相同的无正则条件进行比较。
用于自己的材料时，先确定结构与位移幅度，再生成完整正负配对，保持电子参数一致，按dataset顺序传递对应总力输出，最后检查训练与未拟合数据。增加结构多样性、数据数量和收敛检查后，才进一步讨论可信的材料声子或热力学。这里是固定小体系的谐振流程参考，不是完整非谐力场训练或收敛物性结论。
完整重算的数值与图位于 analysis/CHECKS.json 和 figures/；只分析附带真实数据时，对应为 analysis-reference/CHECKS.json 与 figures-reference/。