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

CP2K约束DFT与片段电荷

用水二聚体实际设置三个片段电子数，说明价电子目标、原子分组与完整输入。

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

想给体系中的某个片段限定电子数，可以使用约束DFT（CDFT）。本例固定水二聚体几何，在整个体系保持电中性的条件下，限定第一个水分子的Becke加权价电子数，让另一个片段通过自洽电子密度作出响应。

下载水二聚体CDFT输入、真实输出、分析脚本与数据图，解压进入run.sh目录。计算环境为CP2K2026.2 P100双精度，GNU13.3、MPICH5.0.1、CUDA12.9、OpenBLAS0.3.33、ScaLAPACK2.2.3；分析使用ASE3.29.0。

1 先做普通自洽，再做电荷约束#

sbatch submit.slurm

模板申请核时P100分区的一张GPU、4个CPU核。每个条件为一个MPI进程、4个OpenMP线程，先运行dimer-unconstrained得到自己这次计算的cdft-RESTART.wfn，再依次运行三个电荷目标，最后自动提取、画图和检查结果。

套餐用户用yskit生成所属分区脚本，在已分配的资源中执行bash run.sh。不要在登录节点运行计算。脚本只使用已分配CUDA_VISIBLE_DEVICES和CPU亲和范围；多卡时每卡一个任务，每任务4核，先完成共同的普通自洽，再分配三个约束目标。

inputs保存新输入，reference保存本文实算原始输入和输出，work保存新计算。已有work请重新解压，避免混用旧结果。下载包不含波函数；脚本现场生成前置文件。参考运行使用原子初猜，公开脚本进一步验证从新生成波函数继续计算，两种初猜均核对最终能量与约束值。

只算一个目标时，附上目录名；脚本仍自动完成其普通自洽前置步骤：

sbatch submit.slurm dimer-Ne8p1

2 TARGET填什么：电子数与净电荷#

本例ATOMS 1 2 3选中XYZ中的第一个水分子，COEFF 1 1 1表示这三个原子的Becke权重相加。XYZ编号从1开始，与Python数组从0开始不同。

使用GTH-PBE-q6氧和GTH-PBE-q1氢赝势，一个H2O片段的价电子基准为6+1+1=8。这里TARGET是权重积分得到的价电子数，并不是要填的净电荷值，也不用额外加上氧的芯电子。

因此本例TARGET 8.0对应片段净电荷约0e，8.1对应约−0.1e，8.2对应约−0.2e；片段电荷按q=8−N计算。ATOMIC_CHARGES TRUE请求打印相应电荷，不改变TARGET的定义。整个二聚体总电子数仍为16，没有给整个体系额外添加电子。

这是选定Becke划分下的片段电子数，不是唯一的原子电荷，也不等同于形式氧化数。更换权重、网格或基组后，应重新检查电荷与能量。

实际固定坐标如下；两个水分子依次为O、H、H，O—O距离2.9Å，盒长16Å，非周期MT边界：

6

O       8.000000000000000      8.000000000000000      6.000000000000000
H       8.000000000000000      8.000000000000000      6.970101302545075
H       7.999999851221545      7.059423260433508      5.762488033741113
O       8.000000000000000      8.000000000000000      8.900000000000000
H       8.000000003400000      8.765345449000000      9.496106434200000
H       7.999999851999999      7.234654621300000      9.496106524400000

3 普通自洽完整输入#

&GLOBAL
 PROJECT cdft
 RUN_TYPE ENERGY
 PRINT_LEVEL MEDIUM
&END GLOBAL
&FORCE_EVAL
 METHOD QS
 &DFT
  BASIS_SET_FILE_NAME ../DATA/BASIS_MOLOPT
  POTENTIAL_FILE_NAME ../DATA/GTH_POTENTIALS
  &QS
   EPS_DEFAULT 1.0E-12
  &END QS
  &MGRID
   CUTOFF 400
   REL_CUTOFF 50
  &END MGRID
  &SCF
   EPS_SCF 1.0E-6
   MAX_SCF 100
   SCF_GUESS ATOMIC
   &OT
    MINIMIZER DIIS
    PRECONDITIONER FULL_SINGLE_INVERSE
   &END OT
  &END SCF
  &XC
   &XC_FUNCTIONAL PBE
   &END XC_FUNCTIONAL
  &END XC
  &POISSON
   PERIODIC NONE
   POISSON_SOLVER MT
  &END POISSON
 &END DFT
 &SUBSYS
  &CELL
   ABC 16 16 16
   PERIODIC NONE
  &END CELL
  &TOPOLOGY
   COORD_FILE_NAME geometry.xyz
   COORD_FILE_FORMAT XYZ
  &END TOPOLOGY
  &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

首先固定几何，用PBE、DZVP-MOLOPT-SR-GTH和400/50Ry网格完成普通自洽。文件中的PROJECT cdft决定输出波函数名cdft-RESTART.wfn。

