FourPhonon四声子热导入门

下载Si真实结构、受力、力常数、8张图和终端脚本，保存为 si-fourphonon.zip。所有本系列文章使用同一个包：

python3 -m zipfile -e si-fourphonon.zip .
cd fourthorder-fourphonon-silicon-wave288
bash preview.sh
ls preview/figures

包中 reference/ 保存本次实算参考；work/ 是您重新生成的文件，results/data/ 是提取后的CSV。受力来自真实VASP计算；下载包不含POTCAR、波函数、完整电子输出或第三方程序。

1 要准备的四个文件#

CONTROL描述原胞、元素、超胞、q网格与求解选项；其余三个文件为 FORCE_CONSTANTS_2ND、FORCE_CONSTANTS_3RD、FORCE_CONSTANTS_4TH。本例来源见Fourthorder重构，不是使用程序自带的材料常数。

CPU模块为 fourphonon/1.3-fix20260924-8173m，固定源码提交5b814b3，oneAPI 2024.2 / Intel MPI 2021.13 / 顺序MKL，需支持AVX512的CPU。版本后缀为集群构建标识。

2 完整CONTROL#

&allocations
 nelements=1, natoms=2, ngrid(:)=3 3 3
&end
&crystal
 lfactor=0.1
 lattvec(:,1)=0 2.715 2.715
 lattvec(:,2)=2.715 0 2.715
 lattvec(:,3)=2.715 2.715 0
 elements="Si", types=1 1
 positions(:,1)=0 0 0
 positions(:,2)=0.25 0.25 0.25
 masses(1)=28.085
 gfactors(1)=0.0
 scell(:)=2 2 2
&end
&parameters
 T=300, scalebroad=1, maxiter=200, eps=1e-6
 num_sample_process_3ph=-1
 num_sample_process_3ph_phase_space=-1
 num_sample_process_4ph=-1
 num_sample_process_4ph_phase_space=-1
&end
&flags
 nonanalytic=.false., espresso=.false., tdep=.false.
 convergence=.false.
 four_phonon=.true.
 four_phonon_iteration=.false.
 isotopes=.false., autoisotopes=.false.
 nanowires=.false., onlyharmonic=.false.
&end

lfactor=0.1把这里以Å填写的晶格长度换成程序所用nm；positions为分数坐标。scell=2 2 2与实际常数的超胞一致。masses=28.085与本例谐波DFT质量一致。

ngrid是声子q网格。four_phonon=.true.启用四阶散射；convergence=.false.在此表示采用RTA，不表示“不检查数值收敛”。isotopes=.false.时同位素g因子不参与求解，本入门例显式写0。

3 执行配套工作流#

bash prepare.sh
bash reconstruct.sh
sbatch cpu.sh
squeue -u "$USER"
# 作业完成后
cat results/data/transport.csv
ls plots/figures

该脚本申请CPU资源后依次运行15个原生小条件，先小网格，再较大网格，不重算DFT。已经完成prepare/reconstruct时直接提交即可。它适用于支持AVX512的Intel CPU队列；队列可用性和资源以门户当前显示为准。

单独计算一个目录时，原生命令如下，在已分配的CPU资源中运行：

module load fourphonon/1.3-fix20260924-8173m
export OMP_NUM_THREADS=8 MKL_NUM_THREADS=1 OPENBLAS_NUM_THREADS=1
cd work/transport/CPU-q3-T300-3plus4ph
mpirun -np 1 fourphonon
cat T300K/BTE.kappa_tensor

配套批处理已经计算过的目录不要再次覆盖；另试参数时复制到新目录，保留原始CONTROL。

4 本例读到什么#

300 K、3³ q网格、RTA下：仅三声子为 26.17719 W/(m·K)，三加四声子为 24.06041 W/(m·K)，降低约8.09%。这是16原子/最近邻教学模型的粗网格结果，不能作为收敛的Si热导预测。

BTE.kappa_tensor第一列是迭代编号，后9列是热导张量，单位W/(m·K)，按Fortran列顺序排列。RTA例只有编号0的一行；本Si模型三个对角元相同、非对角元接近数值零。不要把第一列0误当成温度。

主要输出：

进一步解释频率、散射率与寿命见结果提取。