SIESTA晶体声子：超胞、有限位移与色散图
用硅的54原子超胞串起fcbuild、SIESTA和Vibra，展示完整输入与终端操作，比较两种位移并提取6支声子。
阅读本篇目录
晶体声子需要原胞、超胞力常数与波矢路径。本例从2原子硅原胞生成3×3×3的54原子超胞，分别计算0.01/0.02Å有限位移，再得到Gamma–X–W–L–Gamma的6支声子。它与水分子振动的单Gamma点9个模式不同。
1. 下载并提交完整流程#
unzip siesta-phonons-wave203.zip -d si-phonon
cd si-phonon
sbatch job.slurm
CPU分区可用时提交。环境为SIESTA5.4.2、Intel oneAPI2024.2、Intel MPI2021.13、MKL顺序LP64，配套sisl0.16.4读取结构；fcbuild和vibra来自同一SIESTA安装。两种位移最多同时计算，每组4个MPI进程，每进程1线程。
程序依次执行：准备原胞 → fcbuild生成超胞 → SIESTA有限位移 → vibra计算色散 → 提取CSV并画图。
2. 原胞与波矢路径：phonon.fdf#
SystemLabel si54
NumberOfAtoms 2
LatticeConstant 5.500481763 Ang
%block LatticeVectors
0.0 0.5 0.5
0.5 0.0 0.5
0.5 0.5 0.0
%endblock LatticeVectors
AtomicCoordinatesFormat ScaledCartesian
%block AtomicCoordinatesAndAtomicSpecies
0.0 0.0 0.0 1 28.0855
0.25 0.25 0.25 1 28.0855
%endblock AtomicCoordinatesAndAtomicSpecies
SuperCell_1 1
SuperCell_2 1
SuperCell_3 1
AtomicDispl 0.01 Ang
Eigenvectors true
BandLinesScale ReciprocalLatticeVectors
%block BandLines
1 0.0 0.0 0.0 Gamma
20 0.5 0.0 0.5 X
20 0.5 0.25 0.75 W
20 0.5 0.5 0.5 L
20 0.0 0.0 0.0 Gamma
%endblock BandLines
ScaledCartesian坐标乘以晶格常数，末列是Si质量；不要把它当作以非正交原胞晶格矢量为基的分数坐标。晶格常数来自配套硅结构优化的同类PBE/DZP示例，不等于实验常数。
SuperCell_1/2/3=1代表每个方向包含−1、0、+1三个原胞，因此总原子数是2×3³=54，不是1×1×1。第二组只把AtomicDispl改成0.02Å。
3. fcbuild生成超胞输入#
在作业分配的计算资源内，自动流程中的原生命令是：
fcbuild < phonon.fdf > fcbuild.out
cat FC.fdf
FC.fdf整体描述54原子超胞、晶胞矢量、坐标和位移范围。它自动指定FC.First=27、FC.Last=28，即中心原胞的两个原子。位移中心原胞后，记录对整个超胞的力响应；不需要对54个原子逐个做重复的等价位移。
不要只复制FC.fdf中的坐标而丢掉相配的晶胞、单位和原子编号。工具生成的坐标已写成明确的小数格式。
4. SIESTA计算力常数#
输入的其余部分如下，最后整体包含工具生成的FC.fdf：
SystemLabel si54
NumberOfSpecies 1
%block ChemicalSpeciesLabel
1 14 Si
%endblock ChemicalSpeciesLabel
XC.functional GGA
XC.authors PBE
PAO.BasisSize DZP
PAO.EnergyShift 100 meV
MeshCutoff 150 Ry
ElectronicTemperature 300 K
MaxSCFIterations 300
MinSCFIterations 3
SCF.DM.Converge true
SCF.DM.Tolerance 1.e-5
SCF.FreeE.Converge true
SCF.FreeE.Tolerance 1.e-6 eV
DM.MixingWeight 0.1
DM.NumberPulay 6
DM.UseSaveDM false
%block kgrid_Monkhorst_Pack
1 0 0 0
0 1 0 0
0 0 1 0
%endblock kgrid_Monkhorst_Pack
%include FC.fdf
超胞这里使用Gamma电子k点采样。电子k点与下一步声子色散的波矢不是同一件事。基组、网格和电子k点会影响力常数，实际研究应先检查它们的收敛。
mpirun -np 4 siesta < input.fdf > siesta.out
grep 'SCF cycle converged' siesta.out
ls -lh si54.FC
这些命令已由下载包在作业内自动执行；不要在没有申请计算资源的登录终端重复启动。每个中心原子沿三个方向正负位移，输出si54.FC保留全部力常数。
5. Vibra生成声子色散#
vibra < phonon.fdf > vibra.out
head si54.bands
head si54.vectors
si54.bands的6支频率对应原胞2×3个自由度，而不是54×3条。si54.vectors给出振型。自动流程对两种位移都计算并保留原始输出，然后生成可直接画图的CSV。
head analysis/fc010-dispersion.csv
cat analysis/SUMMARY.json
ls figures
| 位移（Å） | Gamma光学1（cm⁻¹） | 光学2 | 光学3 |
|---|---|---|---|
| 0.01 | 533.024 | 533.024 | 533.029 |
| 0.02 | 532.541 | 532.541 | 532.546 |
两种位移的Gamma点光学频率最大相对差异为0.091%。独立程序从原始力常数重建Gamma点动力学矩阵，与原生光学频率相符；平移求和规则使Gamma声学模式接近零。
6. 如何用于自己的体系#
先优化结构，再替换原胞、物种、质量和匹配赝势；选择合理超胞和波矢路径，运行上述三步，最后比较位移幅度及其他数值参数。原始正负频率都应保留，出现负频时继续检查几何、数值精度和真实不稳定性，不要直接取绝对值。
本例只比较两种位移，没有完成超胞、基组、网格或k点收敛，也没有计算热导率。它提供可复用的实际操作流程，不代替特定材料的收敛研究。官方声子教程也强调超胞与数值精度的影响。