4 CDFT完整输入：8.1个价电子#

&GLOBAL
 PROJECT cdft
 RUN_TYPE ENERGY
 PRINT_LEVEL MEDIUM
&END GLOBAL
&FORCE_EVAL
 METHOD QS
 &DFT
  WFN_RESTART_FILE_NAME ../dimer-unconstrained/cdft-RESTART.wfn
  BASIS_SET_FILE_NAME ../DATA/BASIS_MOLOPT
  POTENTIAL_FILE_NAME ../DATA/GTH_POTENTIALS
  &QS
   &CDFT
    TYPE_OF_CONSTRAINT BECKE
    STRENGTH 0.0
    TARGET 8.1
    ATOMIC_CHARGES TRUE
    &ATOM_GROUP
     ATOMS 1 2 3
     COEFF 1 1 1
     CONSTRAINT_TYPE CHARGE
    &END ATOM_GROUP
    &OUTER_SCF ON
     TYPE CDFT_CONSTRAINT
     OPTIMIZER BISECT
     EPS_SCF 1.0E-5
     MAX_SCF 50
     STEP_SIZE -0.005
     BISECT_TRUST_COUNT 8
    &END OUTER_SCF
    &BECKE_CONSTRAINT
     CUTOFF_TYPE GLOBAL
     GLOBAL_CUTOFF 6.0
     ADJUST_SIZE FALSE
    &END BECKE_CONSTRAINT
    &PROGRAM_RUN_INFO ON
     &EACH
      QS_SCF 1
     &END EACH
     COMMON_ITERATION_LEVELS 2
     ADD_LAST NUMERIC
    &END PROGRAM_RUN_INFO
   &END CDFT
   EPS_DEFAULT 1.0E-12
  &END QS
  &MGRID
   CUTOFF 400
   REL_CUTOFF 50
  &END MGRID
  &SCF
   &OUTER_SCF ON
    EPS_SCF 1.0E-6
    MAX_SCF 10
   &END OUTER_SCF
   EPS_SCF 1.0E-6
   MAX_SCF 100
   SCF_GUESS RESTART
   &OT
    MINIMIZER DIIS
    PRECONDITIONER FULL_SINGLE_INVERSE
   &END OT
  &END SCF
  &XC
   &XC_FUNCTIONAL PBE
   &END XC_FUNCTIONAL
  &END XC
  &POISSON
   PERIODIC NONE
   POISSON_SOLVER MT
  &END POISSON
 &END DFT
 &SUBSYS
  &CELL
   ABC 16 16 16
   PERIODIC NONE
  &END CELL
  &TOPOLOGY
   COORD_FILE_NAME geometry.xyz
   COORD_FILE_FORMAT XYZ
  &END TOPOLOGY
  &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

WFN_RESTART_FILE_NAME指向同一work/cases下普通自洽目录生成的波函数，SCF_GUESS RESTART选择读取它。两处要同时写清；只写RESTART而前置文件路径不对，可能退回其他初猜。

QS/CDFT/OUTER_SCF控制约束乘子的优化，TYPE CDFT_CONSTRAINT和BISECT用于本例的一维约束；DFT/SCF/OUTER_SCF则服务于电子自洽。这两个OUTER_SCF处于不同层级，不能只留其中一处。电子EPS_SCF=1E-6，约束EPS_SCF=1E-5是不同的收敛量。

STRENGTH是约束乘子的初值，后续由优化器更新，不是给体系固定添加的电荷。本例8.0目录初值为0.5，8.1与8.2为0；不同目标的初猜可以不同，最终都要独立检查约束残差。参数不能保证适合任意分子、溶液或电极。

5 本例实际结果#

目录 目标价电子数 实际价电子数 片段净电荷/e 相对普通自洽能量/eV
dimer-Ne8p0 8.0 7.999995274 0.000004726 1.396342
dimer-Ne8p1 8.1 8.100009062 -0.100009062 0.244268
dimer-Ne8p2 8.2 8.199995467 -0.199995467 0.137525

实际固定水二聚体结构，以及三个约束目标对应的最终能量差

普通自洽总能约−34.446368900Ha。这里的能量差是同一固定几何、相同PBE和基组下约束态与普通自洽态的电子总能差，不是电子转移自由能、反应能垒或实验氧化还原电势。

提取与重画#

module load ase/3.29.0
python analyze.py work
python plot_results.py work
python check_results.py work
cat work/analysis/charge-energy.csv
ls work/figures

CSV记录目标电子数、实际电子数、片段净电荷、残差、最终总能和相对普通自洽的能量差。JSON还保存CDFT外循环历史。看本文参考数据时，把work改成reference；比较图需要完整三个目标。

参数定义见CP2K CDFT官方输入说明。接着看约束收敛与能量读取。

需要根据你的体系选择资源？联系讨论计算需求。实际可用软件与权限，以当前账户和运行环境为准。

查看配图

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

100%查看原图 ↗