VASP原生电子–声子矩阵元

电子–声子矩阵元描述原子振动怎样耦合电子态。它与声子DOS、形变势常数或单纯能带曲线不同。本例用VASP原生功能，先计算16原子Si超胞的势导数，再在两原子原胞中输出指定电子态之间的耦合。

1. 先计算超胞势导数#

进入附件si16-elph-potential目录，按POSCAR准备合法Si POTCAR，提交sbatch run.sh。PBE PAW64、ENCUT=550 eV、EDIFF=1E-6、4³Gamma网格，关键设置：

IBRION = 6
POTIM = 0.01
NFREE = 2
ELPH_POT_GENERATE = .TRUE.

程序自动生成有限位移并计算势对位移的响应。这是有限差分步骤，不是原子位置优化；位移结构的力不需要满足优化终止阈值。ELPH_POT_GENERATE

完成后应得到phelel_params.hdf5和CONTCAR_ELPH。后者给出与势导数匹配的原胞，后续不要自行换成方向或原子顺序不同的另一个原胞。

2. 将结果传到原胞计算#

在附件根目录执行：

python prepare_response.py

脚本把CONTCAR_ELPH复制为三个响应目录的POSCAR，同时复制phelel_params.hdf5。每个目录仍需准备同一套Si POTCAR。然后进入si-mels-fullbands，执行sbatch run.sh。

ELPH_RUN = .TRUE.
ELPH_DRIVER = MELS
ELPH_SELFEN_IKPT = 1
ELPH_SELFEN_BAND_START = 1
ELPH_SELFEN_BAND_STOP = 8
ELPH_NBANDS = 8

原胞SCF采用4³Gamma网格；本例输出中的第1个不可约k点已核对为Γ。KPOINTS_ELPH取Γ–X的11个点，原胞倒格坐标终点为(0.5,0,0.5)。不同结构不能直接假定索引1必然是目标k点。

超胞势导数和原胞必须使用一致的泛函、赝势、截断能及兼容FFT网格。这里让VASP根据550 eV自动确定网格，并核对后续读取匹配。原生电子–声子流程

3. 提取矩阵元#

全部响应目录完成后，在附件根目录加载含numpy和h5py的环境，再运行python extract_native_elph.py。集群可使用下列现有环境：

module load py4vasp/0.11.3
python extract_native_elph.py

这里使用现有py4vasp环境中的numpy和h5py读取HDF5，电子–声子矩阵元由VASP计算。

vaspelph.h5保存复数矩阵元，实部与虚部在最后一个维度，单位eV。本例采用默认VDPR的PAW矩阵元定义。输出格式、矩阵元组成

4. 本例实际图像#

图中取Γ侧与路径侧的完整占据价带1–4，计算所有所选带对和声子分支的|g|²之和，除以带对数量后开平方。这样得到所选子空间的耦合范数，便于演示数据读取；它不是散射率，不含能量守恒和占据因子，也不是电子–声子λ。

为避免把接近零频平移模的数值噪声放大，排除|频率|<10⁻⁴ THz的模式，尤其Γ处三个声学平移模。保留模式中单个矩阵元的最大模为0.4472 eV，但单个简并态的矩阵元可随基底混合变化，不能把其大小当成独立可观测量。

当前输出的原始Γ光学频率94.690932与OUTCAR的2PiTHz列对应，除以2π得到15.070530 THz。提取时要核对文件实际单位，不能把同一数值直接标成THz或eV。

5. 下一步能得到什么#

计算散射率和电子自能还需要对中间态、声子和布里渊区积分，不能仅靠这条路径上的矩阵元推导迁移率。可继续阅读电子–声子自能与线宽教程。

本例采用vasp/6.6.1-p100-hdf5（NVHPC25.3、CUDA12.8、HDF5接口），通过Slurm完成。16原子超胞、4³SCF与8条耦合能带用于展示流程，尚未构成参数收敛结论。

下载输入、结果与提取脚本（ZIP）