ALAMODE自洽声子：温度与频率
用SiSW四阶力常数计算0、300、600K声子频率，展示完整输入、量子和经典占据及网格对照。
阅读本篇目录
高阶拟合实例得到FC2、FC3和FC4后，可以继续用ANPHON计算温度相关的自洽声子。本例沿用Si64 经典SW势的力常数，在固定晶胞下计算0、300、600K频率。
“使用经典SW势”说明势能模型；CLASSICAL=0/1选择核振动的量子或经典统计。这两个设置的含义不同。本页展示四阶loop项的SCPH流程，采用 BUBBLE=0、RELAX_STR=0。
1. 下载并运行完整例子#
下载SCPH输入、力常数、脚本与实际结果（ZIP）。在新目录解压，然后进入已分配、支持AVX-512的CPU环境：
unzip alamode-scph-wave189-*.zip -d alamode-scph
cd alamode-scph
module load alamode/1.5.0
bash run.sh
脚本依次准备输入、运行谐波基线和五组SCPH条件、导出有效二阶力常数、复算色散并画图。默认不提交GPU作业，不重新计算DFT或原子力。
| 目录 | 内容 |
|---|---|
model-ifcs.xml |
来自高阶实例的完整力常数，保持原件 |
cases/quantum-q4/ |
4³积分网格、量子统计的SCPH计算 |
effective-fc2/ |
温度相关有效FC2与色散复算 |
analysis/ |
频率CSV、汇总及检查结果 |
figures/ |
三张结果图 |
2. 看清完整输入#
下面是本次量子统计、4³积分网格的实际 input.in：
&general
PREFIX = scph
MODE = SCPH
NKD = 1; KD = Si
MASS = 28.0855
FCSXML = ifcs.xml
CLASSICAL = 0
TMIN = 0; TMAX = 600; DT = 300
RESTART = 0
/
&cell
10.2612128632945208
0 0.5 0.5
0.5 0 0.5
0.5 0.5 0
/
&kpoint
1
G 0 0 0 X 0.5 0.5 0 51
X 0.5 0.5 0 L 0.5 0.5 0.5 51
L 0.5 0.5 0.5 G 0 0 0 51
/
&scph
KMESH_SCPH = 4 4 4
KMESH_INTERPOLATE = 2 2 2
SELF_OFFDIAG = 1
MAXITER = 1000
MIXALPHA = 0.1
TOL_SCPH = 1e-11
WARMSTART = 0
RESTART_SCPH = 0
BUBBLE = 0
RELAX_STR = 0
/
&cell是两原子原胞，晶格缩放单位Bohr；前置力常数来自64原子超胞。这里仍使用同一个a=5.43Å的Si模型。MASS的单位是原子质量单位。
| 参数 | 本例含义 |
|---|---|
MODE=SCPH |
计算自洽声子 |
TMIN/TMAX/DT |
温度范围及间隔，单位K |
KMESH_SCPH |
自洽loop积分所用q网格 |
KMESH_INTERPOLATE |
重整化动力学矩阵的插值网格，本例固定2³ |
SELF_OFFDIAG=1 |
包含本例自能矩阵的非对角部分 |
TOL_SCPH=1e-11 |
自洽迭代的数值阈值 |
MIXALPHA |
混合参数，主要影响迭代过程 |
BUBBLE=0 |
本例四阶loop流程的设置 |
RELAX_STR=0 |
固定晶胞，不进行结构／晶胞热弛豫 |
3. 手动执行时怎么做#
在另一份新解压目录中，可分步操作：
module load alamode/1.5.0
python prepare.py
cd cases/quantum-q4
anphon input.in > anphon.log
grep 'convergence achieved' anphon.log
ls scph.scph_bands scph.scph_dymat scph.scph_dfc2
每个温度都应出现 convergence achieved。程序退出并不能单独证明每个温度都已自洽收敛。若日志存在 not converged，先检查迭代设置和输入，不要直接将对应温度曲线用于结果解释。
完整 run.sh已经执行了这些步骤。单独运行一个条件后可先看该目录输出；若要按步骤完成整批，在根目录依次执行 python prepare.py、python run.py、python effective_fc2.py、python analyze.py。最后的汇总图需要全部条件及有效FC2结果，不必在已完成的目录重复整套操作。
4. 读取实际温度频移#
下表为4³积分网格量子统计下的Γ点光学频率，取三个简并光学模的平均。谐波基线约为595.0767cm⁻¹。
| 温度 / K | Γ光学频率 / cm⁻¹ | 相对谐波频移 / cm⁻¹ |
|---|---|---|
| 0 | 594.0933 | -0.9834 |
| 300 | 592.4350 | -2.6417 |
| 600 | 589.8164 | -5.2603 |
0K的量子SCPH仍含零点振动贡献，因此可以与谐波基线不同。CLASSICAL=1采用经典占据，本例另外计算300、600K作对照。
这些频率来自固定体积的SW模型及本例四阶近似，不能直接替代硅的收敛DFT温度频移。晶胞热膨胀、三阶bubble等需要各自的计算设置。
5. 按温度分组再画图#
SCPH色散的前两列分别是温度和路径横坐标，后面才是六条频率。不要把第一列温度误当成路径横坐标。
head cases/quantum-q4/scph.scph_bands
python - <<'PY'
import numpy as np
a = np.loadtxt('cases/quantum-q4/scph.scph_bands')
b = a[a[:, 0] == 300]
np.savetxt('bands-300K.dat', b[:, 1:])
print('300K点数、路径与六支频率:', b.shape)
PY
回到根目录可按附件脚本重画全部结果：
python analyze.py
ls figures
本例比较2³／4³积分网格，以及混合参数0.1／0.2。自洽迭代收敛与q网格收敛是两件事；更精细研究还要继续检查积分和插值网格。有效FC2的导出见下一篇。
ALAMODE 1.5.0采用Intel oneAPI 2024.2、Intel MPI与MKL，当前CPU构建包含AVX-512；Python操作使用现有ASE环境。ANPHON输入说明