鸢算科技 yuansuan Technology
用户门户
CP2K专题

CP2K ELF、PDOS与电荷分区

用水分子真实网格读取电子密度和ELF，画原子kind投影谱，并比较两种Hirshfeld分区。

更新于 2026-11-08约 9 分钟阅读按任务与环境选择适用配置
阅读本篇目录

用水分子真实网格读取电子密度和ELF，画原子kind投影谱，并比较两种Hirshfeld分区。

1. 下载实例，先提取再画图#

下载电子性质与局域轨道实操包，在新目录解压：

unzip cp2k-electronic-fields-wave294-*.zip
cd cp2k-electronic-fields-wave294
bash preview.sh
bash run.sh

preview/重画三张参考图，results/data/保存重新提取的数值，plots/保存重画结果；随包figures/可直接查看。输出目录受保护，再次运行请使用新解压目录。这里不重新运行DFT，只在已分配资源内使用1个CPU后处理。

原生计算用CP2K2026.2的CPU构建：Intel oneAPI2024.2、icx/icpx/ifx2024.2.1、Intel MPI2021.13、MKL2024.2，MPI＋OpenMP。本例1个MPI进程、4个OpenMP线程。后处理用pymatgen/2026.9.24中的Python3.11、NumPy2.3.5、Matplotlib3.11.2。

reference/是实际CUBE及PDOS导出的自有数值，选定网格完整保留，没有抽样或人为归一化；不含程序、原始二进制WFN或完整运行日志。inputs/保留完整输入与原生PDOS数值，方便核对字段。

2. 请求密度、ELF和投影DOS#

在DFT的PRINT段加入下面内容：

&E_DENSITY_CUBE ON
 STRIDE 2 2 2
&END E_DENSITY_CUBE
&ELF_CUBE ON
 STRIDE 2 2 2
&END ELF_CUBE
&DOS ON
 NLUMO 4
 NDIGITS 10
 &PDOS
  COMPONENTS T
 &END PDOS
&END DOS
&HIRSHFELD ON
 SHAPE_FUNCTION DENSITY
&END HIRSHFELD

2026.2把投影PDOS放在DOS/PDOS内：NLUMO属于外层DOS，COMPONENTS属于内层PDOS。不能把旧版的整个&PDOS直接放到DFT/PRINT下。当前版本层级见CP2K DOS及PDOS。

本例4个占据空间轨道，再请求4个空轨道；不是4个电子。对角化时由SCF/ADDED_MOS提供空轨道。本例使用ADDED_MOS 4，每个占据轨道占据数2，总价电子数8。

3. 读取CUBE的单位和网格#

head -8 water-ELECTRON_DENSITY-1_0.cube
head -8 water-ELF_S1-1_0.cube

本例CUBE为75³输出网格，步长约0.302356 Bohr，即0.1600 Å。STRIDE 2是输出取样间隔，不是把DFT截断减半。电子密度单位e/Bohr³，积分要乘网格体积；若坐标改成Å，密度单位也要一起转换。

本例电子密度积分为 8.00154882 e，接近GTH价电子数8。它不含相同表示下的全部芯电子，不能与全电子水分子的10电子直接混为一谈。小偏差保留原值，没有将整幅数据强行缩放到8。

水分子实际价电子密度与ELF

ELF是0到1的无量纲指标，不是电子密度，也不是直接的成键概率。图中两幅色标含义不同：左边是密度的log10，右边是ELF；采样面实际在x=6.0800 Å，原子位于x=6 Å。低密度真空有数值阈值处理，不应凭真空颜色判断电子成键。参阅官方ELF输出说明。

4. 从原生PDOS文件画出投影谱#

在资料包根目录：

head -8 inputs/water-BERRY-restart/water-k1-1.pdos
head -8 inputs/water-BERRY-restart/water-k2-1.pdos
head -5 results/data/BERRY-PDOS.csv

k1、k2在本例分别是O与H的原子kind，不是k点。H kind包含两个氢原子的投影；不能未经处理就当成单个H的投影。列头依次列出轨道编号、能量、占据、Total以及s/p/d分量，能量单位为a.u.，换成eV乘27.211386。

