SIESTA：原子优化与晶胞优化
分别实际运行固定晶胞与可变晶胞CG，说明力、应力、最终XV和晶格参数的关系。
阅读本篇目录
结构优化应先明确：只移动原子，还是连晶胞一起改变。以下两例都使用Si2、PBE、DZP、150 Ry、4×4×4 k网格；电子自由能差10^-6 eV，原子力0.02 eV/Å，可变晶胞应力0.5 GPa。
1. 固定晶胞，只优化原子位置#
初始a=5.43 Å，第二个Si沿z偏移0.08 Å。除入门教程完整电子输入外，原子优化添加：
MD.TypeOfRun CG
MD.NumCGsteps 80
MD.MaxForceTol 0.02 eV/Ang
MD.MaxStressTol 0.5 GPa
MD.VariableCell false
module load siesta/5.4.2-8173m-mpi sisl/0.16.4
python prepare.py
cat cases/atoms-opt/input.fdf
MD.VariableCell false保留晶胞，优化自由原子；MD前缀沿用软件历史命名，这里运行的是CG几何优化，不是分子动力学。
2. 连同晶胞一起优化#
另一个目录从a=5.30 Å开始，改用：
MD.TypeOfRun CG
MD.NumCGsteps 80
MD.MaxForceTol 0.02 eV/Ang
MD.MaxStressTol 0.5 GPa
MD.VariableCell true
cat cases/cell-opt/input.fdf
sbatch job.slurm
squeue -u "$USER"
可变晶胞允许晶格矢量和原子同时变化，停止时既要满足力，也要满足应力。具体含义见官方结构优化教程。
3. 查看实际收敛结果#
cat cases/atoms-opt/si.FA
cat cases/cell-opt/si.FA
cat analysis/cases.csv
cat analysis/GEOMETRIES.json
| 项目 | 固定晶胞原子优化 | 可变晶胞优化 |
|---|---|---|
| 初始常规a/Å | 5.43 | 5.30 |
| 最终等效常规a/Å | 5.43 | 5.500482 |
| 最大力/eV·Å^-1 | 0.012229 | 0.002408 |
| 最大静态应力/GPa | 3.68180 | 0.08171 |
固定晶胞的最终力满足0.02，但应力仍较大；这不是晶胞优化失败，而是它没有被允许改变晶胞。可变晶胞的力与应力均达到本例要求。
4. 读取XV并继续使用#
python - <<'PY'
import sisl
for name in ["atoms-opt","cell-opt"]:
g=sisl.get_sile("cases/"+name+"/si.XV").read_geometry()
print(name,"cell Angstrom",g.cell,"positions Angstrom",g.xyz)
g.write(name+"-final.xyz")
PY
XV有原生单位，建议用可靠读取器转换为Å，不要直接把文本中的Bohr坐标复制成Å。后续单点应使用最终晶胞和最终原子位置，并按目标性质重新确定基组和网格。
5.500482 Å是本组有限基组与参数下的操作结果，不能视作严格收敛的Si晶格常数；局域基组还可能有Pulay应力，应继续检查基组、网格与能量偏移。
下载Si2实算资料包，在新目录解压，执行sbatch job.slurm申请CPU资源。分区按56核整节点分配；本例小任务每条件2个MPI进程，最多四条件并发。参考数据同样采用2个MPI进程／条件。软件为SIESTA5.4.2 CPU MPI（oneAPI2024.2、IntelMPI2021.13、MKL LP64顺序库、AVX512/O2），后处理sisl0.16.4。包内无软件程序或赝势，准备步骤读取公共PseudoDojo0.4 PBE标量相对论Si PSML，并记录哈希。