VASP计算静电势：从LOCPOT到功函数与电荷分布分析

静电势（Electrostatic Potential）是VASP中描述空间中某一点电势能分布的基本物理量。它不只是一个中间计算结果，而是分析功函数、电荷转移、界面能带排列和局部电场的核心数据。VASP通过LOCPOT文件输出实空间静电势，正确读取和后处理这一数据是本篇的重点。

VASP中静电势的物理定义

在VASP的PAW框架下，总静电势V_eff(r)由三部分构成：

V_eff(r) = V_ion(r) + V_H(r) + V_xc(r)

• V_ion(r)：离子赝势贡献，描述原子核和内层电子的静电势
• V_H(r)：Hartree势，来自价电子的库仑排斥
• V_xc(r)：交换关联势，来自DFT泛函

三者求和后存储在LOCPOT文件中，单位是eV。VASP默认不输出LOCPOT，需要手动设置LVTOT=.TRUE.开启。

需要注意的是，LOCPOT中的静电势包含交换关联势V_xc，严格来说应该叫”总有效势”而非纯静电势。如果只需要纯库仑部分（V_ion + V_H），可以设置LVHAR=.TRUE.，输出到LOCPOT中的是V_ion + V_H，不含V_xc。这个区别在计算功函数时至关重要。

INCAR参数配置

输出完整静电势

LVTOT = .TRUE.
LCHARG = .TRUE.
PREC = Accurate

LVTOT=.TRUE.要求VASP在CHGCAR的基础上额外输出LOCPOT。LCHARG=.TRUE.确保电荷密度也同时输出（默认就是开启的）。PREC=Accurate提升格点精度，避免静电势在原子核附近出现数值噪声。

只输出Hartree+离子势（用于功函数）

LVHAR = .TRUE.
LVTOT = .FALSE.
PREC = Accurate

LVHAR=.TRUE.时，LOCPOT中存储的是V_ion + V_H，去掉了V_xc。这是计算功函数的标准做法，因为功函数的定义是电子从费米能级移动到真空能级的能量差，只涉及静电势，不涉及交换关联。

输出局部势能投影

LOCPROB = .TRUE.
LVTOT = .TRUE.

LOCPROB输出沿特定方向（x/y/z）的势能线积分，适合分析界面势能分布。

功函数计算：完整的实操流程

功函数Φ = V_vacuum – E_Fermi，其中V_vacuum是真空层电势，E_Fermi是费米能级。

第一步：构建slab模型

以金属Cu(111)表面为例：

• 在优化好的bulk Cu基础上切出(111)面
• 构建6-12层slab，底部3-4层固定
• 上方留出15-20Å真空层（必须足够大，否则真空电势不收敛）
• 建议加偶极修正（IDIPOL=3, LDIPOL=.TRUE.）消除slab两侧不对称导致的电场

第二步：自洽计算

PREC = Accurate
EDIFF = 1E-6
ENCUT = 450
ISMEAR = 0
SIGMA = 0.05
LVHAR = .TRUE.
LDIPOL = .TRUE.
IDIPOL = 3
NSW = 0

LDIPOL+IDIPOL=3做z方向偶极修正，这对有极性表面（如氧化物表面）尤其重要。不加偶极修正时，slab两侧电势不对称会导致真空电势偏移0.5-2 eV。

第三步：提取数据

计算完成后从OUTCAR中读费米能级：

grep "E-fermi" OUTCAR

LOCPOT中的势能数据需要用后处理工具提取z方向平均值：

• VASPKIT：vaspkit -task 411，直接输出planar-averaged electrostatic potential沿z方向的曲线
• VESTA：打开LOCPOT，沿z方向做线剖面
• 自写脚本：读取LOCPOT的3D网格数据，沿x-y平面平均后对z作图

第四步：确定真空电势

在z方向势能曲线的真空层区域，电势会趋近一个平台值。取平台中部的平均值作为V_vacuum。如果slab两侧都有真空层，应该取两侧平台值的平均（无偶极修正时两侧不等）。

实际案例：Cu(111)功函数

• 9层Cu(111) slab + 18Å真空层 + 偶极修正
• 费米能级E_F = -5.18 eV
• 真空电势V_vac = -1.32 eV
• 功函数Φ = V_vac – E_F = -1.32 – (-5.18) = 4.86 eV
• 实验值：4.94 eV（Cu多晶），误差1.6%

界面静电势与能带排列

静电势在异质界面研究中用于确定能带排列。以金属/半导体界面为例，关键步骤是分别计算金属和半导体bulk的”宏观平均静电势”V_macro，作为公共参考点：

1. 计算金属bulk的LOCPOT，做面平均得到V_metal(z)
2. 计算半导体bulk的LOCPOT，做面平均得到V_semi(z)
3. 构建界面slab，计算界面LOCPOT
4. 在界面slab的金属侧提取V_metal_offset，半导体侧提取V_semi_offset
5. 用两端的势能偏移校准能带排列

这个方法的核心优势是不依赖绝对能量（DFT中绝对能量无物理意义），而是用静电势的相对偏移确定能带对齐。

常见问题诊断

LOCPOT全为零：检查是否设置了LVTOT=.TRUE.或LVHAR=.TRUE.。如果做了非自洽计算（ICHARG=11），LOCPOT中不会包含新的势能，需要先做自洽计算。

真空电势不收敛（没有平台）：真空层太薄（<15Å），或者k网格在面内不够密（至少12×12才能保证势能在面内充分平均）。增大真空层和k网格密度可以解决。

功函数随slab层数漂移：金属slab的功函数在层数<6层时会显著偏移，原因是量子限制效应。建议至少用8-10层，并检查收敛性。