脚本把能量零点放在HOMO，并用0.15 eV高斯展宽画谱；每个空间轨道权重为1，两个kind的权重之和为1，再乘各自占据才能得到电子数。本例HOMO约−7.2004 eV；这些是Kohn–Sham轨道，不将空轨道能量或轨道差直接称为实验电子亲和能、光学激发能。

实际投影DOS、两种电荷分区与Berry中心

5. 电荷分区先说明参考形状#

默认SHAPE_FUNCTION GAUSSIAN与DENSITY使用不同参考形状，不是同一个结果换了显示方式。后者采用多个高斯展开的原子密度；可选设置见CP2K Hirshfeld。本例分别真实计算得到：

参考形状 O/e H1/e H2/e 原生打印Total/e
默认Gaussian −1.107 +0.559 +0.559 +0.011
Atomic-density −0.305 +0.152 +0.152 −0.000

这是分区电荷，不是原子的整数氧化态。默认Gaussian在本例保留了约0.011 e残差；不要把残差人为摊到各原子后冒充原始结果。逐原子数值只打印三位小数，总和列也独立打印，不能拿三位小数证明更高精度的电荷守恒。

cat results/data/water-Hirshfeld-charges.csv

更换体系时，保持方法、基组、赝势、参考形状一致，再比较相对变化；对带电、金属和强离域体系还要检查分区及网格敏感性。

准备自己的计算并提交#

本例水分子为PBE、DZVP-MOLOPT-SR-GTH/GTH-PBE、12 Å非周期立方盒、MT泊松边界；密度截断400 Ry、REL_CUTOFF=50 Ry、EPS_SCF=1E−6。先优化自己的结构，再请求单点性质；不要把非周期分子的MT设置直接套到周期晶体。

mkdir water-fields
cp inputs/water-BERRY-restart/input.inp water-fields/
cd water-fields

输入中的SCF_GUESS RESTART与WFN_RESTART_FILE_NAME seed.wfn用于读取同一结构、基组和赝势的已收敛波函数。将自己的匹配文件复制为seed.wfn。没有前置WFN时，可在这个新目录改成原子初猜：

sed -i 's/SCF_GUESS RESTART/SCF_GUESS ATOMIC/' input.inp
sed -i '/WFN_RESTART_FILE_NAME/d' input.inp

通过yskit选择CP2K2026.2对应的CPU模板，也可以复制公共模板；CPU分区须处于可用状态。MPI与线程数按模板分配，不混入其他MPI环境：

cp /public/slurm_scripts/hs_user/cp2k/2026.2/8173m/sub_cp2k_2026.2_8173m.sh .
sbatch sub_cp2k_2026.2_8173m.sh
grep 'SCF run converged' cp2k.out
tail -12 cp2k.out
ls *.cube *.pdos

模板的并行配置与本例1×4不同，输入模型不变；小分子不一定随核心数增加而加速。详细提交流程见CP2K入门。

完整输入#

以下原始输入也保留在资料包inputs/中。

water-BERRY-restart 
&GLOBAL
 PROJECT water
 RUN_TYPE ENERGY_FORCE
 PRINT_LEVEL MEDIUM
