GPUMD 控压与晶胞形状

控压之前先决定允许晶胞怎样变化：只改体积、三个长度分别变化，还是连倾斜分量一起变化。对二维材料或固定某一方向的任务，不能把三维各向同性控压直接照搬。本篇用三维Cu、Ar展示本版本的常规控压算法和晶胞选项。

1. 从一个各向同性Cu任务开始#

在已有计算资源中，完成一次共同准备后：

bash run.sh cu-npt-scr-iso
tail -n 1 work/cu-npt-scr-iso/thermo.out
head -n 2 work/cu-npt-scr-iso/restart.xyz

cu-npt-scr-iso / run.in

potential potential.txt
time_step 1
ensemble npt_scr 300 300 100 0 100 1000
dump_thermo 100
dump_xyz 1000 trajectory.xyz precision double mass velocity force potential virial unwrapped_position
dump_restart 20000
run 20000

这里300是温度K，100是温度耦合步数，0是目标压力GPa，后面的100是用于控压数值尺度的模量估计GPa，1000是压力耦合步数。模量估计不是本次预测出的弹性常数，也不是应变率。

1bar约为0.0001GPa，不能把1bar写成1GPa。Ar的标准控压例子使用约1GPa模量估计、2fs时间步；Cu用约100GPa、1fs。选这些数值用于演示稳定操作，不是材料参数推荐表。

2. Berendsen与随机晶胞重标度#

将 npt_scr 换成 npt_ber 是另一种算法，不代表相同系综采样。两者都提供以下形状设置：

三斜例子提供了实际倾斜的初始晶胞，并对原坐标做同样的仿射变换；不是只改Lattice数字而保持坐标不动。三个对角分量和三斜倾角的实际变化见下图。

potential potential.txt
time_step 1
ensemble npt_ber 300 300 100 0 100 1000
dump_thermo 100
dump_xyz 1000 trajectory.xyz precision double mass velocity force potential virial unwrapped_position
dump_restart 20000
run 20000

potential potential.txt
time_step 1
ensemble npt_ber 300 300 100 0 0 0 100 100 100 1000
dump_thermo 100
dump_xyz 1000 trajectory.xyz precision double mass velocity force potential virial unwrapped_position
dump_restart 20000
run 20000

potential potential.txt
time_step 1
ensemble npt_ber 300 300 100 0 0 0 0 0 0 100 100 100 100 100 100 1000
dump_thermo 100
dump_xyz 1000 trajectory.xyz precision double mass velocity force potential virial unwrapped_position
dump_restart 20000
run 20000

potential potential.txt
time_step 1
ensemble npt_scr 300 300 100 0 0 0 100 100 100 1000
dump_thermo 100
dump_xyz 1000 trajectory.xyz precision double mass velocity force potential virial unwrapped_position
dump_restart 20000
run 20000

potential potential.txt
time_step 1
ensemble npt_scr 300 300 100 0 0 0 0 0 0 100 100 100 100 100 100 1000
dump_thermo 100
dump_xyz 1000 trajectory.xyz precision double mass velocity force potential virial unwrapped_position
dump_restart 20000
run 20000

3. MTTK使用另一套关键词语法#

cu-npt-mttk-iso / run.in

potential potential.txt
time_step 1
ensemble npt_mttk temp 300 300 iso 0 0 tperiod 100 pperiod 2000
dump_thermo 100
dump_xyz 1000 trajectory.xyz precision double mass velocity force potential virial unwrapped_position
dump_restart 20000
run 20000

iso 是等比例改体积；aniso 允许三个长度独立变化；tri 允许六个晶胞自由度变化。也可只控制一个方向：

cu-npt-mttk-x / run.in

potential potential.txt
time_step 1
ensemble npt_mttk temp 300 300 x 0 0 tperiod 100 pperiod 2000
dump_thermo 100
dump_xyz 1000 trajectory.xyz precision double mass velocity force potential virial unwrapped_position
dump_restart 20000
run 20000

