Phonopy 格式单位与晶胞核对
用真实结构转换结果检查bohr与Å、LAMMPS晶胞旋转、原子顺序和周期距离，避免把几何文件当完整计算输入。
阅读本篇目录
同一个晶胞可以用Å或bohr写入，也可以整体旋转后写进LAMMPS三斜盒。判断转换是否正确，要把单位统一后比较几何，同时检查原子与元素对应。本篇直接操作配套真实输出。
1. 下载与准备#
下载结构转换和读回实例包，在已分配的CPU资源内操作：
unzip phonopy-interfaces-wave149-*.zip -d phonopy-interfaces
cd phonopy-interfaces
bash run.sh list
模块为 phonopy/4.8.3：Python3.11.16、NumPy2.3.5、SciPy1.16.3、spglib2.7.0；Phonopy C/OpenMP部分由GCC13.3编译，另有phonors0.6.0 Rust后端。CP2K结构解析采用cp2k-input-tools0.9.1的固定393e87bb版本。每次转换、分析使用已分配范围内1个CPU，不需要GPU。
一份新解压目录运行 bash run-all.sh，依次转换28个例子、读回为VASP格式、比较几何并画图；结果在 analysis/roundtrip-errors.csv、analysis/CHECKS.json 和 figures/。另一份新解压目录只提取包内真实结果：
bash analyze.sh
ls analysis-reference figures-reference
同一输出目录不会重复覆盖，继续试另一条路线或重新解压到新目录即可。包内 run.slurm 是CPU分区56核申请模板，实际工具只使用其中1核；按可用资源提交。大型超胞的计算另外申请资源。
2. 比较QE和ABACUS的长度#
bash run.sh si-qe
bash run.sh si-abacus
cat work/si-qe/structure
cat work/si-abacus/structure
Si原始非零晶胞分量是2.715Å。QE输出 CELL_PARAMETERS bohr 后约为5.130606462；ABACUS输出 LATTICE_CONSTANT 1.0 配合相同bohr数值的 LATTICE_VECTORS。两者都不是把晶胞放大了约1.89倍。
原POSCAR： 0.0 2.715000000 2.715000000 Å
QE CELL_PARAMETERS bohr： 0.0 5.130606462 5.130606462 bohr
ABACUS：1.0 bohr × (0.0 5.130606462 5.130606462)
1bohr约0.529177Å。包内分析另外用CODATA的0.529177210903Å检查原始文本，再比较读回几何；不同版本内部常数末位可能略有差异。
3. 结构单位与力单位分开看#
以下是当前Phonopy接口的默认内部约定，供接后续求力流程时查阅。本批仅转换结构，没有计算或验证这14种软件的力。输出软件本身可选择其他单位时，应以对应接口、输入标签和原始输出为准。
| 接口 | 接口默认长度 | 力 | 力常数 |
|---|---|---|---|
| VASP / FHI-aims / CASTEP / PWmat | Š| eV/Š| eV/Ų |
| LAMMPS（metal约定） | Å | eV/Å | eV/Ų |
| QE / Questaal | bohr | Ry/bohr | Ry/bohr² |
| ABACUS / ABINIT / SIESTA | bohr | eV/Å | eV/(Å·bohr) |
| CP2K | Å | Hartree/bohr | Hartree/(Å·bohr) |
| DFTB+ / Elk / exciting / Octopus | bohr | Hartree/bohr | Hartree/bohr² |
不能只看文件名 FORCE_CONSTANTS 就假定全部为eV/Ų。结构转换不转换已有力、应力或力常数文件；改变坐标系后还需要按相同旋转处理矢量/张量，不能给旋转后的结构配原坐标系的力。
4. LAMMPS：晶胞向量为什么变了#
bash run.sh si-lammps
cat work/si-lammps/structure
cat work/si-lammps/roundtrip.vasp
原Si原胞为：
0.000000000 2.715000000 2.715000000
2.715000000 0.000000000 2.715000000
2.715000000 2.715000000 0.000000000
写成受限三斜盒后读回约为：
3.839589822 0.000000000 0.000000000
1.919794911 3.325182326 0.000000000
1.919794911 1.108394109 3.135011962
两者是刚体旋转关系，体积同为40.02575175ų，晶胞边长和夹角保持；原子分数坐标也保持。把两个Cartesian矩阵直接相减会出现大差值，不能据此认定结构错误。
检查时包内计算晶胞度量矩阵 G = cell @ cell.T、体积、按顺序的元素/分数坐标和周期最短原子间距。对应代码片段：
module load phonopy/4.8.3
python - <<'PY'
import numpy as np
from phonopy.interface.vasp import read_vasp
before = read_vasp("inputs/Si.vasp")
after = read_vasp("work/si-lammps/roundtrip.vasp")
print("volume:", before.volume, after.volume)
print("metric error:", np.max(np.abs(
before.cell @ before.cell.T - after.cell @ after.cell.T)))
delta = after.scaled_positions - before.scaled_positions
delta -= np.rint(delta)
print("fractional coordinate error:", np.max(np.abs(delta)))
PY
本例的27个相邻周期像足够检查所选温和三斜晶胞；极端倾斜晶胞应使用可靠的通用最小镜像算法，不能照搬有限邻居枚举。
5. 双元素原子顺序不能省略检查#
bash run.sh nacl-abacus
bash run.sh nacl-qe
cat inputs/NaCl.vasp
cat work/nacl-abacus/roundtrip.vasp
本例顺序为Na、Cl，实际转换后按相同顺序读回。对真实多元素超胞，要同时核对元素列表、每种元素的数量、位置及后续力的排列；格式读写器可能按元素分组，不能只确认总原子数相同。本例两原子不能证明任意交错元素顺序都会保持。
6. 输出与检查范围#
bash run-all.sh
head -5 analysis/roundtrip-errors.csv
python -c 'import json; a=json.load(open("analysis/CHECKS.json")); print(a["conditions"], a["unique_formats"])'
ls figures
结果为28个结构实例、14条格式路线，最大读回边长差2.20×10⁻⁹Å、分数坐标差5.14×10⁻¹⁴、周期距离差5.84×10⁻¹⁰Å，均小于预先选定的几何写读精度范围。图中零误差点只为对数显示放到10⁻¹⁶Å，原始CSV保留真实数值。
WIEN2k、Fleur、CRYSTAL和TURBOMOLE需要各自额外元数据、输入模板或文件组织，本包不把这些接口算作已完成的端到端路线；VASP同格式重写也不单独计为新功能。转换时不传递磁矩、约束、电子参数、基组或势函数，请按目标程序重新检查。
参考Phonopy接口与单位表；跨软件转换操作给出完整POSCAR、14条命令及QE/CP2K步骤。