GPUMDkit专题

GPUMDkit 转换 CP2K 静态与 MD 标签

用8原子Si实际计算区分静态日志和MD三文件路线,匹配帧并补入逐帧virial。

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

下载原始输入输出、转换脚本与绘图资料包。解压进入目录,inputs 是原始参考输入,reference 是实际结果。bash run.sh all 只做转换与检查,不重新启动 DFT;每次在 work 新建目录,LAST-RUN.txt 记录位置。下面独立示例均从资料包根目录开始。

先选对两条路线#

GPUMDkit 1.5.8 的 -cp2k2xyz 短命令读取 静态日志与结构;MD 使用菜单里的 坐标、力、晶胞三文件路线。先区分输入类型,再开始转换。

本篇实际计算了 8 原子周期 Si:一个静态能量/力/应力点,以及 4 步、步长 0.5 fs、共 5 个输出构型的短 NVE。它们用于生成和检查标签,不用于估计扩散或热力学平均。

实际CP2K短轨迹的能量力和逐帧virial

1. 静态计算的完整输入#

静态例子的 input.inp:

&GLOBAL
 PROJECT labels
 RUN_TYPE ENERGY_FORCE
 PRINT_LEVEL MEDIUM
&END GLOBAL
&FORCE_EVAL
 METHOD QS
 STRESS_TENSOR ANALYTICAL
 &DFT
  BASIS_SET_FILE_NAME BASIS_MOLOPT
  POTENTIAL_FILE_NAME GTH_POTENTIALS
  &MGRID
   CUTOFF 400
   REL_CUTOFF 50
  &END MGRID
  &QS
   EPS_DEFAULT 1.0E-10
  &END QS
  &SCF
   EPS_SCF 1.0E-6
   MAX_SCF 250
   SCF_GUESS ATOMIC
   &OT
    MINIMIZER DIIS
    PRECONDITIONER FULL_SINGLE_INVERSE
   &END OT
  &END SCF
  &XC
   &XC_FUNCTIONAL PBE
   &END XC_FUNCTIONAL
  &END XC
 &END DFT
 &SUBSYS
  &CELL
   A 5.43 0 0
   B 0 5.43 0
   C 0 0 5.43
   PERIODIC XYZ
  &END CELL
  &COORD
   Si 0.0200000000 0.0000000000 0.0000000000
   Si 1.3575000000 1.3575000000 1.3575000000
   Si 0.0000000000 2.7150000000 2.7150000000
   Si 1.3575000000 4.0725000000 4.0725000000
   Si 2.7150000000 0.0000000000 2.7150000000
   Si 4.0725000000 1.3575000000 4.0725000000
   Si 2.7150000000 2.7150000000 0.0000000000
   Si 4.0725000000 4.0725000000 1.3575000000
  &END COORD
  &KIND Si
   BASIS_SET DZVP-MOLOPT-SR-GTH
   POTENTIAL GTH-PBE-q4
  &END KIND
 &END SUBSYS
 &PRINT
  &FORCES ON
  &END FORCES
  &STRESS_TENSOR ON
  &END STRESS_TENSOR
 &END PRINT
&END FORCE_EVAL

这是 PBE、DZVP-MOLOPT-SR-GTH/GTH-PBE-q4,CUTOFF=400 Ry、REL_CUTOFF=50 Ry、EPS_SCF=1E-6。第一个 Si 沿 x 位移 0.02 Å,因此有非零力。完整 A/B/C 晶胞与笛卡尔坐标也方便转换器直接对应。

本例使用 cp2k/2026.2-p100:GNU 13.3、MPICH 5.0.1、CUDA 12.9、OpenBLAS 0.3.33。要重算,先分配计算资源,在自己的输入副本中执行:

module purge
module load cp2k/2026.2-p100
export OMP_NUM_THREADS=2
mpiexec -bind-to core:2 -n 1 cp2k.psmp -i input.inp -o output.log

包内另有 compute.sbatch,在新目录重算四个小体系 DFT 条件;run.sbatch 只做已保存数据的后处理。二者用途不同。

2. 转换静态日志#

从资料包根目录准备副本:

module purge
module load gpumdkit/1.5.8
mkdir my-cp2k-scf
cp inputs/cp2k-scf/input.inp inputs/cp2k-scf/output.log my-cp2k-scf/
cd my-cp2k-scf
cp output.log cp2k.log
gpumdkit.sh -cp2k2xyz
head -n 4 cp2k_exyz.xyz

将计算日志命名为 cp2k.log,让工具明确选中它,避免目录里其他 .log 干扰。只放同一次静态计算的结构与日志;本例包括力和解析应力输出。

结果为 1 帧、8 原子,最大力分量约 0.217059 eV/Å。该静态脚本使用自身的 Ha 换算常数,能量约 −846.705405 eV/晶胞。不同脚本的常数精度可能略有差别,批量数据应统一转换路线,不能把这种数值差当成结构物性变化。

3. 准备 MD 的三个文件#

MD 在上述基础上使用 RUN_TYPE MD,并增加 MOTION/MD 的 NVE、4 步和 0.5 fs;逐步输出坐标、力与晶胞。完整输入在 inputs/cp2k-md/input.inp。

从资料包根目录准备独立副本:

mkdir my-cp2k-md
cp inputs/cp2k-md/labels-pos-1.xyz inputs/cp2k-md/labels-frc-1.xyz \
  inputs/cp2k-md/labels-1.cell inputs/cp2k-md/output.log my-cp2k-md/
cp check_cp2k_md.py label_io.py add_cp2k_virial.py my-cp2k-md/
cd my-cp2k-md
python check_cp2k_md.py labels-pos-1.xyz labels-frc-1.xyz labels-1.cell

检查内容包括帧数、原子顺序、步号、时间和晶胞。不能只截短一个文件去凑相同长度;那会把不同构型的力或晶胞配在一起。

4. 使用真正的 MD 入口#

gpumdkit.sh

依次输入 1 格式转换 → 103 CP2K → 2 三文件路线,随后输入:

labels-pos-1.xyz labels-frc-1.xyz labels-1.cell -shifted no

完成后用 0 返回、再用 0 退出。也可直接调用对应原生脚本:

python "$GPUMDkit_path/Scripts/format_conversion/cp2k2xyz.py" \
  labels-pos-1.xyz labels-frc-1.xyz labels-1.cell -shifted no
head -n 4 original.xyz

-shifted no 保留原始总能量,不做原子能量零点拟合。本例得到 5 帧、0–2 fs;力从 Hartree/Bohr 换为 eV/Å,能量从 Hartree 换为 eV。

5. 需要 virial 时,接上对应应力日志#

三文件原生输出包含 E、F、坐标与晶胞,没有自动带出 virial。本例计算日志同时逐帧打印了解析应力,可以继续运行:

python add_cp2k_virial.py original.xyz output.log with-virial.xyz
head -n 4 with-virial.xyz

脚本核对能量顺序、应力块和构型数,再把 CP2K 的 bar 换为 GPa,并乘各帧体积得到 eV/晶胞的 virial。没有对应应力时,不用零矩阵代替缺失标签。

6. 本例最终得到什么#

静态得到 1 个标签构型,MD 得到 5 个标签构型;完整操作得到逐帧能量、力、晶胞,MD 的应力补充路线还提供 virial。数据目录和文件原件均保留。

bash run.sh all

在资料包根目录可一次重现转换、帧对应和图示。继续做 DeepMD 格式衔接前,见标签规范化与字段保留。

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

查看配图

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

100%查看原图 ↗