GPUMDkit 多温度扩散与采样检查

下载输入、实际结果与绘图脚本。以下命令在解压后的目录执行。

module load gpumdkit/1.5.8
bash run.sh
cat LAST-RUN.txt

run.sh在新的work目录中复跑工具，不会重新启动DFT或分子动力学。下面也列出每一步的单独命令。GPUMDkit为1.5.8、Python3.11；GPUMD/NEP为5.9.1，实际GPU计算使用P100单卡。

1. 建立同一套温度目录#

本例为512原子Si，Tersoff Si_Fan_2019模型，固定立方盒长21.724 Å。五个温度为4000、4500、5000、5500、6000 K，各一个种子71–75。先在6000 K运行10 ps，再用10 ps变到目标温度，平衡10 ps，最后采样50 ps。这里只展示高温模型数据处理，不据此给真实Si熔点或室温扩散率。

4000 K实际完整输入如下；换温度时同时改温度与独立随机种子：

potential potential.txt
velocity 6000 seed 71
time_step 0.5
ensemble nvt_nhc 6000 6000 100
run 20000
ensemble nvt_nhc 6000 4000 100
run 20000
ensemble nvt_nhc 4000 4000 100
run 20000
ensemble nvt_nhc 4000 4000 100
compute_msd 10 1000 all_groups 0 save_every 20000
compute_rdf 8 160 500
dump_thermo 1000
dump_xyz 1000 trajectory.xyz
run 100000

步长0.5 fs。compute_msd 10 1000每10步采样，最大相关时间5 ps；save_every 20000每累计10 ps另存一份MSD。all_groups 0按第0套分组输出；本例所有Si都属于组0，因此只有一个组，不是假造多元素数据。

2. 提交、检查与准备绘图#

mkdir md-new
cp reference/md/diffusion-4000/{model.xyz,potential.txt,run.in} md-new/
cd md-new
module load gpumd/5.9.1-p100
gpumd > gpumd.log
tail gpumd.log
ls msd*.out thermo.out
cd ..

也可以在附件根目录直接执行sbatch compute.sbatch diffusion-4000。上述gpumd命令在分配的P100资源内运行；先用附件参考结果练习后处理也可以。默认bash run.sh已把五组目录整理成work/.../diffusion/4000K等结构。进入该diffusion目录：

gpumdkit.sh -plt D save
gpumdkit.sh -plt D_xyz save
cd 4000K
gpumdkit.sh -plt msd_all msd.out Si save
gpumdkit.sh -plt msd_conv save
cd ..

MSD文件首列为相关时间ps，接下来三列为xyz方向MSD，单位Ų；单组后面三列为方向SDC。msd_all后面的Si是组标签，与真实组顺序一致才有意义。

3. 读取实际扩散系数#

工具取MSD相关时间序列40%–80%区间线性拟合；本例约2–4 ps。三维总MSD斜率除以6得到Ų/ps，再乘1E-4得到cm²/s。单方向斜率则除以2。

五点Arrhenius拟合的表观斜率对应 0.838 eV，R²=0.9975。这是这组势模型、固定密度和温区内的结果；高R²并不能消除小盒子、有限采样及势模型误差，也不能据此外推到300 K。

4. 判断是否还需要更长采样#

4000 K累计10、20、30、40、50 ps得到D分别为2.3423、2.4987、2.4368、2.4387、2.4732 Ų/ps。msd_conv横轴来自文件名中的累计步数乘步长，不能把同一个MSD随意改文件名冒充不同长度计算。

这些累计窗口共享轨迹，不能当作五次独立重复来计算置信区间。客户研究应增加独立种子、体系大小和采样长度；若斜率为负、MSD仍平台化，不能强行做正扩散拟合。

5. 分段拟合入口怎么用#

gpumdkit.sh -plt D_PT 5000 save

此命令只是演示在5000 K两侧分别拟合；本例并没有发现5000 K相变。工具把分界点同时纳入两侧：低温三点得到0.820 eV，高温三点0.867 eV。真实研究中必须先用结构或其他物性证据确定分界温度，不能从这张拟合图反推相变。

返回GPUMDkit专题