CP2K专题

CP2K分子动力学系综与时间步

Ar108小模型的NVE、NVT、NPT_I输入与1fs/2fs实际对照。

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

阅读本篇目录

下载完整输入、实际输出、续算文件和五张高清图。在已分配CPU资源内，新解压目录执行 bash run.sh，计算写入 work/cases，分析写入 work/analysis，图片在 work/figures。

1 先准备同一个参考结构#

本例为108个Ar的FCC初始结构，常规胞晶格常数5.26Å，重复3×3×3，盒长15.78Å。LJ截断7.5Å，小于半盒长。初始结构可用ASE重建：

from ase.build import bulk
from ase.io import write
atoms = bulk('Ar', 'fcc', a=5.26, cubic=True) * (3,3,3)
write('ar108.extxyz', atoms)

Ar108初始结构与1ps NPT末态

实际MD输入已包含全部108行坐标：

cat inputs/ar-nve-2fs/input.inp
bash run.sh
ls work/analysis

MD使用同一势函数、初始结构和随机种子20261028。初始温度50K，所有完整轨迹时长1ps：2fs×500步，或1fs×1000步。短轨迹用于检查操作和数值，不作为平衡物性采样。

&MOTION
 &MD
  ENSEMBLE NVE
  STEPS 500
  TIMESTEP [fs] 2
  TEMPERATURE 50
 &END MD
&END MOTION

这里的50K用于初始化速度，之后温度可以变化。对这个初始晶体，动能和势能重新分配，NVE温度一度降到约9.25K；不能因为它没有保持50K，就认定积分程序有问题。应结合总能量守恒判断。

在MD段中改为 ENSEMBLE NVT 并加入：

&THERMOSTAT
 TYPE NOSE
 &NOSE
  TIMECON 100
 &END NOSE
&END THERMOSTAT

本例Nose时间常数100fs，目标50K。恒温器控制的是相应统计分布，不要求每一步温度都等于50K；热化时间、时间常数和系统大小都会影响曲线。

MD段使用 ENSEMBLE NPT_I、同一恒温器，再加：

&BAROSTAT
 PRESSURE [bar] 1.0
 TIMECON 500
&END BAROSTAT

FORCE_EVAL段还需要请求应力：

&FORCE_EVAL
 METHOD FIST
 STRESS_TENSOR ANALYTICAL
 ...
&END FORCE_EVAL

上面省略处是完整输入中已有的MM和SUBSYS段，并非可直接保留的输入字符。NPT依赖维里/应力，不能只写BAROSTAT而不请求相应应力计算。

本例体积从3929.352552Å³变到1ps末态的4402.852026Å³。这是短时间的体积响应，不能当作已经收敛的1bar平衡密度。各向异性晶胞、表面和二维体系需要另选合适的晶胞自由度，不直接套用NPT_I。

head work/analysis/ar-nve-2fs-thermo.csv
head work/analysis/ar-nve-1fs-thermo.csv

NVE设置	总时长	守恒量峰峰变化/eV·atom⁻¹
2fs×500步	1ps	4.89553×10⁻⁷
1fs×1000步	1ps	1.47647×10⁻⁷

时间步、能量守恒和三种系综的实际温度

这个小例减小时间步后能量变化更小。NVE检查物理总能量；NVT/NPT的原生守恒量还包含恒温/恒压扩展自由度，不能把它们与只取动能加势能的列混同。

较大时间步、轻原子、高频键、刚性约束和反应势会改变稳定性要求；不能因为Ar的2fs可用，就套给水分子的高频伸缩。后处理和断点续算见下一篇。

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