CP2K专题

CP2K MP2与RPA解析力

实际RI-MP2和Minimax RI-RPA力计算，完整梯度输入、原生力向量与单位换算。

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

阅读本篇目录

相关能单点之后，可以继续学习当前版本支持的解析力。本篇固定水分子，实际计算RI-MP2和Minimax RI-RPA的原子力，再在下一篇用独立能量位移检查一个分量。

下载解析力、两组位移输入、实际日志与工具脚本，解压进入run.sh目录。CP2K2026.2 P100双精度：GNU13.3、MPICH5.0.1、CUDA12.9、OpenBLAS0.3.33、ScaLAPACK2.2.3；分析ASE3.29.0。

1 从输入生成到提交#

python3 make_displacements.py
sbatch submit.slurm

make_displacements.py从两个基准输入生成8个正负位移输入，只移动指定原子，不重新居中。run.sh也会自动调用它；模板默认完成10次运行：RI-MP2与Minimax RI-RPA各一次解析力，加上每种方法两种步长的正负位移单点。

模板为核时P100分区，一GPU一CPU，每个条件一个MPI与一个OpenMP线程。套餐用户用yskit生成所属分区脚本，在已分配资源中保留bash run.sh。多卡仅使用已经分配的可见GPU/CPU亲和范围，每卡一个任务，空闲卡接下一项；不自动申请资源。不要在登录节点计算，已有work重新解压。

inputs是新计算输入；reference是本文实际原始输出；work是新结果。只运行一个条件可在sbatch后附CASES.json中的名称，例如：

sbatch submit.slurm water-mp2-force

单条件只输出自己的数据；完整力图和差分对照需要整组。

2 模型与方法范围#

同一实际PBE优化水几何，放在14 Å非周期WAVELET盒；轨道基组DZVP-GTH，RI辅助基组RI_DZVP-GTH，300/50 Ry网格、EPS_SCF=1E-7。RI-MP2用HF/GTH-HF参考，RI-RPA用PBE/GTH-PBE参考；两个方法的绝对总能不混作准确度排序。

本次方法与原结构优化方法不同，新势能面上通常存在非零力。若要得到新方法的平衡结构，需要按该方法继续松弛并检查数值参数。

明确固定的坐标为：

3

O       7.000000072300001      6.999999964850000      6.701946737800000
H       7.000000075700000      7.765345413850000      7.298053172000000
H       6.999999924300000      6.234654586150000      7.298053262200000

3 RI-MP2解析力完整输入#

完整input.inp：water-mp2-force#

&GLOBAL
  PRINT_LEVEL MEDIUM
  PROJECT RI_MP2_H2O
  RUN_TYPE ENERGY_FORCE
  &TIMINGS
    THRESHOLD 0.01
  &END TIMINGS
&END GLOBAL

&FORCE_EVAL
  METHOD Quickstep
  &DFT
    BASIS_SET_FILE_NAME ../DATA/HFX_BASIS
    POTENTIAL_FILE_NAME ../DATA/POTENTIAL
    &MGRID
      CUTOFF 300
      REL_CUTOFF 50
    &END MGRID
    &POISSON
      PERIODIC NONE
      POISSON_SOLVER WAVELET
    &END POISSON
    &QS
      EPS_DEFAULT 1.0E-12
      EPS_PGF_ORB 1.0E-18
      METHOD GPW
    &END QS
    &SCF
      EPS_SCF 1.0E-7
      MAX_SCF 100
      SCF_GUESS ATOMIC
    &END SCF
    &XC
      &HF
        FRACTION 1.0000000
        &SCREENING
          EPS_SCHWARZ 1.0E-8
              EPS_SCHWARZ_FORCES 1.0E-8
          SCREEN_ON_INITIAL_P FALSE
        &END SCREENING
      &END HF
      &WF_CORRELATION
        &CANONICAL_GRADIENTS
          EPS_CANONICAL 1.0E-7
          FREE_HFX_BUFFER FALSE
          &CPHF
            EPS_CONV 1.0E-8
            MAX_ITER 50
          &END CPHF
        &END CANONICAL_GRADIENTS
        MEMORY 200.
        NUMBER_PROC 1
        &INTEGRALS
          &WFC_GPW
            CUTOFF 300
            EPS_FILTER 1.0E-8
            EPS_GRID 1.0E-8
            REL_CUTOFF 50
          &END WFC_GPW
        &END INTEGRALS
        &RI_MP2
        &END RI_MP2
      &END WF_CORRELATION
      &XC_FUNCTIONAL NONE
      &END XC_FUNCTIONAL
    &END XC
  &END DFT
  &SUBSYS
    &CELL
      ABC [angstrom] 14.0 14.0 14.0
      PERIODIC NONE
    &END CELL
    &KIND H
      BASIS_SET DZVP-GTH
      BASIS_SET RI_AUX RI_DZVP-GTH
      POTENTIAL GTH-HF-q1
    &END KIND
    &KIND O
      BASIS_SET DZVP-GTH
      BASIS_SET RI_AUX RI_DZVP-GTH
      POTENTIAL GTH-HF-q6
    &END KIND
    &TOPOLOGY
      COORD_FILE_FORMAT xyz
      COORD_FILE_NAME H2O_gas.xyz
    &END TOPOLOGY
  &END SUBSYS
  &PRINT
    &FORCES ON
    &END FORCES
  &END PRINT