这个x方向例子中，另两个晶胞方向保持原值，不能把它当成三个方向都达到零压力。pperiod 是按时间步计的周期参数，本版要求至少200步；本例使用2000步。

potential potential.txt
time_step 1
ensemble npt_mttk temp 300 300 aniso 0 0 tperiod 100 pperiod 2000
dump_thermo 100
dump_xyz 1000 trajectory.xyz precision double mass velocity force potential virial unwrapped_position
dump_restart 20000
run 20000

potential potential.txt
time_step 1
ensemble npt_mttk temp 300 300 tri 0 0 tperiod 100 pperiod 2000
dump_thermo 100
dump_xyz 1000 trajectory.xyz precision double mass velocity force potential virial unwrapped_position
dump_restart 20000
run 20000

4. NPH不再控制温度#

cu-nph-mttk / run.in

potential potential.txt
time_step 1
ensemble nph_mttk iso 0 0 pperiod 2000
dump_thermo 100
dump_xyz 1000 trajectory.xyz precision double mass velocity force potential virial unwrapped_position
dump_restart 20000
run 20000

NPH保留控压与对应的扩展动力学，没有恒温器。温度会随过程变化；判断守恒时不能把18列输出里的物理动能加势能，直接当成包含恒压器自由度的完整扩展守恒量。

5. 读取和使用实际晶胞#

thermo.out 最后九列按 ax ay az bx by bz cx cy cz 给出晶胞矢量。体积要用3×3行列式，三斜盒不能直接把三个对角数字相乘代替一般公式：

module load ase/3.29.0
python - <<'PYBOX'
import numpy as np
a=np.loadtxt('work/cu-npt-mttk-tri/thermo.out')
cell=a[:,9:].reshape(-1,3,3)
volume=np.linalg.det(cell)
print('first/last volume A^3:',volume[0],volume[-1])
PYBOX

本例三斜矩阵恰为下三角，但分析仍按一般行列式处理。液体缺少长时剪切刚度，随意放开所有形状自由度可能产生极端形变，应按物理问题选择约束；不要因为本例20ps运行结束就默认任意液体三斜长跑都合适。QTB控压的温度含义见QTB专题。

资料包与完整运行#

下载完整输入、真实参考数据、末态重算与绘图脚本，解压后进入目录。先重画随包结果：

bash analyze.sh

输出为 reference/analysis/ 下的CSV、数值摘要和 reference/figures/ 中的图片。原生末态重算结果已经随包，分析会先核对其输入，再处理数据。

新解压目录提交准备加一个Cu NVT-BDP实例：

sbatch run.slurm

提交完整46条件、两个材料准备流程和依赖分析：

bash submit.sh

只生成脚本检查，用 bash submit.sh --dry-run。整套完成后：

bash analyze.sh work common

这一步对每个重启末态做零步静态导出，再独立分析和画图，应在已分配的计算资源内执行；生成的分析作业已包含这些命令。只运行单例时，先直接查看其 thermo.out、trajectory.xyz、restart.xyz，统一分析需要46条件齐全。

run.slurm 默认一张P100、一个CPU核、gpu-8p100-hs；套餐用户按权限改为例如 gpu-8p100-tc。整套提交还需修改 submit.py 中的partition。bash prepare.sh 和 bash run.sh 条件名 是已有任务资源内的入口，不在登录或管理节点直接运行。prepare.sh 只在新目录建立 common/cu、common/ar，已有common/work/jobs保留；重跑整套请重新解压。单进程启动，不使用mpirun。

固定环境：GPUMD5.9.1、GCC13.3、CUDA12.9，模块 gpumd/5.9.1-p100；该构建含PLUMED2.10.1和NetCDF-C4.10.1，本批没有施加偏置。分析使用ASE3.29.0环境中的NumPy/Matplotlib。20ps结果用于操作与数值关系参考，正式材料研究仍需检查准备时长、模型、尺寸、时间步和采样；不据此给出收敛热力学量、相界或性能排名。