ABACUS 六种恒温器实操
NHC、Langevin、Anderson、Berendsen及两种重标方法的完整输入和实测。
阅读本篇目录
下载13个条件的完整输入、原生轨迹、续算状态与绘图脚本，解压进入附件根目录。先用已有数据生成图表：
bash analyze.sh
结果在reference/analysis、图片在reference/figures、ASE轨迹在reference/export。从头执行整套流程：
sbatch submit.slurm
脚本在新的work中依次运行13个条件，已有work时请新解压，避免覆盖已完成结果。想先跑一个条件，可在计算资源环境中加载下述模块、运行setup.py后用run_case.py指定目录名；普通用户不要在登录节点运行DFT。提交模板申请1张P100、1个MPI进程、1个线程，套餐用户按实际权限用yskit生成脚本。
本例ABACUS3.10.1模块为abacus/3.10.1-p100，GNU13.3、HPC-X2.50/OpenMPI5.0.10rc2、CUDA12.8环境，实际为P100双精度；后处理为ase/3.29.0。模型采用32原子fcc铝、3价电子LDA赝势Al.pz-vbc.UPF及PZ泛函，40Ry、Γ点、64条带，共96个价电子。不能用同名11价电子赝势替换后仍保留这些带数与精度设置。
这批轨迹只有40fs，用来指导输入、读取和续算，不用于扩散系数、熔点、平衡密度或收敛热力学性质。
1 相同初态，分别设置恒温算法#
先按NVE教程准备结构、质量和k点。本例六种恒温条件全部使用同一份真实初始速度与300K目标，单步0.5fs、共80步；相同温度目标并不保证40fs内已达到平衡。
Nose–Hoover链的完整INPUT：
INPUT_PARAMETERS
suffix tutorial
pseudo_dir ../../PP_ORB
ntype 1
basis_type pw
calculation md
ecutwfc 40
nbands 64
ks_solver cg
symmetry 0
device gpu
precision double
scf_thr 1e-7
scf_nmax 100
smearing_method gaussian
smearing_sigma 0.01
cal_force 1
cal_stress 1
out_chg 0
md_type nvt
md_nstep 80
md_dt 0.5
md_tfirst 300
md_seed 12345
md_dumpfreq 1
md_restartfreq 40
dft_functional pz
init_vel 1
md_tlast 300
dump_force 1
dump_vel 1
dump_virial 1
md_thermostat nhc
md_tfreq 0.01
md_tchain 3
md_type nvt配合md_thermostat nhc，md_tchain 3为链长度，md_tfreq 0.01（fs⁻¹）调节温控耦合。它不是随机碰撞恒温器。NHC的物理势能加动能可以与链交换，不要求它像NVE一样恒定。
2 选对参数，不把不同重标方法混成一种#
| 目录 | 主要设置 | 本例参数含义 |
|---|---|---|
| al-nhc | nvt / nhc | Nose–Hoover链，tfreq=0.01、tchain=3 |
| al-langevin | langevin | 随机阻尼热浴，md_damp=50fs |
| al-anderson | nvt / anderson | 每步每个原子的碰撞概率1/md_nraise，本例1/10 |
| al-berendsen | nvt / berendsen | 按当前与目标温度弱耦合，本例md_nraise=10 |
| al-rescaling | nvt / rescaling | 偏差超过md_tolerance=10K时直接温度重标 |
| al-rescale_v | nvt / rescale_v | 每md_nraise=10步做速度重标 |
rescaling与rescale_v的触发条件不同。程序参数拼写为anderson，请按原生关键字输入。其余五个完整输入如下，结构/KPT与NVE一致：
al-langevin完整INPUT
INPUT_PARAMETERS
suffix tutorial
pseudo_dir ../../PP_ORB
ntype 1
basis_type pw
calculation md
ecutwfc 40
nbands 64
ks_solver cg
symmetry 0
device gpu
precision double
scf_thr 1e-7
scf_nmax 100
smearing_method gaussian
smearing_sigma 0.01
cal_force 1
cal_stress 1
out_chg 0
md_type langevin
md_nstep 80
md_dt 0.5
md_tfirst 300
md_seed 12345
md_dumpfreq 1
md_restartfreq 40
dft_functional pz
init_vel 1
md_tlast 300
dump_force 1
dump_vel 1
dump_virial 1
md_damp 50
al-anderson完整INPUT
INPUT_PARAMETERS
suffix tutorial
pseudo_dir ../../PP_ORB
ntype 1
basis_type pw
calculation md
ecutwfc 40
nbands 64
ks_solver cg
symmetry 0
device gpu
precision double
scf_thr 1e-7
scf_nmax 100
smearing_method gaussian
smearing_sigma 0.01
cal_force 1
cal_stress 1
out_chg 0
md_type nvt
md_nstep 80
md_dt 0.5
md_tfirst 300
md_seed 12345
md_dumpfreq 1
md_restartfreq 40
dft_functional pz
init_vel 1
md_tlast 300
dump_force 1
dump_vel 1
dump_virial 1
md_thermostat anderson
md_nraise 10
md_tolerance 10
al-berendsen完整INPUT
INPUT_PARAMETERS
suffix tutorial
pseudo_dir ../../PP_ORB
ntype 1
basis_type pw
calculation md
ecutwfc 40
nbands 64
ks_solver cg
symmetry 0
device gpu
precision double
scf_thr 1e-7
scf_nmax 100
smearing_method gaussian
smearing_sigma 0.01
cal_force 1
cal_stress 1
out_chg 0
md_type nvt
md_nstep 80
md_dt 0.5
md_tfirst 300
md_seed 12345
md_dumpfreq 1
md_restartfreq 40
dft_functional pz
init_vel 1
md_tlast 300
dump_force 1
dump_vel 1
dump_virial 1
md_thermostat berendsen
md_nraise 10
md_tolerance 10
al-rescaling完整INPUT
INPUT_PARAMETERS
suffix tutorial
pseudo_dir ../../PP_ORB
ntype 1
basis_type pw
calculation md
ecutwfc 40
nbands 64
ks_solver cg
symmetry 0
device gpu
precision double
scf_thr 1e-7
scf_nmax 100
smearing_method gaussian
smearing_sigma 0.01
cal_force 1
cal_stress 1
out_chg 0
md_type nvt
md_nstep 80
md_dt 0.5
md_tfirst 300
md_seed 12345
md_dumpfreq 1
md_restartfreq 40
dft_functional pz
init_vel 1
md_tlast 300
dump_force 1
dump_vel 1
dump_virial 1
md_thermostat rescaling
md_nraise 10
md_tolerance 10
al-rescale_v完整INPUT
INPUT_PARAMETERS
suffix tutorial
pseudo_dir ../../PP_ORB
ntype 1
basis_type pw
calculation md
ecutwfc 40
nbands 64
ks_solver cg
symmetry 0
device gpu
precision double
scf_thr 1e-7
scf_nmax 100
smearing_method gaussian
smearing_sigma 0.01
cal_force 1
cal_stress 1
out_chg 0
md_type nvt
md_nstep 80
md_dt 0.5
md_tfirst 300
md_seed 12345
md_dumpfreq 1
md_restartfreq 40
dft_functional pz
init_vel 1
md_tlast 300
dump_force 1
dump_vel 1
dump_virial 1
md_thermostat rescale_v
md_nraise 10
md_tolerance 10
3 用实际轨迹观察响应#
| 条件 | 40fs平均温度 / K | 末温 / K |
|---|---|---|
| al-nhc | 209.039 | 72.835 |
| al-langevin | 278.859 | 200.155 |
| al-anderson | 293.566 | 242.022 |
| al-berendsen | 269.556 | 236.030 |
| al-rescaling | 296.050 | 291.183 |
| al-rescale_v | 280.612 | 300.000 |
这些数值是短暂响应，不是统计收敛的温控优劣排名。弱耦合或直接重标不保证严格正则系综的涨落分布。强制重标更接近300K，也可能改变真实动力学关联；不要因曲线平就直接拿来算扩散、振动态密度或输运。先确定目标系综与物性，再选择耦合参数和生产段。
python - <<'PY'
import csv
names={'al-nhc','al-langevin','al-anderson','al-berendsen','al-rescaling','al-rescale_v'}
with open('reference/analysis/conditions.csv') as f:
for r in csv.DictReader(f):
if r['name'] in names:
print(r['name'],r['mean_T_K'],r['final_T_K'])
PY
4 检查自由度与速度单位#
本例32原子、无位置冻结，原生温度按3N−3=93个自由度计算。初态已去除质心漂移；随机浴运行后仍应检查质心速度，分母减3不代表每一步都再次清零质心动量。附件使用原生全动能核对原生温度，不擅自删掉漂移能。输出MD_dump的速度是Å/fs；不能直接当作STRU里原子单位的速度重新写回。附件用真实质量和原生速度独立重算动能、温度，并与热力学表核对。