&END FORCE_EVAL

GLOBAL用ENERGY_FORCE请求能量和力。CANONICAL_GRADIENTS控制相关方法的响应求解，CPHF/EPS_CONV控制相应迭代；它们与SCF/EPS_SCF不是同一个收敛量。当前例子固定这些值，先练习完整流程。

4 Minimax RI-RPA解析力完整输入#

完整input.inp：water-rpa-force#

&GLOBAL
  PRINT_LEVEL MEDIUM
  PROJECT RI_RPA_H2O
  RUN_TYPE ENERGY_FORCE
  &TIMINGS
    THRESHOLD 0.01
  &END TIMINGS
&END GLOBAL

&FORCE_EVAL
  METHOD Quickstep
  &DFT
    BASIS_SET_FILE_NAME ../DATA/HFX_BASIS
    POTENTIAL_FILE_NAME ../DATA/GTH_POTENTIALS
    &MGRID
      CUTOFF 300
      REL_CUTOFF 50
    &END MGRID
    &POISSON
      PERIODIC NONE
      POISSON_SOLVER WAVELET
    &END POISSON
    &QS
      EPS_DEFAULT 1.0E-12
      EPS_PGF_ORB 1.0E-18
      METHOD GPW
    &END QS
    &SCF
      EPS_SCF 1.0E-7
      MAX_SCF 100
      SCF_GUESS ATOMIC
      &PRINT
        &RESTART OFF
        &END RESTART
      &END PRINT
    &END SCF
    &XC
      &WF_CORRELATION
        &CANONICAL_GRADIENTS
          EPS_CANONICAL 1.0E-7
          FREE_HFX_BUFFER FALSE
          &CPHF
            EPS_CONV 1.0E-8
            MAX_ITER 50
          &END CPHF
        &END CANONICAL_GRADIENTS
        MEMORY 200.
        NUMBER_PROC 1
        &INTEGRALS
          ERI_METHOD GPW
          &WFC_GPW
            CUTOFF 300
            REL_CUTOFF 50
          &END WFC_GPW
        &END INTEGRALS
        &RI_RPA
          MINIMAX_QUADRATURE TRUE
          RPA_NUM_QUAD_POINTS 8
          &HF
            FRACTION 1.0000000
            &SCREENING
              EPS_SCHWARZ 1.0E-8
              EPS_SCHWARZ_FORCES 1.0E-8
              SCREEN_ON_INITIAL_P FALSE
            &END SCREENING
          &END HF
        &END RI_RPA
      &END WF_CORRELATION
      &XC_FUNCTIONAL PBE
      &END XC_FUNCTIONAL
    &END XC
  &END DFT
  &SUBSYS
    &CELL
      ABC [angstrom] 14.0 14.0 14.0
      PERIODIC NONE
    &END CELL
    &KIND H
      BASIS_SET DZVP-GTH
      BASIS_SET RI_AUX RI_DZVP-GTH
      POTENTIAL GTH-PBE-q1
    &END KIND
    &KIND O
      BASIS_SET DZVP-GTH
      BASIS_SET RI_AUX RI_DZVP-GTH
      POTENTIAL GTH-PBE-q6
    &END KIND
    &TOPOLOGY
      COORD_FILE_FORMAT xyz
      COORD_FILE_NAME H2O_gas.xyz
    &END TOPOLOGY
  &END SUBSYS
  &PRINT
    &FORCES ON
    &END FORCES
  &END PRINT
&END FORCE_EVAL

本版本RPA解析梯度需要Minimax积分。本例用MINIMAX_QUADRATURE TRUE和8个积分点；交换修正教程中的SOSEX、AXK单点输入不直接用作这里的力模板。

5 从原生输出提取力#

grep -E 'ENERGY\||FORCES\|' work/cases/water-mp2-force/output.log
grep -E 'ENERGY\||FORCES\|' work/cases/water-rpa-force/output.log

2026.2这里的输出标题为FORCES| Atomic forces [hartree/bohr]。每行依次是原子编号、Fx、Fy、Fz及模长；本XYZ按O、H1、H2排列，原生编号从1开始，所以H1为第2行。

脚本用ASE的Hartree/Bohr常数将原生力换成eV/Å；不要只乘Hartree到eV而漏掉长度单位。

方法	H1-z解析力/eV·Å⁻¹	总力向量模长/eV·Å⁻¹
RI-MP2	0.029697650	0.0009147
Minimax RI-RPA	0.381219741	0.0083786

当前固定结构原子力的y-z平面投影，箭头采用相同物理比例

当前有限参数下，总力仍有表中所示残余。结构优化等实际任务应继续检查网格、基组和边界精度；下面一个力分量的差分吻合，也不代表全部分量或完整基组已经收敛。

提取与重画#

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

CSV包含方法、位移步长、解析力、差分力和两者之差，力单位eV/Å；JSON还保存三个原子的完整力向量。看本文数据时将work改成reference。

梯度支持条件见MP2官方说明与RPA官方说明。继续看固定坐标的能量差分核对。

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