&END GLOBAL
&FORCE_EVAL
 METHOD QS
 &DFT
  BASIS_SET_FILE_NAME BASIS_MOLOPT
  POTENTIAL_FILE_NAME GTH_POTENTIALS
  WFN_RESTART_FILE_NAME seed.wfn
  &MGRID
   CUTOFF 400
   REL_CUTOFF 50
  &END MGRID
  &SCF
   EPS_SCF 1.0E-6
   MAX_SCF 120
   ADDED_MOS 4
   SCF_GUESS RESTART
   &DIAGONALIZATION
    ALGORITHM STANDARD
   &END DIAGONALIZATION
   &MIXING
    METHOD PULAY_MIXING
    ALPHA 0.3
    NBUFFER 5
   &END MIXING
  &END SCF
  &XC
   &XC_FUNCTIONAL PBE
   &END XC_FUNCTIONAL
  &END XC
  &POISSON
   PERIODIC NONE
   POISSON_SOLVER MT
  &END POISSON
  &LOCALIZE T
   METHOD JACOBI
   OPERATOR BERRY
   EPS_LOCALIZATION 1E-8
   MAX_ITER 5000
   &PRINT
    &WANNIER_CENTERS ON
     IONS+CENTERS T
    &END WANNIER_CENTERS
    &WANNIER_SPREADS ON
    &END WANNIER_SPREADS
    &WANNIER_CUBES ON
     STRIDE 2 2 2
    &END WANNIER_CUBES
    &TOTAL_DIPOLE ON
    &END TOTAL_DIPOLE
   &END PRINT
  &END LOCALIZE
  &PRINT
   &E_DENSITY_CUBE ON
    STRIDE 2 2 2
   &END E_DENSITY_CUBE
   &ELF_CUBE ON
    STRIDE 2 2 2
   &END ELF_CUBE
   &DOS ON
    NLUMO 4
    NDIGITS 10
    &PDOS
     COMPONENTS T
    &END PDOS
   &END DOS
   &HIRSHFELD ON
   &END HIRSHFELD
   &MULLIKEN ON
   &END MULLIKEN
   &MOMENTS ON
    PERIODIC F
    MAX_MOMENT 1
   &END MOMENTS
  &END PRINT
 &END DFT
 &SUBSYS
  &CELL
   ABC 12 12 12
   PERIODIC NONE
  &END CELL
  &COORD
   O 6.000 6.000 6.000
   H 6.000 6.757 6.587
   H 6.000 5.243 6.587
  &END COORD
  &KIND O
   BASIS_SET DZVP-MOLOPT-SR-GTH
   POTENTIAL GTH-PBE-q6
  &END KIND
  &KIND H
   BASIS_SET DZVP-MOLOPT-SR-GTH
   POTENTIAL GTH-PBE-q1
  &END KIND
 &END SUBSYS
 &PRINT
  &FORCES ON
  &END FORCES
 &END PRINT
&END FORCE_EVAL
water-Hirshfeld-density 
&GLOBAL
 PROJECT water
 RUN_TYPE ENERGY_FORCE
 PRINT_LEVEL MEDIUM
&END GLOBAL
&FORCE_EVAL
 METHOD QS
 &DFT
  BASIS_SET_FILE_NAME BASIS_MOLOPT
  POTENTIAL_FILE_NAME GTH_POTENTIALS
  WFN_RESTART_FILE_NAME seed.wfn
  &MGRID
   CUTOFF 400
   REL_CUTOFF 50
  &END MGRID
  &SCF
   EPS_SCF 1.0E-6
   MAX_SCF 120
   ADDED_MOS 4
   SCF_GUESS RESTART
   &DIAGONALIZATION
    ALGORITHM STANDARD
   &END DIAGONALIZATION
   &MIXING
    METHOD PULAY_MIXING
    ALPHA 0.3
    NBUFFER 5
   &END MIXING
  &END SCF
  &XC
   &XC_FUNCTIONAL PBE
   &END XC_FUNCTIONAL
  &END XC
  &POISSON
   PERIODIC NONE
   POISSON_SOLVER MT
  &END POISSON
  &PRINT
   &HIRSHFELD ON
    SHAPE_FUNCTION DENSITY
   &END HIRSHFELD
   &MULLIKEN ON
   &END MULLIKEN
   &MOMENTS ON
    PERIODIC F
    MAX_MOMENT 1
   &END MOMENTS
  &END PRINT
 &END DFT
 &SUBSYS
  &CELL
   ABC 12 12 12
   PERIODIC NONE
  &END CELL
  &COORD
   O 6.000 6.000 6.000
   H 6.000 6.757 6.587
   H 6.000 5.243 6.587
  &END COORD
  &KIND O
   BASIS_SET DZVP-MOLOPT-SR-GTH
   POTENTIAL GTH-PBE-q6
  &END KIND
  &KIND H
   BASIS_SET DZVP-MOLOPT-SR-GTH
   POTENTIAL GTH-PBE-q1
  &END KIND
 &END SUBSYS
 &PRINT
  &FORCES ON
  &END FORCES
 &END PRINT
&END FORCE_EVAL
需要根据你的体系选择资源？联系讨论计算需求。实际可用软件与权限，以当前账户和运行环境为准。

查看配图

放大后可滚动查看细节，按 Esc 关闭。

100%查看原图 ↗