CREST分子动力学实操
从XYZ、完整初始坐标速度和TOML开始,运行水分子短轨迹,读取结构、能量与时间序列。
阅读本篇目录
本篇用一个孤立水分子,演示 CREST 3.0.2 的原生 MD 输入和结果读取。计算采用 GFN2/tblite,步长 0.25 fs,总长 100 fs;这是熟悉命令和输出的小例,不用于液态水性质或收敛热力学统计。
1. 下载、提交、查看结果#
下载完整输入、参考数据与绘图包,解压并进入 crest-md-mtd。用 sbatch job.slurm 提交;CPU 分区按56核整节点分配和计费,脚本最多并发12个单线程小任务。在已分配的CPU作业中,也可以运行:
bash run.sh
脚本依次准备连续MD、短MD、续算和RMSD偏置MD,并做各快照的独立电子能与梯度计算。原生输出在 cases/ 和 snapshot-cases/,CSV在 results/,本次新图在 plots/。软件使用现有 crest/3.0.2;Python与绘图环境由 sisl/0.16.4 提供。
2. 准备XYZ和完整TOML#
inputs/md-full/structure.xyz:
3
Etot= -5.0705444504 g norm= 0.00000890
O -0.0332542193 -0.0429639078 0.0000000000
H 0.9257005935 -0.0656222349 0.0000000000
H -0.2956160149 0.8796812139 0.0000000000
inputs/md-full/input.toml:
input = "structure.xyz"
runtype = "md"
threads = 1
[[calculation.level]]
method = "gfn2"
chrg = 0
uhf = 0
accuracy = 0.0001
[dynamics]
length = 0.1
tstep = 0.25
dump = 0.25
temp = 300.0
hmass = 1.00794075
shake = 0
printstep = 10
restart = "start.mdrestart"
length 的单位是 ps,所以0.1表示100 fs;tstep 和 dump 是 fs。temp 为核运动恒温目标300 K,电子温度另由本版本电子结构默认设置控制。shake=0 不约束O–H键,因而轨迹中能看到内部振动。
start.mdrestart 是包内已有的完整初始状态,第一行为时间,后面每个原子一行:前三列为Bohr坐标,后三列为原子单位速度。此小例保留9个笛卡尔自由度,未去除整体平动和转动。不要把这套固定初速度直接套到其他分子或凝聚态。
自己准备其他体系时,在新目录放好结构、对应的完整初态和TOML,加载模块后运行 crest input.toml > native.out 2>&1。原固定参考只适用于本资料包。
3. 检查原生输出与轨迹#
grep -E 'normal MD termination|MD run completed successfully' cases/md-full/native.out
head -n 10 cases/md-full/crest_dynamics.trj
crest_dynamics.trj 是连续XYZ,每帧注释包含 Epot。这条普通MD中它是电子势能,不是包含核运动动能的总能。程序使用恒温器,不能把它当作严格NVE守恒测试。
本版本400步产生399个轨迹记录,坐标时间为0.25–99.75 fs;最终 .mdrestart 则到100 fs。两者不是缺一步计算,而是轨迹保存位置与最终checkpoint的约定不同。
4. 提取键长、角度和能量#
head -n 6 results/md-timeseries.csv
python -B verify.py
CSV列出每帧的坐标时间、原生势能、独立电子能、偏置能、两条O–H键长和H–O–H角。普通MD这次的O–H范围为 0.92937–0.99223 Å,角度约 101.85–112.45°。

图中结构为实际0.25 fs快照,能量以各轨迹第一帧为零点。reference-plots/ 保留本教程示例,plots/ 是当前新运行生成的图。接着可以看完整MD续算和RMSD偏置。输入格式与具体版本约定见官方TOML文档及3.0.2动力学实现。